org.apache.commons.math3.distribution.NormalDistribution Java Examples

The following examples show how to use org.apache.commons.math3.distribution.NormalDistribution. You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may check out the related API usage on the sidebar.
Example #1
Source File: RegressionSynthesizer.java    From pyramid with Apache License 2.0 6 votes vote down vote up
public  RegDataSet gaussianMixture(){

        NormalDistribution leftGaussian = new NormalDistribution(0.2,0.01);

        NormalDistribution rightGaussian = new NormalDistribution(0.7,0.1);

        RegDataSet dataSet = RegDataSetBuilder.getBuilder()
                .numDataPoints(numDataPoints)
                .numFeatures(1)
                .dense(true)
                .missingValue(false)
                .build();
        for (int i=0;i<numDataPoints;i++){
            double featureValue = Sampling.doubleUniform(0,1);
            double label;
            if (featureValue>0.5){
                label = leftGaussian.sample();
            } else {
                label = rightGaussian.sample();
            }
            dataSet.setFeatureValue(i,0,featureValue);
            dataSet.setLabel(i,label);
        }
        return dataSet;
    }
 
Example #2
Source File: TruncatedNormalTest.java    From pacaya with Apache License 2.0 6 votes vote down vote up
@Test
public void testStandard() {
    assertTrue(TestUtils.checkThrows(() -> new TruncatedNormal(null,  0.0,  1.0,  1E-8, 0.0), IllegalArgumentException.class));
    double mu = 0.0;
    for (double sigma : Arrays.asList(1.0, 0.5, 3.0)) {
        NormalDistribution sn = new NormalDistribution(mu, sigma);
        TruncatedNormal tr = new TruncatedNormal(null, mu, sigma, Double.NEGATIVE_INFINITY,
                Double.POSITIVE_INFINITY);
        assertTrue(tr.isSupportConnected());
        assertFalse(tr.isSupportLowerBoundInclusive());
        assertFalse(tr.isSupportUpperBoundInclusive());
        assertEquals(Double.NEGATIVE_INFINITY, tr.getSupportLowerBound(), tol);
        assertEquals(Double.POSITIVE_INFINITY, tr.getSupportUpperBound(), tol);
        for (double x : Arrays.asList(0.1, 0.5, 2.3, 0.0, -2.8, -10.0, -50.0, 50.0)) {
            assertEquals(sn.density(x), tr.density(x), tol);
            assertEquals(sn.cumulativeProbability(x), tr.cumulativeProbability(x), tol);
            assertEquals(sn.density(x), TruncatedNormal.densityNonTrunc(x, mu, sigma), tol);
            assertEquals(sn.cumulativeProbability(x), TruncatedNormal.cumulativeNonTrunc(x, mu, sigma), tol);
            for (double y : Arrays.asList(x - 1.5, x - 0.5, x, x + 0.5, x + 1.5)) {
                if (y < x) continue;
                assertEquals(sn.probability(x, y), tr.probability(x, y), tol);
                assertEquals(sn.probability(x, y), TruncatedNormal.probabilityNonTrunc(x, y, mu, sigma), tol);
            }
        }
    }
}
 
Example #3
Source File: SchurTransformerTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
@Test
public void testRandomDataNormalDistribution() {
    for (int run = 0; run < 100; run++) {
        Random r = new Random(System.currentTimeMillis());
        NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);

        // matrix size
        int size = r.nextInt(20) + 4;

        double[][] data = new double[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                data[i][j] = dist.sample();
            }
        }

        RealMatrix m = MatrixUtils.createRealMatrix(data);
        RealMatrix s = checkAEqualPTPt(m);
        checkSchurForm(s);
    }
}
 
Example #4
Source File: SchurTransformerTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
@Test
public void testRandomDataNormalDistribution() {
    for (int run = 0; run < 100; run++) {
        Random r = new Random(System.currentTimeMillis());
        NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);

        // matrix size
        int size = r.nextInt(20) + 4;

        double[][] data = new double[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                data[i][j] = dist.sample();
            }
        }

        RealMatrix m = MatrixUtils.createRealMatrix(data);
        RealMatrix s = checkAEqualPTPt(m);
        checkSchurForm(s);
    }
}
 
