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

Example 21 with Vector

use of dr.math.matrixAlgebra.Vector 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(true);
    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: " + treeVariance.length + "\n");
    sb.append("data dim: " + 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 (int i = 0; i < partials.length; ++i) {
                data[offset] = partials[i];
                ++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: " + getLogLikelihood() + " == " + logDensity + "\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() + "\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) {
        for (int j = 0; j < p; ++j) {
            tmp[i][j] = data[i * p + j];
        }
    }
    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 22 with Vector

use of dr.math.matrixAlgebra.Vector 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 23 with Vector

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

the class MarkovModulatedSubstitutionModel method getEigenDecomposition.

//    protected double setupMatrix() {
////        System.err.println("In MM.setupMatrix");
////        setupRelativeRates(relativeRates);
////        double[] pi = freqModel.getFrequencies();
//        setupQMatrix(null, null, q);
////        makeValid(q, stateCount);
//        return 1.0;
//    }
//    public FrequencyModel getFrequencyModel() {
//        return pcFreqModel;
//    }
// TODO Remove
public EigenDecomposition getEigenDecomposition() {
    if (DEBUG) {
        System.err.println("MMSM.getED");
    }
    EigenDecomposition ed = super.getEigenDecomposition();
    if (DEBUG) {
        double[][] q = getQCopy();
        System.err.println(new Matrix(q));
        System.err.println("");
        System.err.println(new dr.math.matrixAlgebra.Vector(ed.getEigenValues()));
        System.err.println("");
        double[] tp = new double[q.length * q.length];
        getTransitionProbabilities(1.0, tp, ed);
        System.err.println(new Vector(tp));
    }
    return ed;
}

Example 24 with Vector

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

the class ComplexSubstitutionModel method printLastProbabilityMatrix.

public void printLastProbabilityMatrix() {
    getLogLikelihood();
    System.err.println((probability == null) ? "Null probability vector" : "Not null probability vector");
    if (probability == null) {
        boolean test = BayesianStochasticSearchVariableSelection.Utils.connectedAndWellConditioned(probability, this);
        System.err.println("BSSVS valid = " + test);
    }
    System.err.println(new Vector(probability));
}

Example 25 with Vector

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

the class HKY method main.

public static void main(String[] args) {
    //        double kappa = 2.0;
    //        double[] pi = new double[]{0.15,0.30,0.20,0.35};
    //        double time = 0.1;
    double kappa = 1.0;
    double[] pi = new double[] { 0.25, 0.25, 0.25, 0.25 };
    double time = 0.1;
    FrequencyModel freqModel = new FrequencyModel(Nucleotides.INSTANCE, pi);
    HKY hky = new HKY(kappa, freqModel);
    EigenDecomposition decomp = hky.getEigenDecomposition();
    //        Matrix evec = new Matrix(decomp.getEigenVectors());
    //        Matrix ivec = new Matrix(decomp.getInverseEigenVectors());
    //        System.out.println("Evec =\n"+evec);
    //        System.out.println("Ivec =\n"+ivec);
    Vector eval = new Vector(decomp.getEigenValues());
    System.out.println("Eval = " + eval);
    double[] probs = new double[16];
    hky.getTransitionProbabilities(time, probs);
    System.out.println("new probs = " + new Vector(probs));
    // check against old implementation
    dr.oldevomodel.substmodel.FrequencyModel oldFreq = new dr.oldevomodel.substmodel.FrequencyModel(Nucleotides.INSTANCE, pi);
    dr.oldevomodel.substmodel.HKY oldHKY = new dr.oldevomodel.substmodel.HKY(kappa, oldFreq);
    oldHKY.setKappa(kappa);
    oldHKY.getTransitionProbabilities(time, probs);
    System.out.println("old probs = " + new Vector(probs));
}