Search in sources :

Example 1 with StructuredCoalescentMultiTypeTree

use of beast.evolution.tree.StructuredCoalescentMultiTypeTree in project MultiTypeTree by tgvaughan.

the class NSR_Test method test.

@Test
public void test() throws Exception {
    System.out.println("NSR test");
    // Fix seed.
    Randomizer.setSeed(42);
    // Assemble migration model:
    RealParameter rateMatrix = new RealParameter("0.1 0.1");
    RealParameter popSizes = new RealParameter("7.0 7.0");
    SCMigrationModel migModel = new SCMigrationModel();
    migModel.initByName("rateMatrix", rateMatrix, "popSizes", popSizes, "typeSet", new TypeSet("A", "B"));
    // Assemble initial MultiTypeTree
    MultiTypeTree mtTree = new StructuredCoalescentMultiTypeTree();
    mtTree.initByName("typeLabel", "deme", "migrationModel", migModel, "leafTypes", "1 0");
    // Set up state:
    State state = new State();
    state.initByName("stateNode", mtTree);
    // Assemble distribution:
    StructuredCoalescentTreeDensity distribution = new StructuredCoalescentTreeDensity();
    distribution.initByName("migrationModel", migModel, "multiTypeTree", mtTree);
    // Set up operators:
    Operator operatorNSR = new NodeShiftRetype();
    operatorNSR.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    // Set up stat analysis logger:
    MultiTypeTreeStatLogger logger = new MultiTypeTreeStatLogger();
    logger.initByName("multiTypeTree", mtTree, "burninFrac", 0.1, "logEvery", 100);
    // Set up MCMC:
    MCMC mcmc = new MCMC();
    mcmc.initByName("chainLength", "1000000", "state", state, "distribution", distribution, "operator", operatorNSR, "logger", logger);
    // Run MCMC:
    mcmc.run();
    System.out.format("height mean = %s\n", logger.getHeightMean());
    System.out.format("height var = %s\n", logger.getHeightVar());
    System.out.format("height ESS = %s\n", logger.getHeightESS());
    // Compare analysis results with truth:
    boolean withinTol = (logger.getHeightESS() > 2000) && (Math.abs(logger.getHeightMean() - 19.0) < 0.5) && (Math.abs(logger.getHeightVar() - 291) < 30);
    Assert.assertTrue(withinTol);
}
Also used : Operator(beast.core.Operator) MultiTypeTreeStatLogger(multitypetree.util.MultiTypeTreeStatLogger) State(beast.core.State) TypeSet(beast.evolution.tree.TypeSet) StructuredCoalescentTreeDensity(multitypetree.distributions.StructuredCoalescentTreeDensity) MCMC(beast.core.MCMC) RealParameter(beast.core.parameter.RealParameter) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree) MultiTypeTree(beast.evolution.tree.MultiTypeTree) SCMigrationModel(beast.evolution.tree.SCMigrationModel) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree) Test(org.junit.Test)

Example 2 with StructuredCoalescentMultiTypeTree

use of beast.evolution.tree.StructuredCoalescentMultiTypeTree in project MultiTypeTree by tgvaughan.

the class STXR_NRR_MTU_TS_Test method test.

@Test
public void test() throws Exception {
    System.out.println("STXR_NRR_MTU_TS test");
    // Fix seed.
    Randomizer.setSeed(42);
    // Assemble migration model:
    RealParameter rateMatrix = new RealParameter("0.1 0.1");
    RealParameter popSizes = new RealParameter("7.0 7.0");
    SCMigrationModel migModel = new SCMigrationModel();
    migModel.initByName("rateMatrix", rateMatrix, "popSizes", popSizes, "typeSet", new TypeSet("A", "B"));
    // Assemble initial MultiTypeTree
    MultiTypeTree mtTree = new StructuredCoalescentMultiTypeTree();
    mtTree.initByName("typeLabel", "deme", "migrationModel", migModel, "leafTypes", "1 1 0 0");
    // Set up state:
    State state = new State();
    state.initByName("stateNode", mtTree);
    // Assemble distribution:
    StructuredCoalescentTreeDensity distribution = new StructuredCoalescentTreeDensity();
    distribution.initByName("migrationModel", migModel, "multiTypeTree", mtTree);
    // Set up operators:
    Operator operatorSTXR = new TypedSubtreeExchangeRandom();
    operatorSTXR.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "mu", 0.2);
    Operator operatorNRR = new NodeRetypeRandom();
    operatorNRR.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "mu", 0.2);
    Operator operatorMTU = new MultiTypeUniform();
    operatorMTU.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    Operator operatorMTTS = new MultiTypeTreeScale();
    operatorMTTS.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "scaleFactor", 1.5, "useOldTreeScaler", false);
    // Set up stat analysis logger:
    MultiTypeTreeStatLogger logger = new MultiTypeTreeStatLogger();
    logger.initByName("multiTypeTree", mtTree, "burninFrac", 0.1, "logEvery", 1000);
    // Set up MCMC:
    MCMC mcmc = new MCMC();
    mcmc.initByName("chainLength", "10000000", "state", state, "distribution", distribution, "operator", operatorSTXR, "operator", operatorNRR, "operator", operatorMTU, "operator", operatorMTTS, "logger", logger);
    // Run MCMC:
    mcmc.run();
    System.out.format("height mean = %s\n", logger.getHeightMean());
    System.out.format("height var = %s\n", logger.getHeightVar());
    System.out.format("height ESS = %s\n", logger.getHeightESS());
    // Compare analysis results with truth:
    boolean withinTol = (logger.getHeightESS() > 2000) && (Math.abs(logger.getHeightMean() - 25.8) < 0.5) && (Math.abs(logger.getHeightVar() - 320) < 30);
    Assert.assertTrue(withinTol);
}
Also used : Operator(beast.core.Operator) MultiTypeTreeStatLogger(multitypetree.util.MultiTypeTreeStatLogger) MCMC(beast.core.MCMC) RealParameter(beast.core.parameter.RealParameter) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree) MultiTypeTree(beast.evolution.tree.MultiTypeTree) SCMigrationModel(beast.evolution.tree.SCMigrationModel) State(beast.core.State) TypeSet(beast.evolution.tree.TypeSet) StructuredCoalescentTreeDensity(multitypetree.distributions.StructuredCoalescentTreeDensity) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree) Test(org.junit.Test)

