Calculates the limiting distribution of the Cramer-von Mises statistic.

    After the second line of equation 1.3 from the Csorgo and Faraway paper.
514640
def cvm_unif_inf(statistic):
    
    args = inf_args / statistic
    return (inf_cs * exp(-args) * kv(.25, args)).sum() / statistic ** .5
738702
def sde(self):
        
        variance = float(self.variance.values)
        lengthscale = float(self.lengthscale)

        F  = np.array(((-1.0/lengthscale,),))
        L  = np.array(((1.0,),))
        Qc = np.array( ((2.0*variance/lengthscale,),) )
        H = np.array(((1.0,),))
        Pinf = np.array(((variance,),))
        P0 = Pinf.copy()

        dF = np.zeros((1,1,2));
        dQc = np.zeros((1,1,2));
        dPinf = np.zeros((1,1,2));

        dF[:,:,0] = 0.0
        dF[:,:,1] = 1.0/lengthscale**2

        dQc[:,:,0] = 2.0/lengthscale
        dQc[:,:,1] = -2.0*variance/lengthscale**2

        dPinf[:,:,0] = 1.0
        dPinf[:,:,1] = 0.0

        dP0 = dPinf.copy()

        return (F, L, Qc, H, Pinf, P0, dF, dQc, dPinf, dP0)
400421
def rvon_mises(mu, kappa, size=None):
    
    # TODO: Just return straight from numpy after release 1.3
    return (np.random.mtrand.vonmises(
        mu, kappa, size) + np.pi) % (2. * np.pi) - np.pi
784740
def cdf(self, x):
        
        t = 0
        N = float(self.n)

        if len(self) == 1:  # only one centroid
            return int(x >= self.C.min_key())

        for i, key in enumerate(self.C.keys()):
            c_i = self.C[key]
            if i == len(self) - 1:
                delta = (c_i.mean - self.C.prev_item(key)[1].mean) / 2.
            else:
                delta = (self.C.succ_item(key)[1].mean - c_i.mean) / 2.
            z = max(-1, (x - c_i.mean) / delta)

            if z < 1:
                return t / N + c_i.count / N * (z + 1) / 2

            t += c_i.count
        return 1
728607
def sse(mean, estimator):
    
    return np.sum((np.asarray(estimator) - np.asarray(mean)) ** 2, axis=0)
1305481
def rho_rec_limit0(self):
        
        return -2 * self.gamma2 * self.expnorm - 0.25 * self.beta2 * (
            1 - self.expnorm)
357605
def probability_density(self, X):
        """"""
        self.check_fit()

        U, V = self.split_matrix(X)

        if self.theta == 1:
            return np.multiply(U, V)

        else:
            a = np.power(np.multiply(U, V), -1)
            tmp = np.power(-np.log(U), self.theta) + np.power(-np.log(V), self.theta)
            b = np.power(tmp, -2 + 2.0 / self.theta)
            c = np.power(np.multiply(np.log(U), np.log(V)), self.theta - 1)
            d = 1 + (self.theta - 1) * np.power(tmp, -1.0 / self.theta)
            return self.cumulative_distribution(X) * a * b * c * d
613068
private void doMutation(double probability, DoubleSolution solution) {
    double rnd, delta1, delta2, mutPow, deltaq;
    double y, yl, yu, val, xy;

    for (int i = 0; i < solution.getNumberOfVariables(); i++) {
      if (randomGenerator.getRandomValue() <= probability) {
        y = solution.getVariableValue(i);
        yl = solution.getLowerBound(i) ;
        yu = solution.getUpperBound(i) ;
        if (yl == yu) {
          y = yl ;
        } else {
          delta1 = (y - yl) / (yu - yl);
          delta2 = (yu - y) / (yu - yl);
          rnd = randomGenerator.getRandomValue();
          mutPow = 1.0 / (distributionIndex + 1.0);
          if (rnd <= 0.5) {
            xy = 1.0 - delta1;
            val = 2.0 * rnd + (1.0 - 2.0 * rnd) * (Math.pow(xy, distributionIndex + 1.0));
            deltaq = Math.pow(val, mutPow) - 1.0;
          } else {
            xy = 1.0 - delta2;
            val = 2.0 * (1.0 - rnd) + 2.0 * (rnd - 0.5) * (Math.pow(xy, distributionIndex + 1.0));
            deltaq = 1.0 - Math.pow(val, mutPow);
          }
          y = y + deltaq * (yu - yl);
          y = solutionRepair.repairSolutionVariableValue(y, yl, yu);
        }
        solution.setVariableValue(i, y);
      }
    }
  }
