Search in sources :

Example 6 with TypeSet

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

the class MigrationModelLogger method init.

@Override
public void init(PrintStream out) {
    String outName;
    TypeSet typeSet = migModel.getTypeSet();
    if (migModel.getID() == null || migModel.getID().matches("\\s*"))
        outName = "migModel";
    else
        outName = migModel.getID();
    out.print(outName + ".popSizeScaleFactor\t");
    for (int i = 0; i < migModel.getNTypes(); i++) {
        if (mtTree != null)
            out.print(outName + ".popSize_" + typeSet.getTypeName(i) + "\t");
        else
            out.print(outName + ".popSize_" + i + "\t");
    }
    out.print(outName + ".rateMatrixScaleFactor\t");
    for (int i = 0; i < migModel.getNTypes(); i++) {
        for (int j = 0; j < migModel.getNTypes(); j++) {
            if (i == j)
                continue;
            if (mtTree != null)
                out.format("%s.rateMatrix_%s_%s\t", outName, typeSet.getTypeName(i), typeSet.getTypeName(j));
            else
                out.format("%s.rateMatrix_%d_%d\t", outName, i, j);
        }
    }
    if (migModel.rateMatrixFlagsInput.get() != null) {
        for (int i = 0; i < migModel.getNTypes(); i++) {
            for (int j = 0; j < migModel.getNTypes(); j++) {
                if (i == j)
                    continue;
                if (mtTree != null)
                    out.format("%s.rateMatrixFlag_%s_%s\t", outName, typeSet.getTypeName(i), typeSet.getTypeName(j));
                else
                    out.format("%s.rateMatrixFlag_%d_%d\t", outName, i, j);
            }
        }
    }
}
Also used : TypeSet(beast.evolution.tree.TypeSet)

Example 7 with TypeSet

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

the class Ewing_Test method test2.

@Test
public void test2() throws Exception {
    System.out.println("Ewing_test 2");
    // Fix seed.
    Randomizer.setSeed(42);
    // Assemble initial MultiTypeTree
    String newickStr = "(((1[&deme=1]:0.5)[&deme=0]:0.5,2[&deme=0]:1)[&deme=0]:1," + "3[&deme=0]:2)[&deme=0]:0;";
    MultiTypeTreeFromNewick mtTree = new MultiTypeTreeFromNewick();
    mtTree.initByName("value", newickStr, "typeLabel", "deme");
    // 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 distribution:
    StructuredCoalescentTreeDensity distribution = new StructuredCoalescentTreeDensity();
    distribution.initByName("migrationModel", migModel, "multiTypeTree", mtTree, "checkValidity", true);
    // Set up state:
    State state = new State();
    state.initByName("stateNode", mtTree);
    // Set up operators:
    TypedWilsonBaldingEasy twbOperator = new TypedWilsonBaldingEasy();
    twbOperator.initByName("weight", 1.0, "migrationModel", migModel, "multiTypeTree", mtTree);
    TypedSubtreeExchangeEasy tsxOperator = new TypedSubtreeExchangeEasy();
    tsxOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "isNarrow", true);
    MultiTypeUniform mtuOperator = new MultiTypeUniform();
    mtuOperator.initByName("weight", 1.0, "migrationModel", migModel, "multiTypeTree", mtTree);
    MultiTypeTreeScale mttsOperator = new MultiTypeTreeScale();
    mttsOperator.initByName("weight", 1.0, "scaleFactor", 0.8, "useOldTreeScaler", true, "migrationModel", migModel, "multiTypeTree", mtTree);
    TypePairBirthDeath tpbdOperator = new TypePairBirthDeath();
    tpbdOperator.initByName("weight", 1.0, "migrationModel", migModel, "multiTypeTree", mtTree);
    TypeMergeSplit tmsOperator = new TypeMergeSplit();
    tmsOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "includeRoot", true);
    // 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", twbOperator, "operator", tsxOperator, "operator", mtuOperator, "operator", mttsOperator, "operator", tpbdOperator, "operator", tmsOperator, "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() > 1000) && (Math.abs(logger.getHeightMean() - 23) < 0.5) && (Math.abs(logger.getHeightVar() - 300) < 30.0);
    Assert.assertTrue(withinTol);
}
Also used : MultiTypeTreeStatLogger(multitypetree.util.MultiTypeTreeStatLogger) MCMC(beast.core.MCMC) RealParameter(beast.core.parameter.RealParameter) SCMigrationModel(beast.evolution.tree.SCMigrationModel) State(beast.core.State) MultiTypeTreeFromNewick(beast.evolution.tree.MultiTypeTreeFromNewick) TypeSet(beast.evolution.tree.TypeSet) StructuredCoalescentTreeDensity(multitypetree.distributions.StructuredCoalescentTreeDensity) Test(org.junit.Test)