Example #5
Source File: CurveEstimationTest.java    From finmath-lib with Apache License 2.0 6 votes vote down vote up
/**
 * Regression matrix (currently no test, just for inspection)
 */
@Test
public void testRegressionMatrix() {
	final double[] X = { 0.0 , 1.0 };
	final double[] Y = { 1.0 , 0.8 };

	if(isJBLASPresent) {
		// The following code only works if JBlas is present

		final AbstractRealDistribution kernel=new NormalDistribution();
		final double K=kernel.density(0.0);
		final DoubleMatrix R=new DoubleMatrix(new double[] {K*(Y[0]+Y[1]),Y[1]*(X[1]-X[0])*K} );
		final DoubleMatrix M=new DoubleMatrix(new double[][] {{2*K,K*(X[1]-X[0])},{K*(X[1]-X[0]),K*(X[1]-X[0])*(X[1]-X[0])}} );
		final double detM= M.get(0,0)*M.get(1, 1)-M.get(1,0)*M.get(1,0);
		DoubleMatrix MInv=new DoubleMatrix(new double[][] {{M.get(1, 1),-M.get(1,0)},{-M.get(1,0),M.get(0,0)}} );
		MInv=MInv.mul(1/detM);
		final DoubleMatrix a=MInv.mmul(R);
		System.out.println(R.toString());
		System.out.println(M.toString());
		System.out.println(a.toString());
	}
}
 
Example #6
Source File: FragPileupGen.java    From HMMRATAC with GNU General Public License v3.0 6 votes vote down vote up
/**
 * Returns a new Distribution 
 * @param p a double specifying which distribution to use
 * @param m a double representing the mean for the desired distribution
 * @param l a double representing the standard deviation, or more generally, the lambda of the desired distribution
 * @return A new AbstractRealDistribution
 */
private AbstractRealDistribution getDist(double p,double m, double l){
	if (p == 1){
		//return new LaplaceDistribution(m,l);
		return new ExponentialDistribution(m);
	}
	if (p == 2){
		return new NormalDistribution(m,l);
	}
	if (p == 0.5 || p == 3){
		return new ExponentialDistribution(m);
	}
	if (p == 0){
		//return new ModifiedLaplaceDistribution(m,l);
		return new ExponentialDistribution(m);
	}
	else{
		return null;
	}
	
}
 
Example #7
Source File: CurveEstimationTest.java    From finmath-lib with Apache License 2.0 6 votes vote down vote up
/**
 * Regression matrix (currently no test, just for inspection)
 */
@Test
public void testRegressionMatrix() {
	final double[] X = { 0.0 , 1.0 };
	final double[] Y = { 1.0 , 0.8 };

	if(isJBLASPresent) {
		// The following code only works if JBlas is present

		final AbstractRealDistribution kernel=new NormalDistribution();
		final double K=kernel.density(0.0);
		final DoubleMatrix R=new DoubleMatrix(new double[] {K*(Y[0]+Y[1]),Y[1]*(X[1]-X[0])*K} );
		final DoubleMatrix M=new DoubleMatrix(new double[][] {{2*K,K*(X[1]-X[0])},{K*(X[1]-X[0]),K*(X[1]-X[0])*(X[1]-X[0])}} );
		final double detM= M.get(0,0)*M.get(1, 1)-M.get(1,0)*M.get(1,0);
		DoubleMatrix MInv=new DoubleMatrix(new double[][] {{M.get(1, 1),-M.get(1,0)},{-M.get(1,0),M.get(0,0)}} );
		MInv=MInv.mul(1/detM);
		final DoubleMatrix a=MInv.mmul(R);
		System.out.println(R.toString());
		System.out.println(M.toString());
		System.out.println(a.toString());
	}
}
 
Example #8
Source File: MannWhitneyUTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Umin smallest Mann-Whitney U value
 * @param n1 number of subjects in first sample
 * @param n2 number of subjects in second sample
 * @return two-sided asymptotic p-value
 * @throws ConvergenceException if the p-value can not be computed
 * due to a convergence error
 * @throws MaxCountExceededException if the maximum number of
 * iterations is exceeded
 */
