Examples with Matrix dr.math.matrixAlgebra.Matrix used on opensource projects

Example 46 with Matrix

use of dr.math.matrixAlgebra.Matrix in project beast-mcmc by beast-dev.

the class FullyConjugateMultivariateTraitLikelihood method getReport.

@Override
public String getReport() {
    StringBuilder sb = new StringBuilder();
    // sb.append(this.g)
    // System.err.println("Hello");
    sb.append("Tree:\n");
    sb.append(getId()).append("\t");
    sb.append(treeModel.toString());
    sb.append("\n\n");
    double[][] treeVariance = computeTreeVariance();
    double[][] traitPrecision = getDiffusionModel().getPrecisionmatrix();
    Matrix traitVariance = new Matrix(traitPrecision).inverse();
    double[][] jointVariance = KroneckerOperation.product(treeVariance, traitVariance.toComponents());
    sb.append("Tree variance:\n");
    sb.append(new Matrix(treeVariance));
    sb.append(matrixMin(treeVariance)).append("\t").append(matrixMax(treeVariance)).append("\t").append(matrixSum(treeVariance));
    sb.append("\n\n");
    sb.append("Trait variance:\n");
    sb.append(traitVariance);
    sb.append("\n\n");
    // sb.append("Joint variance:\n");
    // sb.append(new Matrix(jointVariance));
    // sb.append("\n\n");
    sb.append("Tree dim: ").append(treeVariance.length).append("\n");
    sb.append("data dim: ").append(jointVariance.length);
    sb.append("\n\n");
    double[] data = new double[jointVariance.length];
    System.arraycopy(meanCache, 0, data, 0, jointVariance.length);
    if (nodeToClampMap != null) {
        int offset = treeModel.getExternalNodeCount() * getDimTrait();
        for (Map.Entry<NodeRef, RestrictedPartials> clamps : nodeToClampMap.entrySet()) {
            double[] partials = clamps.getValue().getPartials();
            for (double partial : partials) {
                data[offset] = partial;
                ++offset;
            }
        }
    }
    sb.append("Data:\n");
    sb.append(new Vector(data)).append("\n");
    sb.append(data.length).append("\t").append(vectorMin(data)).append("\t").append(vectorMax(data)).append("\t").append(vectorSum(data));
    sb.append(treeModel.getNodeTaxon(treeModel.getExternalNode(0)).getId());
    sb.append("\n\n");
    MultivariateNormalDistribution mvn = new MultivariateNormalDistribution(new double[data.length], new Matrix(jointVariance).inverse().toComponents());
    double logDensity = mvn.logPdf(data);
    sb.append("logLikelihood: ").append(getLogLikelihood()).append(" == ").append(logDensity).append("\n\n");
    final WishartSufficientStatistics sufficientStatistics = getWishartStatistics();
    final double[] outerProducts = sufficientStatistics.getScaleMatrix();
    sb.append("Outer-products (DP):\n");
    sb.append(new Vector(outerProducts));
    sb.append(sufficientStatistics.getDf()).append("\n");
    Matrix treePrecision = new Matrix(treeVariance).inverse();
    final int n = data.length / traitPrecision.length;
    final int p = traitPrecision.length;
    double[][] tmp = new double[n][p];
    for (int i = 0; i < n; ++i) {
        System.arraycopy(data, i * p, tmp[i], 0, p);
    }
    Matrix y = new Matrix(tmp);
    Matrix S = null;
    try {
        // Using Matrix-Normal form
        S = y.transpose().product(treePrecision).product(y);
    } catch (IllegalDimension illegalDimension) {
        illegalDimension.printStackTrace();
    }
    sb.append("Outer-products (from tree variance:\n");
    sb.append(S);
    sb.append("\n\n");
    return sb.toString();
}

Example 47 with Matrix

use of dr.math.matrixAlgebra.Matrix in project beast-mcmc by beast-dev.

the class ApproximateFactorAnalysisPrecisionMatrix method getParameterAsMatrix.

@Override
public double[][] getParameterAsMatrix() {
    computeValues();
    double[][] matrix = new double[dim][dim];
    for (int i = 0; i < dim; ++i) {
        System.arraycopy(values, i * dim, matrix[i], 0, dim);
    }
    if (DEBUG) {
        System.err.println("vec:");
        System.err.println(new Vector(values));
        System.err.println(new Matrix(matrix));
        System.err.println("");
    }
    return matrix;
}

Example 48 with Matrix

use of dr.math.matrixAlgebra.Matrix in project beast-mcmc by beast-dev.

the class PrecisionMatrixGibbsOperator method doOperationDontFireChange.

public void doOperationDontFireChange() {
    if (wishartIsModel) {
        // TODO Deprecate
        setupWishartStatistics(priorModel);
        priorStatistics = setupStatistics(priorModel);
    }
    final double[][] scaleMatrix = getOperationScaleMatrixAndSetObservationCount();
    final double treeDf = numberObservations;
    final double df = priorDf + treeDf * pathWeight;
    double[][] draw = WishartDistribution.nextWishart(df, scaleMatrix);
    if (DEBUG) {
        System.err.println("draw = " + new Matrix(draw));
    }
    for (int i = 0; i < dim; i++) {
        Parameter column = precisionParam.getParameter(i);
        for (int j = 0; j < dim; j++) column.setParameterValueQuietly(j, draw[j][i]);
    }
}

Example 49 with Matrix

use of dr.math.matrixAlgebra.Matrix in project beast-mcmc by beast-dev.

the class VarianceProportionStatistic method getReport.

@Override
public String getReport() {
    Matrix mat = new Matrix(dimTrait, dimTrait);
    for (int i = 0; i < dimTrait; i++) {
        int offset = dimTrait * i;
        for (int j = 0; j < dimTrait; j++) {
            mat.set(i, j, getStatisticValue(offset + j));
        }
    }
    StringBuilder sb = new StringBuilder();
    sb.append("Variance proportion statistic: " + ratio.name());
    sb.append("\n");
    sb.append("stat value = ");
    sb.append(mat);
    sb.append("\n\n");
    return sb.toString();
}

Example 50 with Matrix

use of dr.math.matrixAlgebra.Matrix in project beast-mcmc by beast-dev.

the class VarianceProportionStatistic method getMatrixSqrt.

// TODO: Move method below to a different class
// TODO: implement this for DenseMatrix64F rather than Matrix
private static Matrix getMatrixSqrt(Matrix M, Boolean invert) {
    DoubleMatrix2D S = new DenseDoubleMatrix2D(M.toComponents());
    RobustEigenDecomposition eigenDecomp = new RobustEigenDecomposition(S, 100);
    DoubleMatrix1D eigenValues = eigenDecomp.getRealEigenvalues();
    int dim = eigenValues.size();
    for (int i = 0; i < dim; i++) {
        double value = sqrt(eigenValues.get(i));
        if (invert) {
            value = 1 / value;
        }
        eigenValues.set(i, value);
    }
    DoubleMatrix2D eigenVectors = eigenDecomp.getV();
    for (int i = 0; i < dim; i++) {
        for (int j = 0; j < dim; j++) {
            eigenVectors.set(i, j, eigenVectors.get(i, j) * eigenValues.get(j));
        }
    }
    DoubleMatrix2D storageMatrix = new DenseDoubleMatrix2D(dim, dim);
    eigenVectors.zMult(eigenDecomp.getV(), storageMatrix, 1, 0, false, true);
    return new Matrix(storageMatrix.toArray());
}