Examples with NormalDistribution dr.math.distributions.NormalDistribution used on opensource projects

Example 1 with NormalDistribution

use of dr.math.distributions.NormalDistribution in project beast-mcmc by beast-dev.

the class SliceOperator method main.

public static void main(String[] arg) {
    // Define normal model
    // Starting value
    Parameter meanParameter = new Parameter.Default(1.0);
    // Fixed value
    Variable<Double> stdev = new Variable.D(1.0, 1);
    ParametricDistributionModel densityModel = new NormalDistributionModel(meanParameter, stdev);
    DistributionLikelihood likelihood = new DistributionLikelihood(densityModel);
    // Define prior
    // Hyper-priors
    DistributionLikelihood prior = new DistributionLikelihood(new NormalDistribution(0.0, 1.0));
    prior.addData(meanParameter);
    // Define data
    likelihood.addData(new Attribute.Default<double[]>("Data", new double[] { 0.0, 1.0, 2.0 }));
    List<Likelihood> list = new ArrayList<Likelihood>();
    list.add(likelihood);
    list.add(prior);
    CompoundLikelihood posterior = new CompoundLikelihood(0, list);
    SliceOperator sliceSampler = new SliceOperator(meanParameter);
    final int length = 10000;
    double mean = 0;
    double variance = 0;
    for (int i = 0; i < length; i++) {
        sliceSampler.doOperation(posterior);
        double x = meanParameter.getValue(0);
        mean += x;
        variance += x * x;
    }
    mean /= length;
    variance /= length;
    variance -= mean * mean;
    System.out.println("E(x) = " + mean);
    System.out.println("V(x) = " + variance);
}

Example 2 with NormalDistribution

use of dr.math.distributions.NormalDistribution in project beast-mcmc by beast-dev.

the class RiemannApproximation method main.

public static void main(String[] args) {
    UnivariateFunction normalPDF = new NormalDistribution(0.0, 1.0).getProbabilityDensityFunction();
    UnivariateFunction normalPDF2 = new NormalDistribution(0.0, 1.0).getProbabilityDensityFunction();
    UnivariateFunction normalPDF3 = new NormalDistribution(0.0, 1.0).getProbabilityDensityFunction();
    double Z = 1.0;
    //double Z = Math.sqrt(2*Math.PI);
    CompoundFunction threeNormals = new CompoundFunction(new UnivariateFunction[] { normalPDF, normalPDF2, normalPDF3 }, Z);
    System.out.println("Riemann approximation to the integral of a three normal distribution:");
    RiemannApproximation integrator = new RiemannApproximation(100000);
    System.out.println("integral(N(0.0, 1.0))=" + integrator.integrate(normalPDF, -4.0, 4.0));
    System.out.println("integral(N(1.0, 2.0))=" + integrator.integrate(normalPDF2, -8.0, 8.0));
    System.out.println("integral(N(2.0, 3.0))=" + integrator.integrate(normalPDF3, -16.0, 16.0));
    double integral = integrator.integrate(threeNormals, -16.0, 16.0);
    System.out.println("Riemann approximation to the integral of the compound of three normal distribution:");
    System.out.println("integral(N(0.0, 1.0)*N(1.0, 2.0)*N(2.0, 3.0))=" + integral);
    System.out.println("Estimate normalizing constant is " + (1.0 / integral));
/*System.out.println("Ten monte carlo approximations to the integral of a normal distribution:");
          MonteCarloIntegral integrator2 = new MonteCarloIntegral(10000);
          for (int i = 0; i < 10; i++) {
              System.out.println(integrator2.integrate(normalPDF, -4.0, 4.0));
          }*/
}

Example 3 with NormalDistribution

use of dr.math.distributions.NormalDistribution in project beast-mcmc by beast-dev.

the class DirichletProcessPriorLogger method getNew.

// END: getColumns
private void getNew() {
    this.categoryProbabilities = getCategoryProbs();
    this.newCategoryIndex = Utils.sample(categoryProbabilities);
    this.meanForCategory = uniquelyRealizedParameters.getParameterValue(newCategoryIndex);
    //TODO: generalize
    double sd = precisionParameter.getParameterValue(0);
    //		System.out.println("FUBAR:" + sd);
    NormalDistribution nd = new NormalDistribution(meanForCategory, sd);
    this.newX = (Double) nd.nextRandom();
}