private double calculateAsymptoticPValue(final double Umin,
                                         final int n1,
                                         final int n2)
    throws ConvergenceException, MaxCountExceededException {

    /* long multiplication to avoid overflow (double not used due to efficiency
     * and to avoid precision loss)
     */
    final long n1n2prod = (long) n1 * n2;

    // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
    final double EU = n1n2prod / 2.0;
    final double VarU = n1n2prod * (n1 + n2 + 1) / 12.0;

    final double z = (Umin - EU) / FastMath.sqrt(VarU);

    // No try-catch or advertised exception because args are valid
    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2 * standardNormal.cumulativeProbability(z);
}
 
Example #9
Source File: RandomCirclePointGenerator.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param x Abscissa of the circle center.
 * @param y Ordinate of the circle center.
 * @param radius Radius of the circle.
 * @param xSigma Error on the x-coordinate of the circumference points.
 * @param ySigma Error on the y-coordinate of the circumference points.
 * @param seed RNG seed.
 */
public RandomCirclePointGenerator(double x,
                                  double y,
                                  double radius,
                                  double xSigma,
                                  double ySigma,
                                  long seed) {
    final RandomGenerator rng = new Well44497b(seed);
    this.radius = radius;
    cX = new NormalDistribution(rng, x, xSigma,
                                NormalDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
    cY = new NormalDistribution(rng, y, ySigma,
                                NormalDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
    tP = new UniformRealDistribution(rng, 0, MathUtils.TWO_PI,
                                     UniformRealDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
}
 
Example #10
Source File: EigenDecompositionTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
@Test
@Ignore
public void testNormalDistributionUnsymmetricMatrix() {
    for (int run = 0; run < 100; run++) {
        Random r = new Random(System.currentTimeMillis());
        NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);

        // matrix size
        int size = r.nextInt(20) + 4;

        double[][] data = new double[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                data[i][j] = dist.sample();
            }
        }

        RealMatrix m = MatrixUtils.createRealMatrix(data);
        checkUnsymmetricMatrix(m);
    }
}
 
Example #11
Source File: CNLOHCaller.java    From gatk-protected with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
/**
 * Calculate the likelihood for the copy ratio.  Models the credible interval as a Gaussian.
 *
 * @param rho -- see class docs
 * @param m -- see class docs
 * @param n -- see class docs
 * @param lambda -- see class docs
 * @param credibleMode mode of the credible interval from for the copy ratio.
 * @param credibleLow  low point of the credible interval from for the copy ratio.
 * @param credibleHigh  high point of the credible interval from for the copy ratio.
 * @param additionalVariance Additional variance that should be included in the segment mean pdf
 * @param normalPloidy ploidy for the normal cells.  For humans, this is 2.
 * @return likelihood, with a minimum value of {@link CNLOHCaller#MIN_L}
 */
@VisibleForTesting
static double calculateFcr(final double rho, final int m, final int n, final double lambda, final double credibleMode,
                     final double credibleLow, final double credibleHigh, final double additionalVariance, final double normalPloidy) {
    ParamUtils.inRange(rho, 0.0, 1.0, "Invalid rho value: " + rho + ".  Must be [0,1]");
    ParamUtils.isPositiveOrZero(m, "M must be positive.");
    ParamUtils.isPositiveOrZero(n, "N must be positive.");
    ParamUtils.isPositiveOrZero(additionalVariance, "additional variance must be positive.");
    if (Double.isNaN(credibleMode) || Double.isNaN(credibleLow) || Double.isNaN(credibleHigh)){
        return 1;
    }
    final double cr = calculateCopyRatio(rho, m, n, lambda, normalPloidy);
    final double avgDist = ((credibleMode - credibleLow) + (credibleHigh - credibleMode))/2.0;
    return Math.max(new NormalDistribution(null, credibleMode,
            (avgDist / NUM_STD_95_CONF_INTERVAL_GAUSSIAN) + Math.sqrt(additionalVariance)).density(cr), MIN_L);
}
 
Example #12
Source File: MannWhitneyUTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Umin smallest Mann-Whitney U value
 * @param n1 number of subjects in first sample
 * @param n2 number of subjects in second sample
 * @return two-sided asymptotic p-value
 * @throws ConvergenceException if the p-value can not be computed
 * due to a convergence error
 * @throws MaxCountExceededException if the maximum number of
 * iterations is exceeded
 */
private double calculateAsymptoticPValue(final double Umin,
                                         final int n1,
                                         final int n2)
    throws ConvergenceException, MaxCountExceededException {

    /* long multiplication to avoid overflow (double not used due to efficiency
     * and to avoid precision loss)
     */
    final long n1n2prod = (long) n1 * n2;

    // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
    final double EU = n1n2prod / 2.0;
    final double VarU = n1n2prod * (n1 + n2 + 1) / 12.0;

    final double z = (Umin - EU) / FastMath.sqrt(VarU);

    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2 * standardNormal.cumulativeProbability(z);
}
 
Example #13
Source File: SliceSamplerUnitTest.java    From gatk-protected with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
@Test(expectedExceptions = IllegalArgumentException.class)
public void testInitialPointOutOfRange() {
    rng.setSeed(RANDOM_SEED);

    final double mean = 5.;
    final double standardDeviation = 0.75;
    final NormalDistribution normalDistribution = new NormalDistribution(mean, standardDeviation);
    final Function<Double, Double> normalLogPDF = normalDistribution::logDensity;

    final double xInitial = -10.;
    final double xMin = 0.;
    final double xMax = 1.;
    final double width = 0.5;
    final SliceSampler normalSampler = new SliceSampler(rng, normalLogPDF, xMin, xMax, width);
    normalSampler.sample(xInitial);
}
 
Example #14
Source File: CoverageDropoutDetectorTest.java    From gatk-protected with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
private Object[][] getUnivariateGaussianTargets(final double sigma) {
    Random rng = new Random(337);
    final RandomGenerator randomGenerator = RandomGeneratorFactory.createRandomGenerator(rng);
    NormalDistribution n = new NormalDistribution(randomGenerator, 1, sigma);
    final int numDataPoints = 10000;
    final int numEventPoints = 2000;

    final List<ReadCountRecord.SingleSampleRecord> targetList = new ArrayList<>();
    for (int i = 0; i < (numDataPoints - numEventPoints); i++){
        targetList.add(new ReadCountRecord.SingleSampleRecord(new Target("arbitrary_name", new SimpleInterval("chr1", 100 + 2*i, 101 + 2 * i)), n.sample()));
    }
    for (int i = (numDataPoints - numEventPoints); i < numDataPoints; i++){
        targetList.add(new ReadCountRecord.SingleSampleRecord(new Target("arbitrary_name", new SimpleInterval("chr1", 100 + 2 * i, 101 + 2 * i)), 0.5 + n.sample()));
    }

    HashedListTargetCollection<ReadCountRecord.SingleSampleRecord> targets = new HashedListTargetCollection<>(targetList);

    List<ModeledSegment> segments = new ArrayList<>();
    segments.add(new ModeledSegment(new SimpleInterval("chr1", 100, 16050), 8000, 1));
    segments.add(new ModeledSegment(new SimpleInterval("chr1", 16100, 20200), 2000, 1.5));

    return new Object [] []{ {targets, segments}};
}
 
Example #15
Source File: MannWhitneyUTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Umin smallest Mann-Whitney U value
 * @param n1 number of subjects in first sample
 * @param n2 number of subjects in second sample
 * @return two-sided asymptotic p-value
 * @throws ConvergenceException if the p-value can not be computed
 * due to a convergence error
 * @throws MaxCountExceededException if the maximum number of
 * iterations is exceeded
 */
private double calculateAsymptoticPValue(final double Umin,
                                         final int n1,
                                         final int n2)
    throws ConvergenceException, MaxCountExceededException {

    final int n1n2prod = n1 * n2;

    // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
    final double EU = (double) n1n2prod / 2.0;
    final double VarU = (double) (n1n2prod * (n1 + n2 + 1)) / 12.0;

    final double z = (Umin - EU) / FastMath.sqrt(VarU);

    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2 * standardNormal.cumulativeProbability(z);
}
 
Example #16
Source File: EigenDecompositionTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
@Test
public void testNormalDistributionUnsymmetricMatrix() {
    for (int run = 0; run < 100; run++) {
        Random r = new Random(System.currentTimeMillis());
        NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);

        // matrix size
        int size = r.nextInt(20) + 4;

        double[][] data = new double[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                data[i][j] = dist.sample();
            }
        }

        RealMatrix m = MatrixUtils.createRealMatrix(data);
        checkUnsymmetricMatrix(m);
    }
}
 
