Examples with RandomTree beast.evolution.tree.RandomTree used on opensource projects

Example 1 with RandomTree

use of beast.evolution.tree.RandomTree in project beast2 by CompEvol.

the class WilsonBaldingTest method topologyDistribution.

/**
 * Test topology distribution.
 * @throws Exception
 */
@Test
public void topologyDistribution() throws Exception {
    // Fix seed: will hopefully ensure success of test unless something
    // goes terribly wrong.
    Randomizer.setSeed(42);
    // Assemble model:
    ConstantPopulation constantPop = new ConstantPopulation();
    constantPop.initByName("popSize", new RealParameter("10000.0"));
    List<Object> alignmentInitArgs = new ArrayList<Object>();
    for (int i = 0; i < 4; i++) {
        Sequence thisSeq = new Sequence();
        thisSeq.initByName("taxon", String.valueOf(i), "value", "?");
        alignmentInitArgs.add("sequence");
        alignmentInitArgs.add(thisSeq);
    }
    Alignment alignment = new Alignment();
    alignment.initByName(alignmentInitArgs.toArray());
    Tree tree = new RandomTree();
    tree.initByName("taxa", alignment, "populationModel", constantPop);
    TreeIntervals treeIntervals = new TreeIntervals();
    treeIntervals.initByName("tree", tree);
    Coalescent coalescentDistrib = new Coalescent();
    coalescentDistrib.initByName("treeIntervals", treeIntervals, "populationModel", constantPop);
    // Set up state:
    State state = new State();
    state.initByName("stateNode", tree);
    // Set up operator:
    WilsonBalding wilsonBalding = new WilsonBalding();
    wilsonBalding.initByName("weight", "1", "tree", tree);
    // Set up logger:
    TreeReport treeReport = new TreeReport();
    treeReport.initByName("logEvery", "100", "burnin", "200000", "credibleSetPercentage", "95.0", "log", tree, "silent", true);
    // Set up MCMC:
    MCMC mcmc = new MCMC();
    mcmc.initByName("chainLength", "2000000", "state", state, "distribution", coalescentDistrib, "operator", wilsonBalding, "logger", treeReport);
    // Run MCMC:
    mcmc.run();
    // Obtain analysis results:
    TreeTraceAnalysis analysis = treeReport.getAnalysis();
    Map<String, Integer> topologyCounts = analysis.getTopologyCounts();
    int totalTreesUsed = analysis.getNTrees();
    // Test topology distribution against ideal:
    double tol = 0.005;
    for (int i = 0; i < topologies.length; i++) {
        double thisProb = topologyCounts.get(topologies[i]) / (double) totalTreesUsed;
        boolean withinTol = (thisProb > probs[i] - tol && thisProb < probs[i] + tol);
        Assert.assertTrue(withinTol);
        System.err.format("Topology %s rel. freq. %.3f", topologies[i], thisProb);
        if (withinTol)
            System.err.println(" (Within tolerance " + tol + " of " + String.valueOf(probs[i]) + ")");
        else
            System.err.println(" (FAILURE: outside tolerance " + tol + " of " + String.valueOf(probs[i]) + ")");
    }
}

Example 2 with RandomTree

use of beast.evolution.tree.RandomTree in project beast2 by CompEvol.

the class RandomTreeTest method testCoalescentTimes.

@Test
public void testCoalescentTimes() throws Exception {
    Randomizer.setSeed(53);
    int Nleaves = 10;
    int Niter = 5000;
    // (Serially sampled) coalescent time means and variances
    // estimated from 50000 trees simulated using MASTER
    double[] coalTimeMeansTruth = { 1.754662, 2.833337, 3.843532, 4.850805, 5.849542, 6.847016, 7.8482, 8.855137, 10.15442 };
    double[] coalTimeVarsTruth = { 0.2751625, 0.2727121, 0.2685172, 0.2705117, 0.2678611, 0.2671793, 0.2686952, 0.2828477, 1.076874 };
    // Assemble BEASTObjects needed by RandomTree
    StringBuilder traitSB = new StringBuilder();
    List<Sequence> seqList = new ArrayList<Sequence>();
    for (int i = 0; i < Nleaves; i++) {
        String taxonID = "t " + i;
        seqList.add(new Sequence(taxonID, "?"));
        if (i > 0)
            traitSB.append(",");
        traitSB.append(taxonID).append("=").append(i);
    }
    Alignment alignment = new Alignment(seqList, "nucleotide");
    TaxonSet taxonSet = new TaxonSet(alignment);
    TraitSet timeTrait = new TraitSet();
    timeTrait.initByName("traitname", "date-backward", "taxa", taxonSet, "value", traitSB.toString());
    ConstantPopulation popFunc = new ConstantPopulation();
    popFunc.initByName("popSize", new RealParameter("1.0"));
    // Create RandomTree and TreeInterval instances
    RandomTree tree = new RandomTree();
    TreeIntervals intervals = new TreeIntervals();
    // Estimate coalescence time moments
    double[] coalTimeMeans = new double[Nleaves - 1];
    double[] coalTimeVars = new double[Nleaves - 1];
    double[] coalTimes = new double[Nleaves - 1];
    for (int i = 0; i < Niter; i++) {
        tree.initByName("taxa", alignment, "populationModel", popFunc, "trait", timeTrait);
        intervals.initByName("tree", tree);
        intervals.getCoalescentTimes(coalTimes);
        for (int j = 0; j < Nleaves - 1; j++) {
            coalTimeMeans[j] += coalTimes[j];
            coalTimeVars[j] += coalTimes[j] * coalTimes[j];
        }
    }
    // Normalise means and variances
    for (int j = 0; j < Nleaves - 1; j++) {
        coalTimeMeans[j] /= Niter;
        coalTimeVars[j] /= Niter;
        coalTimeVars[j] -= coalTimeMeans[j] * coalTimeMeans[j];
    }
    // Test means and variances against independently estimated values
    for (int j = 0; j < Nleaves - 1; j++) {
        assert (relError(coalTimeMeans[j], coalTimeMeansTruth[j]) < 5e-3);
        assert (relError(coalTimeVars[j], coalTimeVarsTruth[j]) < 5e-2);
    }
}

Example 3 with RandomTree

use of beast.evolution.tree.RandomTree in project beast2 by CompEvol.

the class UnorderedAlignmentsTest method getTree.

public static Tree getTree(TaxonSet taxa) throws Exception {
    Tree tree = new RandomTree();
    TraitSet dates = getDates(taxa);
    ConstantPopulation constant = new ConstantPopulation();
    constant.initByName("popSize", new RealParameter("5.0"));
    tree.initByName("taxonset", taxa, "populationModel", constant, "trait", dates);
    return tree;
}

Example 4 with RandomTree

use of beast.evolution.tree.RandomTree in project beast2 by CompEvol.

the class StarBeastStartState method initWithMRCACalibrations.

private void initWithMRCACalibrations(List<MRCAPrior> calibrations) {
    final Tree spTree = speciesTreeInput.get();
    final RandomTree rnd = new RandomTree();
    rnd.setInputValue("taxonset", spTree.getTaxonset());
    for (final MRCAPrior cal : calibrations) {
        rnd.setInputValue("constraint", cal);
    }
    ConstantPopulation pf = new ConstantPopulation();
    pf.setInputValue("popSize", new RealParameter("1.0"));
    rnd.setInputValue("populationModel", pf);
    rnd.initAndValidate();
    spTree.assignFromWithoutID((Tree) rnd);
    final double rootHeight = spTree.getRoot().getHeight();
    randomInitGeneTrees(rootHeight);
}