1005720
def calc_pvalue(self, study_count, study_n, pop_count, pop_n):
        """"""
        # k, n = study_true, study_tot,
        # K, N = population_true, population_tot
        # def pvalue_population(int k, int n, int K, int N): ...
        return self.pval_fnc(study_count, study_n, pop_count, pop_n).two_tail
67728
def density(self, r, rho0, Rs):
        
        rho = rho0 / (r/Rs * (1 + (r/Rs))**3)
        return rho
602121
def _compute_mle(self):
        
        min_y = np.min(self._interp.y)
        if self._interp.y[0] == min_y:
            self._mle = self._interp.x[0]
        elif self._interp.y[-1] == min_y:
            self._mle = self._interp.x[-1]
        else:
            argmin_y = np.argmin(self._interp.y)
            ix0 = max(argmin_y - 4, 0)
            ix1 = min(argmin_y + 4, len(self._interp.x) - 1)

            while np.sign(self._interp.derivative(self._interp.x[ix0])) == \
                    np.sign(self._interp.derivative(self._interp.x[ix1])):
                ix0 += 1

            self._mle = scipy.optimize.brentq(self._interp.derivative,
                                              self._interp.x[ix0],
                                              self._interp.x[ix1],
                                              xtol=1e-10 *
                                              np.median(self._interp.x))
686924
public static double cdf(double val, double mu, double sigma, double k) {
    final double x = (val - mu) / sigma;
    if(k > 0 || k < 0) {
      if(k * x > 1) {
        return k > 0 ? 1 : 0;
      }
      return FastMath.exp(-FastMath.exp(FastMath.log(1 - k * x) / k));
    }
    else { // Gumbel case:
      return FastMath.exp(-FastMath.exp(-x));
    }
  }
1034431
public static double cdf(double val, double k, double lambda, double theta) {
    return (val > theta) ? //
        (1.0 - FastMath.exp(-FastMath.pow((val - theta) / lambda, k))) //
        : val == val ? 0.0 : Double.NaN;
  }
532363
def erfinv(x, a=.147):
    
    lnx = log(1 - x * x)
    part1 = (2 / (a * pi) + lnx / 2)
    part2 = lnx / a
    sgn = 1 if x > 0 else -1
    return sgn * sqrt(sqrt(part1 * part1 - part2) - part1)
771400
def _minimize_c(c, z=0, a_tilde=1, b_tilde=-1,
                Ascaling=900, omega_M_0=0.25, omega_lambda_0=0.75):
    

    # Fn 1 (LHS of Eqn 18)

    Y1 = np.log(2) - 0.5
    Yc = np.log(1+c) - c/(1+c)
    f1 = Y1/Yc

    # Fn 2 (RHS of Eqn 18)

    # Eqn 14 - Define the mean inner density
    rho_2 = 200 * c**3 * Y1 / Yc

    # Eqn 17 rearranged to solve for Formation Redshift
    # essentially when universe had rho_2 density
    zf = (((1 + z)**3 + omega_lambda_0/omega_M_0) *
          (rho_2/Ascaling) - omega_lambda_0/omega_M_0)**(1/3) - 1

    # RHS of Eqn 19
    f2 = ((1 + zf - z)**a_tilde) * np.exp((zf - z) * b_tilde)

    # LHS - RHS should be zero for the correct concentration
    return(f1-f2)
699394
private static double snorm(double[] x, int itol) {
        int n = x.length;

        if (itol <= 3) {
            double ans = 0.0;
            for (int i = 0; i < n; i++) {
                ans += x[i] * x[i];
            }
            return Math.sqrt(ans);
        } else {
            int isamax = 0;
            for (int i = 0; i < n; i++) {
                if (Math.abs(x[i]) > Math.abs(x[isamax])) {
                    isamax = i;
                }
            }

            return Math.abs(x[isamax]);
        }
    }
883858
def kurtosis(data):
    

    if not data:
        raise StatisticsError('kurtosis requires at least one data point')

    size = len(data)
    sd = stdev(data) ** 4

    if not sd:
        return 0.0

    mn = mean(data)
    return sum(map(lambda x: ((x - mn) ** 4 / sd), data)) / size - 3