Example #17
Source File: WilcoxonSignedRankTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Wmin smallest Wilcoxon signed rank value
 * @param N number of subjects (corresponding to x.length)
 * @return two-sided asymptotic p-value
 */
private double calculateAsymptoticPValue(final double Wmin, final int N) {

    final double ES = (double) (N * (N + 1)) / 4.0;

    /* Same as (but saves computations):
     * final double VarW = ((double) (N * (N + 1) * (2*N + 1))) / 24;
     */
    final double VarS = ES * ((double) (2 * N + 1) / 6.0);

    // - 0.5 is a continuity correction
    final double z = (Wmin - ES - 0.5) / FastMath.sqrt(VarS);

    // No try-catch or advertised exception because args are valid
    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2*standardNormal.cumulativeProbability(z);
}
 
Example #18
Source File: buggyMannWhitneyUTest.java    From coming with MIT License 6 votes vote down vote up
/**
 * @param Umin smallest Mann-Whitney U value
 * @param n1 number of subjects in first sample
 * @param n2 number of subjects in second sample
 * @return two-sided asymptotic p-value
 * @throws ConvergenceException if the p-value can not be computed
 * due to a convergence error
 * @throws MaxCountExceededException if the maximum number of
 * iterations is exceeded
 */
private double calculateAsymptoticPValue(final double Umin,
                                         final int n1,
                                         final int n2)
    throws ConvergenceException, MaxCountExceededException {

    final int n1n2prod = n1 * n2;

    // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
    final double EU = n1n2prod / 2.0;
    final double VarU = n1n2prod * (n1 + n2 + 1) / 12.0;

    final double z = (Umin - EU) / FastMath.sqrt(VarU);

    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2 * standardNormal.cumulativeProbability(z);
}
 
