Test multi-column convergence with varying amount of proximal/distal sampling

  :return:
10997
def experimentVaryingSynapseSampling(expParams,
                                     sampleSizeDistalList,
                                     sampleSizeProximalList):
  
  numRpts = 20
  df = None
  args = []
  for sampleSizeProximal in sampleSizeProximalList:
    for sampleSizeDistal in sampleSizeDistalList:

      for rpt in range(numRpts):
        l4Params = getL4Params()
        l2Params = getL2Params()
        l2Params["sampleSizeProximal"] = sampleSizeProximal
        l2Params["minThresholdProximal"] = sampleSizeProximal
        l2Params["sampleSizeDistal"] = sampleSizeDistal
        l2Params["activationThresholdDistal"] = sampleSizeDistal

        args.append(
          {
            "numObjects": expParams['numObjects'],
            "numLocations": expParams['numLocations'],
            "numFeatures": expParams['numFeatures'],
            "numColumns": expParams['numColumns'],
            "trialNum": rpt,
            "l4Params": l4Params,
            "l2Params": l2Params,
            "profile": True,
            "objectSeed": rpt,
          }
        )

  pool = Pool(processes=expParams['numWorkers'])
  result = pool.map(runExperiment, args)

  #
  #       if df is None:
  #         df = pd.DataFrame.from_dict(result, orient='index')
  #       else:
  #         df = pd.concat([df, pd.DataFrame.from_dict(result, orient='index')], axis=1)
  #
  # df = df.transpose()
  return result
919926
def check_constraint_convergence(X, L, LX, Z, U, R, S, step_f, step_g, e_rel, e_abs):
    

    if isinstance(L, list):
        M = len(L)
        convergence = True
        errors = []
        # recursive call
        for i in range(M):
            c, e = check_constraint_convergence(X, L[i], LX[i], Z[i], U[i], R[i], S[i],
                                                step_f, step_g[i], e_rel, e_abs)
            convergence &= c
            errors.append(e)
        return convergence, errors
    else:
        # check convergence of prime residual R and dual residual S
        e_pri, e_dual = get_variable_errors(X, L, LX, Z, U, step_g, e_rel, e_abs)
        lR = l2(R)
        lS = l2(S)
        convergence = (lR <= e_pri) and (lS <= e_dual)
        return convergence, (e_pri, e_dual, lR, lS)
948373
def _locateConvergencePoint(stats, minOverlap, maxOverlap):
    
    for i, v in enumerate(stats[::-1]):
      if not (v >= minOverlap and v <= maxOverlap):
        return len(stats) - i + 1

    # Never differs - converged in one iteration
    return 1
257893
def veq_samples(R_dist,Prot_dist,N=1e4,alpha=0.23,l0=20,sigl=20):
    
    ls = stats.norm(l0,sigl).rvs(N)
    Prots = Prot_dist.rvs(N)
    Prots *= diff_Prot_factor(ls,alpha)
    return R_dist.rvs(N)*2*np.pi*RSUN/(Prots*DAY)/1e5
