Output reliability curves for each machine learning model, size threshold, and time window.

        :param prob_thresholds:
        :param dilation_radius:
        :return:
359652
def reliability_curves(self, prob_thresholds):
        
        all_rel_curves = {}
        for model_name in self.model_names:
            all_rel_curves[model_name] = {}
            for size_threshold in self.size_thresholds:
                all_rel_curves[model_name][size_threshold] = {}
                for h, hour_window in enumerate(self.hour_windows):
                    hour_range = (hour_window.start, hour_window.stop)
                    all_rel_curves[model_name][size_threshold][hour_range] = \
                        DistributedReliability(prob_thresholds, 1)
                    if self.obs_mask:
                        all_rel_curves[model_name][size_threshold][hour_range].update(
                            self.window_forecasts[model_name][size_threshold][h][
                                self.window_obs[self.mask_variable][h] > 0],
                            self.dilated_obs[size_threshold][h][self.window_obs[self.mask_variable][h] > 0]
                        )
                    else:
                        all_rel_curves[model_name][size_threshold][hour_range].update(
                            self.window_forecasts[model_name][size_threshold][h],
                            self.dilated_obs[size_threshold][h]
                        )
        return all_rel_curves
688655
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        # extracting dictionary of coefficients specific to required
        # intensity measure type.
        C = self.COEFFS[imt]

        mean = (self._get_magnitude_scaling_term(C, rup.mag) +
                self._get_distance_scaling_term(C, rup.mag, dists.rrup) +
                self._get_style_of_faulting_term(C, rup.rake) +
                self._get_site_scaling_term(C, sites.vs30))

        stddevs = self._get_stddevs(imt,
                                    rup.mag,
                                    len(dists.rrup),
                                    stddev_types)
        return mean, stddevs
538394
def precision_curve(self, delta_tau=0.01):
        
        # compute thresholds based on the sorted probabilities
        orig_thresh = self.threshold
        sorted_labels, sorted_probs = self.sorted_values

        scores = []
        taus = []
        tau = 0
        for k in range(len(sorted_labels)):
            # compute new accuracy
            self.threshold = tau
            scores.append(self.precision)
            taus.append(tau)

            # update threshold
            tau = sorted_probs[k]

        # add last datapoint
        tau = 1.0
        self.threshold = tau
        scores.append(1.0)
        taus.append(tau)

        self.threshold = orig_thresh
        return scores, taus
527754
def summary(self):
        
        ci = 1 - self.alpha
        lower_upper_bounds = self._compute_confidence_bounds_of_parameters()
        df = pd.DataFrame(index=self._fitted_parameter_names)
        df["coef"] = self._fitted_parameters_
        df["se(coef)"] = self._compute_standard_errors().loc["se"]
        df["lower %g" % ci] = lower_upper_bounds["lower-bound"]
        df["upper %g" % ci] = lower_upper_bounds["upper-bound"]
        df["p"] = self._compute_p_values()
        with np.errstate(invalid="ignore", divide="ignore"):
            df["-log2(p)"] = -np.log2(df["p"])
        return df
290469
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        # extracting dictionary of coefficients specific to required
        # intensity measure type.
        C = self.COEFFS[imt]
        mean = (self._compute_style_of_faulting_term(rup, C) +
                self._compute_magnitude_scaling(rup.mag, C) +
                self._compute_distance_scaling(dists.rjb, C) +
                self._compute_site_term(sites.vs30, C))

        stddevs = self._get_stddevs(C, stddev_types, num_sites=len(sites.vs30))

        return mean, stddevs
778281
def convenience_calc_probs(self, params):
        
        orig_nest_coefs, betas = self.convenience_split_params(params)
        natural_nest_coefs = nc.naturalize_nest_coefs(orig_nest_coefs)

        args = [natural_nest_coefs,
                betas,
                self.design,
                self.rows_to_obs,
                self.rows_to_nests]
        kwargs = {"chosen_row_to_obs": self.chosen_row_to_obs,
                  "return_type": "long_and_chosen_probs"}

        probability_results = general_calc_probabilities(*args, **kwargs)

        return probability_results
313181
def bounds(self):
        
        return [scipy.stats.lognorm.interval(self.i, s, loc=0, scale=em) for s, em in zip(self.sigma, self.emu)]
