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

Example 6 with NormalDistribution

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

the class LoadingsGibbsTruncatedOperator method drawI.

private void drawI(int i, int column) {
    double[] draws = null;
    precisionArray = new double[loadings.getColumnDimension()][loadings.getColumnDimension()];
    double[][] variance;
    meanMidArray = new double[loadings.getColumnDimension()];
    meanArray = new double[loadings.getColumnDimension()];
    double[][] cholesky = null;
    NormalDistribution conditioned;
    getPrecision(i, precisionArray);
    if (LFM.getLoadings().getParameterValue(i, column) != 0) {
        variance = (new SymmetricMatrix(precisionArray)).inverse().toComponents();
        //        try {
        //            cholesky = new CholeskyDecomposition(variance).getL();
        //        } catch (IllegalDimension illegalDimension) {
        //            illegalDimension.printStackTrace();
        //        }
        getMean(i, variance, meanMidArray, meanArray);
        if (LFM.getFactorDimension() != 1)
            conditioned = getConditionalDistribution(meanArray, variance, column, i);
        else
            conditioned = new NormalDistribution(meanArray[0], Math.sqrt(variance[0][0]));
    } else
        //TODO generify
        conditioned = new NormalDistribution(0, Math.sqrt(1 / priorPrecision));
    if (prior instanceof MomentDistributionModel) {
        if (MathUtils.nextDouble() < .5) {
            getTruncatedDraw(i, column, conditioned);
            getCutoffDraw(i, column, conditioned);
        } else {
            getCutoffDraw(i, column, conditioned);
            getTruncatedDraw(i, column, conditioned);
        }
    } else {
        loadings.setParameterValue(i, column, conditioned.quantile(MathUtils.nextDouble()));
    }
}

Example 7 with NormalDistribution

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

the class LatentLiabilityGibbs method sampleNode.

public double sampleNode(NodeRef node) {
    final int thisNumber = node.getNumber();
    double[] traitValue = getNodeTrait(node);
    double[] mean = new double[dim];
    for (int i = 0; i < dim; i++) {
        mean[i] = preMeans[thisNumber][i];
    }
    double p = preP[thisNumber];
    double[][] thisP = new double[dim][dim];
    for (int i = 0; i < dim; i++) {
        for (int j = 0; j < dim; j++) {
            thisP[i][j] = p * precisionParam.getParameterValue(i, j);
        }
    }
    /*	  Sample it all traits together as one multivariate normal
            *
               MultivariateNormalDistribution distribution  = new MultivariateNormalDistribution(mean, thisP);


 		double[] oldValue = getNodeTrait(node);
 		double[] value = distribution.nextMultivariateNormal();

 		setNodeTrait(node,value);



		double pOld = distribution.logPdf(oldValue[]);

		double pNew = distribution.logPdf(value[]);


 	//	printInformation(oldValue[0]);

 	//	printInformation(value[0]);

  */
    // double[] newTraitValue = getNodeTrait(node);
    //double pNew = distribution.logPdf(newTraitValue);
    ///////////  Individually gibbs sample each entry in the vector
    //  for(int entry=0;entry<dim; entry++){
    int entry = MathUtils.nextInt(dim);
    double thisMean = getConditionalMean(entry, thisP, traitValue, mean);
    double SD = Math.sqrt(1 / thisP[entry][entry]);
    double oldValue = getNodeTrait(node, entry);
    double value = MathUtils.nextGaussian();
    value *= SD;
    value += thisMean;
    //	 printInformation(oldValue);
    //    printInformation(value);
    NormalDistribution distribution = new NormalDistribution(thisMean, SD);
    double pOld = distribution.logPdf(oldValue);
    double pNew = distribution.logPdf(value);
    setNodeTrait(node, entry, value);
    double logq = pOld - pNew;
    traitModel.getTraitParameter().getParameter(thisNumber).fireParameterChangedEvent();
    return logq;
}