1013509
def _multitaper_cross_spectrum(self, clm, slm, k, convention='power',
                                   unit='per_l', clat=None, clon=None,
                                   coord_degrees=True, lmax=None,
                                   taper_wt=None):
        
        if lmax is None:
            lmax = min(clm.lmax, slm.lmax)

        if (clat is not None and clon is not None and clat == self.clat and
                clon == self.clon and coord_degrees is self.coord_degrees and
                k <= self.nwinrot):
            # use the already stored coeffs
            pass
        elif (clat is None and clon is None) and \
                (self.clat is not None and self.clon is not None and
                 k <= self.nwinrot):
            # use the already stored coeffs
            pass
        else:
            if clat is None:
                clat = self.clat
            if clon is None:
                clon = self.clon
            if (clat is None and clon is not None) or \
                    (clat is not None and clon is None):
                raise ValueError('clat and clon must both be input. ' +
                                 'clat = {:s}, clon = {:s}'
                                 .format(repr(clat), repr(clon)))
            if clat is None and clon is None:
                self.rotate(clat=90., clon=0., coord_degrees=True, nwinrot=k)
            else:
                self.rotate(clat=clat, clon=clon, coord_degrees=coord_degrees,
                            nwinrot=k)

        sh1 = clm.to_array(normalization='4pi', csphase=1, lmax=lmax)
        sh2 = slm.to_array(normalization='4pi', csphase=1, lmax=lmax)

        if taper_wt is None:
            mtse, sd = _shtools.SHMultiTaperMaskCSE(sh1, sh2, self.coeffs,
                                                    lmax1=lmax, lmax2=lmax,
                                                    k=k)
        else:
            mtse, sd = _shtools.SHMultiTaperMaskCSE(sh1, sh2, self.coeffs,
                                                    lmax1=lmax, lmax2=lmax,
                                                    k=k, taper_wt=taper_wt)

        if (unit == 'per_l'):
            pass
        elif (unit == 'per_lm'):
            degree_l = _np.arange(len(mtse))
            mtse /= (2.0 * degree_l + 1.0)
            sd /= (2.0 * degree_l + 1.0)
        else:
            raise ValueError(
                "unit must be 'per_l' or 'per_lm'." +
                "Input value was {:s}".format(repr(unit)))

        if (convention == 'power'):
            return mtse, sd
        elif (convention == 'energy'):
            return mtse * 4.0 * _np.pi, sd * 4.0 * _np.pi
        else:
            raise ValueError(
                "convention must be 'power' or 'energy'." +
                "Input value was {:s}".format(repr(convention)))
707795
def _fit_m(D, a0, logp, tol=1e-7, maxiter=1000):
    ''''''
    N,K = D.shape
    s = a0.sum()

    for i in xrange(maxiter):
        m = a0 / s
        a1 = _ipsi(logp + (m*(psi(a0) - logp)).sum())
        a1 = a1/a1.sum() * s

        if norm(a1 - a0) < tol:
            return a1
        a0 = a1

    raise Exception('Failed to converge after {} iterations, s is {}'
            .format(maxiter, s))
833806
def _single_multichannel_mp_run(X, Phi, bound, selection_rule, stop_crit,
                                max_iter, verbose=False, pad=0,
                                random_state=None, memory=Memory(None)):
    
    rng = check_random_state(random_state)

    # padding as v stak
    pad = int(pad)
    n_channels = X.shape[0]
    X = np.hstack((np.zeros((n_channels, pad)), X,
                   np.zeros((n_channels, pad))))
    n_samples = X.shape[1]
    n_projs = Phi.doth(X).shape[1]
    err_mse = {}

    # Initialisation
    residual = np.hstack((X.copy(), np.zeros((n_channels,
                                              max(Phi.sizes) / 2))))

    s_rep = np.zeros((n_channels, n_projs))
    X_est = np.zeros((n_channels, n_samples))
    # Main algorithm
    coeffs = np.zeros((n_channels, n_projs))

    it_number = 0
    current_lambda = 1
    for c_idx in range(n_channels):
        err_mse[c_idx] = []
        err_mse[c_idx].append(linalg.norm(residual[c_idx, :]))

    # Decomposition loop: stopping criteria is either SNR or iteration number
    while (current_lambda > bound) & (it_number < max_iter):

        # pick a shift at random : in each size
        rndshifts = {}
        for c_idx in range(n_channels):
            rndshifts[c_idx] = []
        for s_idx, L in enumerate(Phi.sizes):
            shift = rng.randint(low=0, high=L / 4)
            for c_idx in range(n_channels):
                coeffs[c_idx, s_idx * n_samples:(s_idx + 1) * n_samples] = \
                    mdct(residual[c_idx, shift:shift + n_samples], L).ravel()
                rndshifts[c_idx].append(shift)

        # Multichannel mode : we combine projections
        combined = selection_rule(coeffs ** 2)

        # Select a new element
        idx = np.argmax(np.abs(combined))
        # find scale and frequency bin of selected atom
        s_idx = idx // n_samples
        L = Phi.sizes[s_idx]
        F = n_samples // (L // 2)
        frame = (idx - (s_idx * n_samples)) % F
        freq_bin = ((idx - (s_idx * n_samples))) // F

        mdct_wf = memory.cache(mdct_waveform)

        # Update coefficients and residual
        current_lambda_array = np.zeros(n_channels)
        for c_idx in range(n_channels):
            s_rep[c_idx, idx] += coeffs[c_idx, idx]

            # Only one method now : local update via a cached waveform
            pos = (frame * L / 2) - L / 4 + rndshifts[c_idx][s_idx]
            residual[c_idx, pos:pos + L] -= coeffs[c_idx, idx] * \
                mdct_wf(L, freq_bin)

            # also add it to the reconstruction
            X_est[c_idx, pos:pos + L] += coeffs[c_idx, idx] * \
                mdct_wf(L, freq_bin)

            # error computation (err_mse)
            err_mse[c_idx].append(linalg.norm(residual[c_idx, :]))

            current_lambda_array[c_idx] = np.sqrt(
                1. - err_mse[c_idx][-1] / err_mse[c_idx][-2])

        current_lambda = stop_crit(current_lambda_array)

        if verbose:
            print("Iteration %d : Current lambda of %1.4f" % (
                  it_number, current_lambda))
        it_number += 1

    return X_est[:, pad: -pad], err_mse