Example 8 with TypeSet

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

the class Ewing_Test method test1.

@Test
public void test1() throws Exception {
    System.out.println("Ewing_test1");
    // Fix seed.
    Randomizer.setSeed(1);
    // Assemble initial MultiTypeTree
    String newickStr = "((1[&deme=0]:1,2[&deme=0]:1)[&deme=0]:1," + "3[&deme=0]:2)[&deme=0]:0;";
    MultiTypeTreeFromNewick mtTree = new MultiTypeTreeFromNewick();
    mtTree.initByName("value", newickStr, "typeLabel", "deme");
    // 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 distribution:
    StructuredCoalescentTreeDensity distribution = new StructuredCoalescentTreeDensity();
    distribution.initByName("migrationModel", migModel, "multiTypeTree", mtTree);
    // Set up state:
    State state = new State();
    state.initByName("stateNode", mtTree);
    // Set up operators:
    TypedWilsonBaldingEasy twbOperator = new TypedWilsonBaldingEasy();
    twbOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    TypedSubtreeExchangeEasy tsxOperator = new TypedSubtreeExchangeEasy();
    tsxOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "isNarrow", true);
    MultiTypeUniform mtuOperator = new MultiTypeUniform();
    mtuOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    MultiTypeTreeScale mttsOperator = new MultiTypeTreeScale();
    mttsOperator.initByName("weight", 1.0, "scaleFactor", 0.8, "useOldTreeScaler", true, "multiTypeTree", mtTree, "migrationModel", migModel);
    TypePairBirthDeath tpbdOperator = new TypePairBirthDeath();
    tpbdOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel);
    TypeMergeSplit tmsOperator = new TypeMergeSplit();
    tmsOperator.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "includeRoot", true);
    // Set up stat analysis logger:
    MultiTypeTreeStatLogger logger = new MultiTypeTreeStatLogger();
    logger.initByName("multiTypeTree", mtTree, "burninFrac", 0.1, "logEvery", 1000);
    // Logger fileLogger = new Logger();
    // TreeHeightLogger thLogger = new TreeHeightLogger();
    // thLogger.initByName("tree", mtTree);
    // fileLogger.initByName(
    // "fileName", "test.log",
    // "logEvery", 1000,
    // "model", distribution,
    // "log", thLogger);
    // Set up MCMC:
    MCMC mcmc = new MCMC();
    mcmc.initByName("chainLength", "10000000", "state", state, "distribution", distribution, "operator", twbOperator, "operator", tsxOperator, "operator", mtuOperator, "operator", mttsOperator, "operator", tpbdOperator, "operator", tmsOperator, "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() > 3000) && (Math.abs(logger.getHeightMean() - 19.15) < 0.5) && (Math.abs(logger.getHeightVar() - 310) < 20);
    Assert.assertTrue(withinTol);
}
Also used : MultiTypeTreeStatLogger(multitypetree.util.MultiTypeTreeStatLogger) MCMC(beast.core.MCMC) RealParameter(beast.core.parameter.RealParameter) SCMigrationModel(beast.evolution.tree.SCMigrationModel) State(beast.core.State) MultiTypeTreeFromNewick(beast.evolution.tree.MultiTypeTreeFromNewick) TypeSet(beast.evolution.tree.TypeSet) StructuredCoalescentTreeDensity(multitypetree.distributions.StructuredCoalescentTreeDensity) Test(org.junit.Test)