Example #19
Source File: MannWhitneyUTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Umin smallest Mann-Whitney U value
 * @param n1 number of subjects in first sample
 * @param n2 number of subjects in second sample
 * @return two-sided asymptotic p-value
 * @throws ConvergenceException if the p-value can not be computed
 * due to a convergence error
 * @throws MaxCountExceededException if the maximum number of
 * iterations is exceeded
 */
private double calculateAsymptoticPValue(final double Umin,
                                         final int n1,
                                         final int n2)
    throws ConvergenceException, MaxCountExceededException {

    /* long multiplication to avoid overflow (double not used due to efficiency
     * and to avoid precision loss)
     */
    final long n1n2prod = (long) n1 * n2;

    // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
    final double EU = n1n2prod / 2.0;
    final double VarU = n1n2prod * (n1 + n2 + 1) / 12.0;

    final double z = (Umin - EU) / FastMath.sqrt(VarU);

    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2 * standardNormal.cumulativeProbability(z);
}
 
Example #20
Source File: NormalDistTest.java    From macrobase with Apache License 2.0 6 votes vote down vote up
@Test
public void checkLUT() {
    int numElements =
            (int)((NormalDist.MAXZSCORE - NormalDist.MINZSCORE)/NormalDist.GRANULARITY) + 1;
    double[] LUT = new double[numElements];
    NormalDistribution dist = new NormalDistribution();
    int minKey = (int) Math.round(NormalDist.MINZSCORE / NormalDist.GRANULARITY);
    int maxKey = (int) Math.round(NormalDist.MAXZSCORE / NormalDist.GRANULARITY);
    int offset = -minKey;
    for (int i = minKey; i <= maxKey; i ++) {
        double zscore = i * NormalDist.GRANULARITY;
        LUT[i+offset] = dist.cumulativeProbability(zscore);
    }

    assertEquals(NormalDist.CDF_LUT.length, LUT.length);
    for (int i = 0; i < LUT.length; i++) {
        assertEquals(NormalDist.CDF_LUT[i], LUT[i], 0.001);
    }
}
 
Example #21
Source File: WilcoxonSignedRankTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Wmin smallest Wilcoxon signed rank value
 * @param N number of subjects (corresponding to x.length)
 * @return two-sided asymptotic p-value
 */
private double calculateAsymptoticPValue(final double Wmin, final int N) {

    final double ES = (double) (N * (N + 1)) / 4.0;

    /* Same as (but saves computations):
     * final double VarW = ((double) (N * (N + 1) * (2*N + 1))) / 24;
     */
    final double VarS = ES * ((double) (2 * N + 1) / 6.0);

    // - 0.5 is a continuity correction
    final double z = (Wmin - ES - 0.5) / FastMath.sqrt(VarS);

    final NormalDistribution standardNormal = new NormalDistribution(0, 1);

    return 2*standardNormal.cumulativeProbability(z);
}
 