1027436
def _sample_rows(self, X, Y, sample_shape, seed):
        
        if sample_shape[0] is None or X.shape[0] <= sample_shape[0]:
            X_sample, Y_sample = X, Y
        elif Y is None:
            np.random.seed(seed)
            row_indices = np.random.choice(
                X.shape[0], size=sample_shape[0], replace=False
            )
            X_sample, Y_sample = X.iloc[row_indices], Y
        else:
            drop_size = X.shape[0] - sample_shape[0]
            sample_size = sample_shape[0]
            sss = StratifiedShuffleSplit(
                n_splits=2, test_size=drop_size, train_size=sample_size, random_state=seed
            )
            row_indices, _ = next(sss.split(X, Y))
            X_sample, Y_sample = X.iloc[row_indices], Y.iloc[row_indices]
        return (X_sample, Y_sample)
1166025
def compute_samples(self):
        
        eps = tf.random_normal([self.batch_size, self.eq_samples, self.iw_samples, self.num_latent])
        z = tf.reshape(eps * self.std + self.mu, [-1, self.num_latent])
        return z
470664
def sample_given_context(self, c, c_dims):
        '''
        
        '''
        index = self.discrete_progress.sample_given_context(c, c_dims, self.space)
        return self.space.rand_value(index).flatten()[list(set(range(len(self.space.cardinalities))) - set(c_dims))]
258536
def prox_gradf_lim(xy, step, boundary=None):
    """"""
    return prox_lim(prox_gradf(xy,step), step, boundary=boundary)
406028
def matches_sample(output, target, threshold, is_correct, actual_output):
        
        return (bool(output > threshold) == bool(target)) == is_correct and actual_output == bool(output > threshold)
27580
def plotConvergenceByObject(results, objectRange, featureRange, numTrials,
                            linestyle='-'):
  
  ########################################################################
  #
  # Accumulate all the results per column in a convergence array.
  #
  # Convergence[f,o] = how long it took it to converge with f unique features
  # and o objects.

  convergence = numpy.zeros((max(featureRange), max(objectRange) + 1))
  for r in results:
    if r["numFeatures"] in featureRange:
      convergence[r["numFeatures"] - 1, r["numObjects"]] += r["convergencePoint"]

  convergence /= numTrials

  ########################################################################
  #
  # Create the plot. x-axis=

  # Plot each curve
  legendList = []
  colorList = ['r', 'b', 'g', 'm', 'c', 'k', 'y']

  for i in range(len(featureRange)):
    f = featureRange[i]
    print "features={} objectRange={} convergence={}".format(
      f,objectRange, convergence[f-1,objectRange])
    legendList.append('Unique features={}'.format(f))
    plt.plot(objectRange, convergence[f-1, objectRange],
             color=colorList[i], linestyle=linestyle)

  # format
  plt.legend(legendList, loc="lower right", prop={'size':10})
  plt.xlabel("Number of objects in training set")
  plt.xticks(range(0,max(objectRange)+1,10))
  plt.yticks(range(0,int(convergence.max())+2))
  plt.ylabel("Average number of touches")
  plt.title("Number of touches to recognize one object (single column)")
1059943
def samples_multidimensional_uniform(bounds,num_data):
    '''
    

    '''
    dim = len(bounds)
    Z_rand = np.zeros(shape=(num_data,dim))
    for k in range(0,dim): Z_rand[:,k] = np.random.uniform(low=bounds[k][0],high=bounds[k][1],size=num_data)
    return Z_rand