Example 9 with TypeSet

use of beast.evolution.tree.TypeSet 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)

Example 10 with TypeSet

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

the class TWB_TS_Test method test1.

@Test
public void test1() throws Exception {
    System.out.println("TWB_test 1");
    // Fix seed.
    Randomizer.setSeed(42);
    // Assemble initial MultiTypeTree
    String newickStr = "((1[&deme=0]:1,2[&deme=0]:1)[&deme=0]:1," + "3[&deme=0]:2)[&deme=0]:0;";
    MultiTypeTreeFromNewick mtTree = new MultiTypeTreeFromNewick();
    mtTree.initByName("value", newickStr, "typeLabel", "deme");
    // 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 distribution:
    StructuredCoalescentTreeDensity distribution = new StructuredCoalescentTreeDensity();
    distribution.initByName("migrationModel", migModel, "multiTypeTree", mtTree);
    // Set up state:
    State state = new State();
    state.initByName("stateNode", mtTree);
    // Set up operator:
    TypedWilsonBalding operatorTWB = new TypedWilsonBalding();
    operatorTWB.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "alpha", 0.2);
    Operator operatorMTTS = new MultiTypeTreeScale();
    operatorMTTS.initByName("weight", 1.0, "multiTypeTree", mtTree, "migrationModel", migModel, "scaleFactor", 0.8, "useOldTreeScaler", false);
    // Set up stat analysis logger:
    MultiTypeTreeStatLogger logger = new MultiTypeTreeStatLogger();
    logger.initByName("multiTypeTree", mtTree, "burninFrac", 0.2, "logEvery", 1000);
    // Set up MCMC:
    MCMC mcmc = new MCMC();
    mcmc.initByName("chainLength", "1000000", "state", state, "distribution", distribution, "operator", operatorTWB, "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("0 0 0"));
    double simHeightMean = DiscreteStatistics.mean(heights);
    double simHeightVar = DiscreteStatistics.variance(heights);
    // Compare results with simulation results:
    boolean withinTol = (logger.getHeightESS() > 400) && (Math.abs(logger.getHeightMean() - simHeightMean) < 1.0) && (Math.abs(logger.getHeightVar() - simHeightVar) < 30);
    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) SCMigrationModel(beast.evolution.tree.SCMigrationModel) State(beast.core.State) MultiTypeTreeFromNewick(beast.evolution.tree.MultiTypeTreeFromNewick) TypeSet(beast.evolution.tree.TypeSet) StructuredCoalescentTreeDensity(multitypetree.distributions.StructuredCoalescentTreeDensity) Test(org.junit.Test)

Aggregations

TypeSet (beast.evolution.tree.TypeSet)12 RealParameter (beast.core.parameter.RealParameter)11 SCMigrationModel (beast.evolution.tree.SCMigrationModel)11 StructuredCoalescentTreeDensity (multitypetree.distributions.StructuredCoalescentTreeDensity)10 Test (org.junit.Test)10 MCMC (beast.core.MCMC)9 State (beast.core.State)9 MultiTypeTreeStatLogger (multitypetree.util.MultiTypeTreeStatLogger)9 Operator (beast.core.Operator)7 MultiTypeTreeFromNewick (beast.evolution.tree.MultiTypeTreeFromNewick)7 IntegerParameter (beast.core.parameter.IntegerParameter)4 MultiTypeTree (beast.evolution.tree.MultiTypeTree)3 StructuredCoalescentMultiTypeTree (beast.evolution.tree.StructuredCoalescentMultiTypeTree)3