Example 3 with StructuredCoalescentMultiTypeTree

use of beast.evolution.tree.StructuredCoalescentMultiTypeTree in project MultiTypeTree by tgvaughan.

the class InitMigrationModelConnector method customConnector.

public static boolean customConnector(BeautiDoc doc) {
    for (BEASTInterface p : doc.getPartitions("Tree")) {
        TreeLikelihood treeLikelihood = (TreeLikelihood) p;
        StructuredCoalescentMultiTypeTree tree = (StructuredCoalescentMultiTypeTree) treeLikelihood.treeInput.get();
        String pID = BeautiDoc.parsePartition(tree.getID());
        SCMigrationModel migModel = (SCMigrationModel) doc.pluginmap.get("migModel.t:" + pID);
        SCMigrationModel migModelInit = (SCMigrationModel) doc.pluginmap.get("migModelInit.t:" + pID);
        String rateMatrixStr = getParameterString((RealParameter) migModel.rateMatrixInput.get());
        String popSizesStr = getParameterString((RealParameter) migModel.popSizesInput.get());
        // Ensure model has appropriate number of demes
        int uniqueTraitCount = uniqueTraitsInData(tree).size();
        StringBuilder rateMatrixStrBuilder = new StringBuilder();
        StringBuilder popSizesStrBuilder = new StringBuilder();
        migModel.getTypeSet().initAndValidate();
        if (migModel.popSizesInput.get().getDimension() != migModel.getNTypes()) {
            for (int i = 0; i < migModel.getNTypes(); i++) {
                popSizesStrBuilder.append(" 1.0");
                for (int j = 0; j < migModel.getNTypes(); j++) {
                    if (j == i)
                        continue;
                    rateMatrixStrBuilder.append(" 1.0");
                }
            }
            popSizesStr = popSizesStrBuilder.toString();
            rateMatrixStr = rateMatrixStrBuilder.toString();
            ((RealParameter) migModel.popSizesInput.get()).setDimension(migModel.getNTypes());
            ((RealParameter) migModel.popSizesInput.get()).valuesInput.setValue(popSizesStr, (RealParameter) migModel.popSizesInput.get());
            ((RealParameter) migModel.rateMatrixInput.get()).setDimension(migModel.getNTypes() * (migModel.getNTypes() - 1));
            ((RealParameter) migModel.rateMatrixInput.get()).valuesInput.setValue(rateMatrixStr, (RealParameter) migModel.rateMatrixInput.get());
            ((RealParameter) migModel.popSizesInput.get()).initAndValidate();
            ((RealParameter) migModel.rateMatrixInput.get()).initAndValidate();
            migModel.initAndValidate();
        }
        ((RealParameter) migModelInit.popSizesInput.get()).setDimension(migModel.getNTypes());
        ((RealParameter) migModelInit.popSizesInput.get()).valuesInput.setValue(popSizesStr, (RealParameter) migModelInit.popSizesInput.get());
        ((RealParameter) migModelInit.rateMatrixInput.get()).setDimension(migModel.getNTypes() * (migModel.getNTypes() - 1));
        ((RealParameter) migModelInit.rateMatrixInput.get()).valuesInput.setValue(rateMatrixStr, (RealParameter) migModelInit.rateMatrixInput.get());
        try {
            ((RealParameter) migModelInit.popSizesInput.get()).initAndValidate();
            ((RealParameter) migModelInit.rateMatrixInput.get()).initAndValidate();
            migModelInit.initAndValidate();
        } catch (Exception ex) {
            System.err.println("Error configuring initial migration model.");
        }
    }
    return false;
}
Also used : TreeLikelihood(beast.evolution.likelihood.TreeLikelihood) RealParameter(beast.core.parameter.RealParameter) BEASTInterface(beast.core.BEASTInterface) SCMigrationModel(beast.evolution.tree.SCMigrationModel) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree)