288156
public static boolean isConvergenced(List<KmeansPoint> basePoints, List<KmeansPoint> newPoints,
            double convergenceThres)
    {
        boolean result = true;

        for (int pointIndex = 0; pointIndex < basePoints.size(); pointIndex++)
        {
            double distance = MathUtils.distance(basePoints.get(pointIndex).getDataPoint(),
                    newPoints.get(pointIndex).getDataPoint());
            if (distance > convergenceThres)
            {
                result = false;
                break;
            }
        }

        return result;
    }
118804
def _ep_pairwise(
        n_items, comparisons, alpha, match_moments, max_iter, initial_state):
    
    # Static variable that allows to check the # of iterations after the call.
    _ep_pairwise.iterations = 0
    m = len(comparisons)
    prior_inv = alpha * np.eye(n_items)
    if initial_state is None:
        # Initially, mean and covariance come from the prior.
        mean = np.zeros(n_items)
        cov = (1 / alpha) * np.eye(n_items)
        # Initialize the natural params in the function space.
        tau = np.zeros(m)
        nu = np.zeros(m)
        # Initialize the natural params in the space of thetas.
        prec = np.zeros((n_items, n_items))
        xs = np.zeros(n_items)
    else:
        tau, nu = initial_state
        mean, cov, xs, prec = _init_ws(
                n_items, comparisons, prior_inv, tau, nu)
    for _ in range(max_iter):
        _ep_pairwise.iterations += 1
        # Keep a copy of the old parameters for convergence testing.
        tau_old = np.array(tau, copy=True)
        nu_old = np.array(nu, copy=True)
        for i in nprand.permutation(m):
            a, b = comparisons[i]
            # Update mean and variance in function space.
            f_var = cov[a,a] + cov[b,b] - 2 * cov[a,b]
            f_mean = mean[a] - mean[b]
            # Cavity distribution.
            tau_tot = 1.0 / f_var
            nu_tot = tau_tot * f_mean
            tau_cav = tau_tot - tau[i]
            nu_cav = nu_tot - nu[i]
            cov_cav = 1.0 / tau_cav
            mean_cav = cov_cav * nu_cav
            # Moment matching.
            logpart, dlogpart, d2logpart = match_moments(mean_cav, cov_cav)
            # Update factor params in the function space.
            tau[i] = -d2logpart / (1 + d2logpart / tau_cav)
            delta_tau = tau[i] - tau_old[i]
            nu[i] = ((dlogpart - (nu_cav / tau_cav) * d2logpart)
                     / (1 + d2logpart / tau_cav))
            delta_nu = nu[i] - nu_old[i]
            # Update factor params in the weight space.
            prec[(a, a, b, b), (a, b, a, b)] += delta_tau * MAT_ONE_FLAT
            xs[a] += delta_nu
            xs[b] -= delta_nu
            # Update mean and covariance.
            if abs(delta_tau) > 0:
                phi = -1.0 / ((1.0 / delta_tau) + f_var) * MAT_ONE
                upd_mat = cov.take([a, b], axis=0)
                cov = cov + upd_mat.T.dot(phi).dot(upd_mat)
            mean = cov.dot(xs)
        # Recompute the global parameters for stability.
        cov = inv_posdef(prior_inv + prec)
        mean = cov.dot(xs)
        if _converged((tau, nu), (tau_old, nu_old)):
            return mean, cov
    raise RuntimeError(
            "EP did not converge after {} iterations".format(max_iter))
1324649
def estimate_column_scales(
            self,
            X_centered,
            row_scales):
        
        n_rows, n_cols = X_centered.shape
        row_scales = np.asarray(row_scales)

        if len(row_scales) != n_rows:
            raise ValueError("Expected length %s, got shape %s" % (
                n_rows, row_scales.shape,))

        column_variances = np.nanmean(
            X_centered ** 2 / (row_scales ** 2).reshape((n_rows, 1)),
            axis=0)
        column_variances[column_variances == 0] = 1.0
        assert len(column_variances) == n_cols, "%d != %d" % (
            len(column_variances),
            n_cols)
        return np.sqrt(column_variances)