355971
def filter_threshold(df, inst_rc, threshold, num_occur=1):
  '''
  
  '''
  from copy import deepcopy

  inst_df = deepcopy(df['mat'])

  if inst_rc == 'col':
    inst_df = inst_df.transpose()

  inst_df = inst_df.abs()

  ini_rows = inst_df.index.values.tolist()

  inst_df[inst_df < threshold] = 0
  inst_df[inst_df >= threshold] = 1

  tmp_sum = inst_df.sum(axis=1)

  tmp_sum = tmp_sum[tmp_sum >= num_occur]

  keep_names = tmp_sum.index.values.tolist()

  if inst_rc == 'row':
    if len(keep_names) < len(ini_rows):
      df['mat'] = grab_df_subset(df['mat'], keep_rows=keep_names)

      if 'mat_up' in df:
        df['mat_up'] = grab_df_subset(df['mat_up'], keep_rows=keep_names)
        df['mat_dn'] = grab_df_subset(df['mat_dn'], keep_rows=keep_names)

      if 'mat_orig' in df:
        df['mat_orig'] = grab_df_subset(df['mat_orig'], keep_rows=keep_names)

  elif inst_rc == 'col':
    inst_df = inst_df.transpose()

    inst_rows = inst_df.index.values.tolist()
    inst_cols = keep_names

    df['mat'] = grab_df_subset(df['mat'], inst_rows, inst_cols)

    if 'mat_up' in df:
      df['mat_up'] = grab_df_subset(df['mat_up'], inst_rows, inst_cols)
      df['mat_dn'] = grab_df_subset(df['mat_dn'], inst_rows, inst_cols)

    if 'mat_orig' in df:
      df['mat_orig'] = grab_df_subset(df['mat_orig'], inst_rows, inst_cols)

  return df
289781
def get_thresholds(points=100, power=3) -> list:
    """"""
    return [(i / (points + 1)) ** power for i in range(1, points + 1)]
619486
def runExperimentS(dirName):
  
  # Results are put into a pkl file which can be used to generate the plots.
  # dirName is the absolute path where the pkl file will be placed.
  resultsFilename = os.path.join(dirName, "superimposed_training.pkl")
  results = runExperiment(
    {
      "numSequences": 50,
      "numObjects": 50,
      "seqLength": 10,
      "numFeatures": 100,
      "trialNum": 8,
      "numLocations": 100,
      "numRepetitions": 30,
      "sampleSizeProximalL2": 15,
      "minThresholdProximalL2": 10,
      "figure": "S",
      "stripStats": False,
    }
  )

  # Pickle results for plotting and possible later debugging
  with open(resultsFilename, "wb") as f:
    cPickle.dump(results, f)

  # Debugging
  with open(resultsFilename, "rb") as f:
    r = cPickle.load(f)

    r.pop("objects")
    r.pop("sequences")
    stat = r.pop("statistics")
    pprint.pprint(r)
    sObject = 0
    sSequence = 0
    for i in range(0, 50):
      sObject += sum(stat[i]['L4 PredictedActive C0'])
    for i in range(50, 100):
      sSequence += sum(stat[i]['L4 PredictedActive C0'])
    print sObject, sSequence
121494
def scores_to_preds(self, threshold, use_probs = True):
        
        self.threshold = threshold

        if use_probs:
            if self.probs is None:
                raise DataError("Probabilities are not available to make "
                                "predictions.")
            else:
                word = "probabilities"
                scores = self.probs
        else:
            if self.scores is None:
                raise DataError("Scores are not available to make predictions.")
            else:
                word = "scores"
                scores = self.scores

        if threshold > np.max(scores) or threshold < np.min(scores):
            warnings.warn("Threshold {} is outside the range of the "
                          "{}.".format(self.threshold, word))

        if self.preds is not None:
            warnings.warn("Overwriting predictions")
        self.preds = (scores >= threshold)*1
281540
def _thresholds_init(self):
        
        thresholds = getattr(self._py3status_module, "thresholds", [])
        self._thresholds = {}
        if isinstance(thresholds, list):
            try:
                thresholds.sort()
            except TypeError:
                pass
            self._thresholds[None] = [(x[0], self._get_color(x[1])) for x in thresholds]
            return
        elif isinstance(thresholds, dict):
            for key, value in thresholds.items():
                if isinstance(value, list):
                    try:
                        value.sort()
                    except TypeError:
                        pass
                    self._thresholds[key] = [
                        (x[0], self._get_color(x[1])) for x in value
                    ]
490685
def precision_pct_pred_pos_curve(self, interval=False, delta_tau=0.001):
        
        # compute thresholds based on the sorted probabilities
        orig_thresh = self.threshold
        sorted_labels, sorted_probs = self.sorted_values

        precisions = []
        pct_pred_pos = []
        taus = []
        tau = 0
        if not interval:
            for k in range(len(sorted_labels)):
                # compute new accuracy
                self.threshold = tau
                precisions.append(self.precision)
                pct_pred_pos.append(self.pct_pred_pos)
                taus.append(tau)
                
                # update threshold
                tau = sorted_probs[k]

        else:
            while tau < 1.0:
                # compute new accuracy
                self.threshold = tau
                precisions.append(self.precision)
                pct_pred_pos.append(self.pct_pred_pos)
                taus.append(tau)
                
                # update threshold
                tau += delta_tau

        # add last datapoint
        tau = 1.0
        self.threshold = tau
        precisions.append(self.precision)
        pct_pred_pos.append(self.pct_pred_pos)
        taus.append(tau)

        precisions.append(1.0)
        pct_pred_pos.append(0.0)
        taus.append(1.0 + 1e-12)

        self.threshold = orig_thresh
        return precisions, pct_pred_pos, taus