Example 4 with NormalDistribution

use of dr.math.distributions.NormalDistribution in project beast-mcmc by beast-dev.

the class LatentFactorLiabilityGibbsOperator method drawTruncatedNormalDistribution.

double drawTruncatedNormalDistribution(double mean, double precision, double lower, double upper) {
    double sd = Math.sqrt(1 / precision);
    NormalDistribution normal = new NormalDistribution(mean, sd);
    double newLower = normal.cdf(lower);
    double newUpper = normal.cdf(upper);
    double cdfDraw = 1.0;
    int iterator = 0;
    boolean invalid = true;
    double draw = 0;
    while (iterator < 10000 && invalid) {
        cdfDraw = MathUtils.nextDouble() * (newUpper - newLower) + newLower;
        draw = normal.quantile(cdfDraw);
        if (!Double.isNaN(draw) && draw > lower && draw < upper) {
            invalid = false;
        }
        iterator++;
    }
    //            System.out.println(upper);}
    if (Double.isNaN(draw) || Double.isInfinite(draw)) {
        if (Double.isInfinite(lower)) {
            //                System.out.println(upper);
            return upper;
        } else if (Double.isInfinite(upper)) {
            //                System.out.println(lower);
            return lower;
        } else
            return (lower + upper) / 2;
    } else
        return draw;
}

Example 5 with NormalDistribution

use of dr.math.distributions.NormalDistribution in project beast-mcmc by beast-dev.

the class LoadingsGibbsTruncatedOperator method getConditionalDistribution.

private NormalDistribution getConditionalDistribution(double[] meanArray, double[][] variance, int column, int row) {
    double[][] newVariance = new double[meanArray.length - 1][meanArray.length - 1];
    for (int i = 0; i < meanArray.length; i++) {
        for (int j = 0; j < meanArray.length; j++) {
            if (i < column && j < column) {
                newVariance[i][j] = variance[i][j];
            } else if (i < column && j > column) {
                newVariance[i][j - 1] = variance[i][j];
            } else if (i > column && j < column) {
                newVariance[i - 1][j] = variance[i][j];
            } else if (i > column && j > column) {
                newVariance[i - 1][j - 1] = variance[i][j];
            } else {
            }
        }
    }
    double[][] smallPrecision = (new SymmetricMatrix(newVariance)).inverse().toComponents();
    double[] meanStore1 = new double[meanArray.length - 1];
    double[] meanStore2 = new double[meanArray.length - 1];
    double[] precStore = new double[meanArray.length - 1];
    for (int i = 0; i < meanArray.length; i++) {
        if (i < column) {
            meanStore1[i] = LFM.getLoadings().getParameterValue(row, i) - meanArray[i];
        } else if (i > column) {
            meanStore1[i - 1] = LFM.getLoadings().getParameterValue(row, i) - meanArray[i];
        } else {
        }
    }
    for (int i = 0; i < meanArray.length - 1; i++) {
        for (int j = 0; j < meanArray.length - 1; j++) {
            meanStore2[i] += smallPrecision[i][j] * meanStore1[j];
        }
    }
    double mean = meanArray[column];
    for (int i = 0; i < meanArray.length - 1; i++) {
        if (i < column) {
            mean += meanStore2[i] * variance[i][column];
        } else {
            mean += meanStore2[i] * variance[i + 1][column];
        }
    }
    for (int i = 0; i < meanArray.length - 1; i++) {
        for (int j = 0; j < meanArray.length - 1; j++) {
            if (i < column)
                precStore[i] += smallPrecision[i][j] * variance[j][column];
            else
                precStore[i] += smallPrecision[i][j] * variance[j + 1][column];
        }
    }
    double varianceElement = variance[column][column];
    for (int i = 0; i < meanArray.length - 1; i++) {
        if (i < column)
            varianceElement -= precStore[i] * variance[i][column];
        else
            varianceElement -= precStore[i] * variance[i + 1][column];
    }
    return new NormalDistribution(mean, Math.sqrt(varianceElement));
}