Example #22
Source File: QQPlotExample.java    From tablesaw with Apache License 2.0 6 votes vote down vote up
public static void main(String[] args) throws Exception {
  Table baseball = Table.read().csv("../data/baseball.csv");
  Plot.show(QQPlot.create("batting averages and On-base percent", baseball, "BA", "SLG"));

  // example with different sized arrays;
  double[] first = new NormalDistribution().sample(100);
  double[] second = new NormalDistribution().sample(200);
  Plot.show(
      QQPlot.create(
          "Test of different sized arrays", "short array", first, "long array", second));

  // example with different sized arrays;
  first = new NormalDistribution().sample(20000);
  second = new NormalDistribution().sample(19990);
  Plot.show(
      QQPlot.create(
          "Test of different sized arrays", "long array", first, "short array", second));
}
 
Example #23
Source File: WilcoxonSignedRankTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Wmin smallest Wilcoxon signed rank value
 * @param N number of subjects (corresponding to x.length)
 * @return two-sided asymptotic p-value
 */
private double calculateAsymptoticPValue(final double Wmin, final int N) {

    final double ES = (double) (N * (N + 1)) / 4.0;

    /* Same as (but saves computations):
     * final double VarW = ((double) (N * (N + 1) * (2*N + 1))) / 24;
     */
    final double VarS = ES * ((double) (2 * N + 1) / 6.0);

    // - 0.5 is a continuity correction
    final double z = (Wmin - ES - 0.5) / FastMath.sqrt(VarS);

    // No try-catch or advertised exception because args are valid
    // pass a null rng to avoid unneeded overhead as we will not sample from this distribution
    final NormalDistribution standardNormal = new NormalDistribution(null, 0, 1);

    return 2*standardNormal.cumulativeProbability(z);
}
 
Example #24
Source File: WilsonScoreInterval.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/** {@inheritDoc} */
public ConfidenceInterval createInterval(int numberOfTrials, int numberOfSuccesses, double confidenceLevel) {
    IntervalUtils.checkParameters(numberOfTrials, numberOfSuccesses, confidenceLevel);
    final double alpha = (1.0 - confidenceLevel) / 2;
    final NormalDistribution normalDistribution = new NormalDistribution();
    final double z = normalDistribution.inverseCumulativeProbability(1 - alpha);
    final double zSquared = FastMath.pow(z, 2);
    final double mean = (double) numberOfSuccesses / (double) numberOfTrials;

    final double factor = 1.0 / (1 + (1.0 / numberOfTrials) * zSquared);
    final double modifiedSuccessRatio = mean + (1.0 / (2 * numberOfTrials)) * zSquared;
    final double difference = z *
                              FastMath.sqrt(1.0 / numberOfTrials * mean * (1 - mean) +
                                            (1.0 / (4 * FastMath.pow(numberOfTrials, 2)) * zSquared));

    final double lowerBound = factor * (modifiedSuccessRatio - difference);
    final double upperBound = factor * (modifiedSuccessRatio + difference);
    return new ConfidenceInterval(lowerBound, upperBound, confidenceLevel);
}
 
Example #25
Source File: ClusterAlgorithmComparison.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
public static List<Vector2D> makeCircles(int samples, boolean shuffle, double noise, double factor, final RandomGenerator random) {
    if (factor < 0 || factor > 1) {
        throw new IllegalArgumentException();
    }
    
    NormalDistribution dist = new NormalDistribution(random, 0.0, noise, 1e-9);

    List<Vector2D> points = new ArrayList<Vector2D>();
    double range = 2.0 * FastMath.PI;
    double step = range / (samples / 2.0 + 1);
    for (double angle = 0; angle < range; angle += step) {
        Vector2D outerCircle = new Vector2D(FastMath.cos(angle), FastMath.sin(angle));
        Vector2D innerCircle = outerCircle.scalarMultiply(factor);
        
        points.add(outerCircle.add(generateNoiseVector(dist)));
        points.add(innerCircle.add(generateNoiseVector(dist)));
    }
    
    if (shuffle) {
        Collections.shuffle(points, new RandomAdaptor(random));
    }

    return points;
}
 