747971
def evaluate_hourly_forecasts(self):
        
        score_columns = ["Run_Date", "Forecast_Hour", "Ensemble Name", "Model_Name", "Forecast_Variable",
                         "Neighbor_Radius", "Smoothing_Radius", "Size_Threshold", "ROC", "Reliability"]
        all_scores = pd.DataFrame(columns=score_columns)
        for h, hour in enumerate(range(self.start_hour, self.end_hour + 1)):
            for neighbor_radius in self.neighbor_radii:
                n_filter = disk(neighbor_radius)
                for s, size_threshold in enumerate(self.size_thresholds):
                    print("Eval hourly forecast {0:02d} {1} {2} {3} {4:d} {5:d}".format(hour, self.model_name,
                                                                                        self.forecast_variable,
                                                                                        self.run_date, neighbor_radius,
                                                                                        size_threshold))
                    hour_obs = fftconvolve(self.raw_obs[self.mrms_variable][h] >= self.obs_thresholds[s],
                                           n_filter, mode="same")
                    hour_obs[hour_obs > 1] = 1
                    hour_obs[hour_obs < 1] = 0
                    if self.obs_mask:
                        hour_obs = hour_obs[self.raw_obs[self.mask_variable][h] > 0]
                    for smoothing_radius in self.smoothing_radii:
                        hour_var = "neighbor_prob_r_{0:d}_s_{1:d}_{2}_{3:0.2f}".format(neighbor_radius,
                                                                                       smoothing_radius,
                                                                                       self.forecast_variable,
                                                                                       size_threshold)
                        if self.obs_mask:
                            hour_forecast = self.hourly_forecasts[hour_var][h][self.raw_obs[self.mask_variable][h] > 0]
                        else:
                            hour_forecast = self.hourly_forecasts[hour_var][h]
                        roc = DistributedROC(thresholds=self.probability_levels, obs_threshold=0.5)
                        roc.update(hour_forecast, hour_obs)
                        rel = DistributedReliability(thresholds=self.probability_levels, obs_threshold=0.5)
                        rel.update(hour_forecast, hour_obs)
                        row = [self.run_date, hour, self.ensemble_name, self.model_name, self.forecast_variable,
                               neighbor_radius,
                               smoothing_radius, size_threshold, roc, rel]
                        all_scores.loc[hour_var + "_{0:d}".format(hour)] = row
        return all_scores
184329
def crps(self):
        
        return np.sum(self.errors["F_2"].values - self.errors["F_O"].values * 2.0 + self.errors["O_2"].values) / \
            (self.thresholds.size * self.num_forecasts)
1213762
def get_rates(self, mmin, mmax=np.inf):
        
        nsrcs = self.number_sources()
        for iloc, source in enumerate(self.source_model):
            print("Source Number %s of %s, Name = %s, Typology = %s" % (
                iloc + 1,
                nsrcs,
                source.name,
                source.__class__.__name__))
            if isinstance(source, CharacteristicFaultSource):
                self._get_fault_rates(source, mmin, mmax)
            elif isinstance(source, ComplexFaultSource):
                self._get_fault_rates(source, mmin, mmax)
            elif isinstance(source, SimpleFaultSource):
                self._get_fault_rates(source, mmin, mmax)
            elif isinstance(source, AreaSource):
                self._get_area_rates(source, mmin, mmax)
            elif isinstance(source, PointSource):
                self._get_point_rates(source, mmin, mmax)
            else:
                print("Source type %s not recognised - skipping!" % source)
                continue
325343
def confidence_interval_min_width(mu, post, alpha=0.9):
    '''
    
    '''
    if not 0 < alpha < 1:
        raise ValueError("Confidence level must be in (0,1).")

    # choose a step size for the sliding confidence window
    alpha_step = 0.01

    # initialize the lower and upper limits
    mu_low = numpy.min(mu)
    mu_high = numpy.max(mu)

    # find the smallest window (by delta-mu) stepping by dalpha
    for ai in numpy.arange(0, 1 - alpha, alpha_step):
        ml = compute_lower_limit(mu, post, 1 - ai)
        mh = compute_upper_limit(mu, post, alpha + ai)
        if mh - ml < mu_high - mu_low:
            mu_low = ml
            mu_high = mh

    return mu_low, mu_high