Examples with WilsonBalding beast.evolution.operators.WilsonBalding used on opensource projects

Example 1 with WilsonBalding

use of beast.evolution.operators.WilsonBalding 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]) + ")");
    }
}