Example #26
Source File: WilcoxonSignedRankTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param Wmin smallest Wilcoxon signed rank value
 * @param N number of subjects (corresponding to x.length)
 * @return two-sided asymptotic p-value
 */
private double calculateAsymptoticPValue(final double Wmin, final int N) {

    final double ES = (double) (N * (N + 1)) / 4.0;

    /* Same as (but saves computations):
     * final double VarW = ((double) (N * (N + 1) * (2*N + 1))) / 24;
     */
    final double VarS = ES * ((double) (2 * N + 1) / 6.0);

    // - 0.5 is a continuity correction
    final double z = (Wmin - ES - 0.5) / FastMath.sqrt(VarS);

    // No try-catch or advertised exception because args are valid
    // pass a null rng to avoid unneeded overhead as we will not sample from this distribution
    final NormalDistribution standardNormal = new NormalDistribution(null, 0, 1);

    return 2*standardNormal.cumulativeProbability(z);
}
 
Example #27
Source File: ClusterAlgorithmComparison.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
public static List<Vector2D> makeMoons(int samples, boolean shuffle, double noise, RandomGenerator random) {
    NormalDistribution dist = new NormalDistribution(random, 0.0, noise, 1e-9);

    int nSamplesOut = samples / 2;
    int nSamplesIn = samples - nSamplesOut;
    
    List<Vector2D> points = new ArrayList<Vector2D>();
    double range = FastMath.PI;
    double step = range / (nSamplesOut / 2.0);
    for (double angle = 0; angle < range; angle += step) {
        Vector2D outerCircle = new Vector2D(FastMath.cos(angle), FastMath.sin(angle));
        points.add(outerCircle.add(generateNoiseVector(dist)));
    }

    step = range / (nSamplesIn / 2.0);
    for (double angle = 0; angle < range; angle += step) {
        Vector2D innerCircle = new Vector2D(1 - FastMath.cos(angle), 1 - FastMath.sin(angle) - 0.5);
        points.add(innerCircle.add(generateNoiseVector(dist)));
    }
    
    if (shuffle) {
        Collections.shuffle(points, new RandomAdaptor(random));
    }

    return points;
}
 
Example #28
Source File: HessenbergTransformerTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
@Test
public void testRandomDataNormalDistribution() {
    for (int run = 0; run < 100; run++) {
        Random r = new Random(System.currentTimeMillis());
        NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);

        // matrix size
        int size = r.nextInt(20) + 4;

        double[][] data = new double[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                data[i][j] = dist.sample();
            }
        }

        RealMatrix m = MatrixUtils.createRealMatrix(data);
        RealMatrix h = checkAEqualPHPt(m);
        checkHessenbergForm(h);
    }
}
 
Example #29
Source File: RandomCirclePointGenerator.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
/**
 * @param x Abscissa of the circle center.
 * @param y Ordinate of the circle center.
 * @param radius Radius of the circle.
 * @param xSigma Error on the x-coordinate of the circumference points.
 * @param ySigma Error on the y-coordinate of the circumference points.
 * @param seed RNG seed.
 */
public RandomCirclePointGenerator(double x,
                                  double y,
                                  double radius,
                                  double xSigma,
                                  double ySigma,
                                  long seed) {
    final RandomGenerator rng = new Well44497b(seed);
    this.radius = radius;
    cX = new NormalDistribution(rng, x, xSigma,
                                NormalDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
    cY = new NormalDistribution(rng, y, ySigma,
                                NormalDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
    tP = new UniformRealDistribution(rng, 0, MathUtils.TWO_PI,
                                     UniformRealDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
}
 
Example #30
Source File: EigenDecompositionTest.java    From astor with GNU General Public License v2.0 6 votes vote down vote up
@Test
@Ignore
public void testNormalDistributionUnsymmetricMatrix() {
    for (int run = 0; run < 100; run++) {
        Random r = new Random(System.currentTimeMillis());
        NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);

        // matrix size
        int size = r.nextInt(20) + 4;

        double[][] data = new double[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                data[i][j] = dist.sample();
            }
        }

        RealMatrix m = MatrixUtils.createRealMatrix(data);
        checkUnsymmetricMatrix(m);
    }
}