Examples with MultivariateTraitTree dr.evolution.tree.MultivariateTraitTree used on opensource projects

Example 11 with MultivariateTraitTree

use of dr.evolution.tree.MultivariateTraitTree in project beast-mcmc by beast-dev.

the class MultivariateTraitUtils method computeTreeTraitMean.

public static double[] computeTreeTraitMean(FullyConjugateMultivariateTraitLikelihood trait, double[] rootValue, boolean conditionOnRoot) {
    double[] root = trait.getPriorMean();
    if (conditionOnRoot) {
        System.err.println("WARNING: Not yet fully implemented (conditioning on root in simulator)");
        //root = new double[root.length];
        root = rootValue;
    }
    final int nTaxa = trait.getTreeModel().getExternalNodeCount();
    double[] mean = new double[root.length * nTaxa];
    for (int i = 0; i < nTaxa; ++i) {
        System.arraycopy(root, 0, mean, i * root.length, root.length);
    }
    if (trait.driftModels != null) {
        MultivariateTraitTree myTreeModel = trait.getTreeModel();
        for (int i = 0; i < nTaxa; ++i) {
            double[] shiftContribution = getShiftContributionToMean(myTreeModel.getExternalNode(i), trait);
            for (int j = 0; j < trait.dimTrait; ++j) {
                mean[i * trait.dimTrait + j] = mean[i * trait.dimTrait + j] + shiftContribution[j];
            }
        }
    }
    return mean;
}

Example 12 with MultivariateTraitTree

use of dr.evolution.tree.MultivariateTraitTree in project beast-mcmc by beast-dev.

the class MultivariateTraitUtils method computeTreeVariance.

public static double[][] computeTreeVariance(FullyConjugateMultivariateTraitLikelihood trait, boolean conditionOnRoot) {
    MultivariateTraitTree treeModel = trait.getTreeModel();
    final int tipCount = treeModel.getExternalNodeCount();
    double[][] variance = new double[tipCount][tipCount];
    for (int i = 0; i < tipCount; i++) {
        // Fill in diagonal
        double marginalTime = trait.getRescaledLengthToRoot(treeModel.getExternalNode(i));
        variance[i][i] = marginalTime;
        for (int j = i + 1; j < tipCount; j++) {
            NodeRef mrca = findMRCA(trait, i, j);
            if (DEBUG) {
                System.err.println(trait.getTreeModel().getRoot().getNumber());
                System.err.print("Taxa pair: " + i + " : " + j + " (" + mrca.getNumber() + ") = ");
            }
            double length = trait.getRescaledLengthToRoot(mrca);
            if (DEBUG) {
                System.err.println(length);
            }
            variance[i][j] = length;
        }
    }
    // Make symmetric
    for (int i = 0; i < tipCount; i++) {
        for (int j = i + 1; j < tipCount; j++) {
            variance[j][i] = variance[i][j];
        }
    }
    if (DEBUG) {
        System.err.println("");
        System.err.println("New tree conditional variance:\n" + new Matrix(variance));
    }
    if (!conditionOnRoot) {
        double priorSampleSize = trait.getPriorSampleSize();
        for (int i = 0; i < tipCount; ++i) {
            for (int j = 0; j < tipCount; ++j) {
                variance[i][j] += 1.0 / priorSampleSize;
            }
        }
        if (DEBUG) {
            System.err.println("");
            System.err.println("New tree unconditional variance:\n" + new Matrix(variance));
        }
    }
    return variance;
}