Example 4 with StructuredCoalescentMultiTypeTree

use of beast.evolution.tree.StructuredCoalescentMultiTypeTree in project MultiTypeTree by tgvaughan.

the class STX_NR_MTU_TS_Test method test.

@Test
public void test() throws Exception {
    System.out.println("STX_NR_MTU_TS test");
    // Test passing locally, not on Travis.  WHY!?
    // Fix seed.
    Randomizer.setSeed(53);
    // Assemble migration model:
    RealParameter rateMatrix = new RealParameter("0.1 0.1");
    RealParameter popSizes = new RealParameter("7.0 7.0");
    SCMigrationModel migModel = new SCMigrationModel();
    migModel.initByName("rateMatrix", rateMatrix, "popSizes", popSizes, "typeSet", new TypeSet("A", "B"));
    // Assemble initial MultiTypeTree
    MultiTypeTree mtTree = new StructuredCoalescentMultiTypeTree();
    mtTree.initByName("typeLabel", "deme", "migrationModel", migModel, "leafTypes", "1 1 0 0");
    // Set up state:
    State state = new State();
    state.initByName("stateNode", mtTree);
    // Assemble distribution:
    StructuredCoalescentTreeDensity distribution = new StructuredCoalescentTreeDensity();
    distribution.initByName("migrationModel", migModel, "multiTypeTree", mtTree);
    // Set up operators:
    Operator operatorSTX = new TypedSubtreeExchange();
    operatorSTX.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    Operator operatorNR = new NodeRetype();
    operatorNR.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    Operator operatorMTU = new MultiTypeUniform();
    operatorMTU.initByName("weight", 1.0, "migrationModel", migModel, "multiTypeTree", mtTree);
    Operator operatorMTTS = new MultiTypeTreeScale();
    operatorMTTS.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "scaleFactor", 1.5, "useOldTreeScaler", false);
    // Set up stat analysis logger:
    MultiTypeTreeStatLogger logger = new MultiTypeTreeStatLogger();
    logger.initByName("multiTypeTree", mtTree, "burninFrac", 0.1, "logEvery", 1000);
    // Set up MCMC:
    MCMC mcmc = new MCMC();
    mcmc.initByName("chainLength", "1000000", "state", state, "distribution", distribution, "operator", operatorSTX, "operator", operatorNR, "operator", operatorMTU, "operator", operatorMTTS, "logger", logger);
    // Run MCMC:
    mcmc.run();
    System.out.format("height mean = %s\n", logger.getHeightMean());
    System.out.format("height var = %s\n", logger.getHeightVar());
    System.out.format("height ESS = %s\n", logger.getHeightESS());
    // Direct simulation:
    double[] heights = UtilMethods.getSimulatedHeights(migModel, new IntegerParameter("1 1 0 0"));
    double simHeightMean = DiscreteStatistics.mean(heights);
    double simHeightVar = DiscreteStatistics.variance(heights);
    // Compare analysis results with truth:
    boolean withinTol = (logger.getHeightESS() > 500) && (Math.abs(logger.getHeightMean() - simHeightMean) < 2.0) && (Math.abs(logger.getHeightVar() - simHeightVar) < 50);
    Assert.assertTrue(withinTol);
}
Also used : Operator(beast.core.Operator) MultiTypeTreeStatLogger(multitypetree.util.MultiTypeTreeStatLogger) IntegerParameter(beast.core.parameter.IntegerParameter) MCMC(beast.core.MCMC) RealParameter(beast.core.parameter.RealParameter) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree) MultiTypeTree(beast.evolution.tree.MultiTypeTree) SCMigrationModel(beast.evolution.tree.SCMigrationModel) State(beast.core.State) TypeSet(beast.evolution.tree.TypeSet) StructuredCoalescentTreeDensity(multitypetree.distributions.StructuredCoalescentTreeDensity) StructuredCoalescentMultiTypeTree(beast.evolution.tree.StructuredCoalescentMultiTypeTree) Test(org.junit.Test)

Aggregations

RealParameter (beast.core.parameter.RealParameter)4 SCMigrationModel (beast.evolution.tree.SCMigrationModel)4 StructuredCoalescentMultiTypeTree (beast.evolution.tree.StructuredCoalescentMultiTypeTree)4 MCMC (beast.core.MCMC)3 Operator (beast.core.Operator)3 State (beast.core.State)3 MultiTypeTree (beast.evolution.tree.MultiTypeTree)3 TypeSet (beast.evolution.tree.TypeSet)3 StructuredCoalescentTreeDensity (multitypetree.distributions.StructuredCoalescentTreeDensity)3 MultiTypeTreeStatLogger (multitypetree.util.MultiTypeTreeStatLogger)3 Test (org.junit.Test)3 BEASTInterface (beast.core.BEASTInterface)1 IntegerParameter (beast.core.parameter.IntegerParameter)1 TreeLikelihood (beast.evolution.likelihood.TreeLikelihood)1