Search in sources :

Example 16 with TreeLikelihood

use of beast.evolution.likelihood.TreeLikelihood in project beast2 by CompEvol.

the class SequenceSimulator method main.

// printUsageAndExit
@SuppressWarnings("unchecked")
public static void main(String[] args) {
    try {
        // parse arguments
        if (args.length < 2) {
            printUsageAndExit();
        }
        String fileName = args[0];
        int replications = Integer.parseInt(args[1]);
        PrintStream out = System.out;
        if (args.length == 3) {
            File file = new File(args[2]);
            out = new PrintStream(file);
        }
        // grab the file
        String xml = "";
        BufferedReader fin = new BufferedReader(new FileReader(fileName));
        while (fin.ready()) {
            xml += fin.readLine();
        }
        fin.close();
        // parse the xml
        XMLParser parser = new XMLParser();
        BEASTInterface beastObject = parser.parseFragment(xml, true);
        // find relevant objects from the model
        TreeLikelihood treeLikelihood = getTreeLikelihood(beastObject);
        if (treeLikelihood == null) {
            throw new IllegalArgumentException("No treelikelihood found in file. Giving up now.");
        }
        Alignment data = ((Input<Alignment>) treeLikelihood.getInput("data")).get();
        Tree tree = ((Input<Tree>) treeLikelihood.getInput("tree")).get();
        SiteModel pSiteModel = ((Input<SiteModel>) treeLikelihood.getInput("siteModel")).get();
        BranchRateModel pBranchRateModel = ((Input<BranchRateModel>) treeLikelihood.getInput("branchRateModel")).get();
        // feed to sequence simulator and generate leaves
        SequenceSimulator treeSimulator = new SequenceSimulator();
        treeSimulator.init(data, tree, pSiteModel, pBranchRateModel, replications);
        XMLProducer producer = new XMLProducer();
        Alignment alignment = treeSimulator.simulate();
        xml = producer.toRawXML(alignment);
        out.println("<beast version='2.0'>");
        out.println(xml);
        out.println("</beast>");
    } catch (Exception e) {
        e.printStackTrace();
    }
}
Also used : PrintStream(java.io.PrintStream) XMLProducer(beast.util.XMLProducer) TreeLikelihood(beast.evolution.likelihood.TreeLikelihood) SiteModel(beast.evolution.sitemodel.SiteModel) XMLParserException(beast.util.XMLParserException) IOException(java.io.IOException) FileNotFoundException(java.io.FileNotFoundException) Alignment(beast.evolution.alignment.Alignment) Input(beast.core.Input) BranchRateModel(beast.evolution.branchratemodel.BranchRateModel) BufferedReader(java.io.BufferedReader) Tree(beast.evolution.tree.Tree) FileReader(java.io.FileReader) BEASTInterface(beast.core.BEASTInterface) XMLParser(beast.util.XMLParser) File(java.io.File)

Example 17 with TreeLikelihood

use of beast.evolution.likelihood.TreeLikelihood in project beast2 by CompEvol.

the class TreeLikelihoodTest method testHKY85GILikelihood.

@Test
public void testHKY85GILikelihood() throws Exception {
    // Set up HKY85+G+I model: estimated freqs, kappa = 39.464538, 4 gamma categories, shape = 0.587649, prop invariant = 0.486548
    Alignment data = BEASTTestCase.getAlignment();
    Tree tree = BEASTTestCase.getTree(data);
    Frequencies freqs = new Frequencies();
    freqs.initByName("data", data);
    HKY hky = new HKY();
    hky.initByName("kappa", "39.464538", "frequencies", freqs);
    SiteModel siteModel = new SiteModel();
    siteModel.initByName("mutationRate", "1.0", "gammaCategoryCount", 4, "shape", "0.587649", "proportionInvariant", "0.486548", "substModel", hky);
    TreeLikelihood likelihood = newTreeLikelihood();
    likelihood.initByName("data", data, "tree", tree, "siteModel", siteModel);
    double logP = 0;
    logP = likelihood.calculateLogP();
    assertEquals(logP, -1789.639227747059, BEASTTestCase.PRECISION);
    likelihood.initByName("useAmbiguities", true, "data", data, "tree", tree, "siteModel", siteModel);
    logP = likelihood.calculateLogP();
    assertEquals(logP, -1789.639227747059, BEASTTestCase.PRECISION);
}
Also used : Alignment(beast.evolution.alignment.Alignment) HKY(beast.evolution.substitutionmodel.HKY) BeagleTreeLikelihood(beast.evolution.likelihood.BeagleTreeLikelihood) TreeLikelihood(beast.evolution.likelihood.TreeLikelihood) Tree(beast.evolution.tree.Tree) SiteModel(beast.evolution.sitemodel.SiteModel) Frequencies(beast.evolution.substitutionmodel.Frequencies) UncertainAlignmentTest(test.beast.evolution.alignment.UncertainAlignmentTest) Test(org.junit.Test)

Example 18 with TreeLikelihood

use of beast.evolution.likelihood.TreeLikelihood in project beast2 by CompEvol.

the class TreeLikelihoodTest method testMarginalisationOfLikelihoodBinary.

@Test
public void testMarginalisationOfLikelihoodBinary() throws Exception {
    // test summation over all patterns adds to 1 for binary data
    Sequence German_ST = new Sequence("German_ST", "           10110010");
    Sequence Dutch_List = new Sequence("Dutch_List", "          11010100");
    Sequence English_ST = new Sequence("English_ST", "          11101000");
    Alignment data = new Alignment();
    data.initByName("sequence", German_ST, "sequence", Dutch_List, "sequence", English_ST, "dataType", "binary");
    Tree tree = BEASTTestCase.getTree(data, "(English_ST:0.22743347188019544,(German_ST:0.10557648379843088,Dutch_List:0.10557648379843088):0.12185698808176457):0.0;");
    RealParameter frequencies = new RealParameter("0.683 0.317");
    Frequencies freqs = new Frequencies();
    freqs.initByName("frequencies", frequencies);
    GeneralSubstitutionModel covarion = new GeneralSubstitutionModel();
    covarion.initByName("frequencies", freqs, "rates", new RealParameter("1.0 1.0"));
    SiteModel siteModel = new SiteModel();
    siteModel.initByName("mutationRate", "1.0", "gammaCategoryCount", 1, "substModel", covarion);
    TreeLikelihood likelihood = newTreeLikelihood();
    likelihood.initByName("data", data, "tree", tree, "siteModel", siteModel);
    likelihood.initByName("useAmbiguities", false, "data", data, "tree", tree, "siteModel", siteModel);
    likelihood.calculateLogP();
    double[] logPs = likelihood.getPatternLogLikelihoods();
    double P = 0;
    for (double d : logPs) {
        P += Math.exp(d);
    }
    assertEquals(P, 1.0, BEASTTestCase.PRECISION);
}
Also used : Alignment(beast.evolution.alignment.Alignment) BeagleTreeLikelihood(beast.evolution.likelihood.BeagleTreeLikelihood) TreeLikelihood(beast.evolution.likelihood.TreeLikelihood) Tree(beast.evolution.tree.Tree) RealParameter(beast.core.parameter.RealParameter) GeneralSubstitutionModel(beast.evolution.substitutionmodel.GeneralSubstitutionModel) SiteModel(beast.evolution.sitemodel.SiteModel) Sequence(beast.evolution.alignment.Sequence) Frequencies(beast.evolution.substitutionmodel.Frequencies) UncertainAlignmentTest(test.beast.evolution.alignment.UncertainAlignmentTest) Test(org.junit.Test)

Example 19 with TreeLikelihood

use of beast.evolution.likelihood.TreeLikelihood in project beast2 by CompEvol.

the class TreeLikelihoodTest method testBeagleRNALikelihood.

/**
 * Test only effective when BEAGLE installed - otherwise it will always
 * pass since BEAGLE code will never be run.
 *
 * @throws Exception
 */
@Test
public void testBeagleRNALikelihood() throws Exception {
    Sequence seq1 = new Sequence("t1", "GUACGUACGUAC");
    Sequence seq2 = new Sequence("t2", "UACGUACGUACG");
    Sequence seq3 = new Sequence("t3", "ACGUACGUACGU");
    Alignment data = new Alignment();
    data.initByName("sequence", seq1, "sequence", seq2, "sequence", seq3, "dataType", "nucleotide");
    Tree tree = BEASTTestCase.getTree(data, "((t1:0.5,t2:0.5):0.5,t3:1.0):0.0;");
    SiteModel siteModel = new SiteModel();
    siteModel.initByName("gammaCategoryCount", 1, "substModel", new JukesCantor());
    TreeLikelihood likelihoodNoBeagle = newTreeLikelihood();
    likelihoodNoBeagle.initByName("data", data, "tree", tree, "siteModel", siteModel);
    double logLnoBeagle = likelihoodNoBeagle.calculateLogP();
    System.setProperty("java.only", "false");
    TreeLikelihood likelihoodBeagle = new TreeLikelihood();
    likelihoodBeagle.initByName("data", data, "tree", tree, "siteModel", siteModel);
    double logLBeagle = likelihoodBeagle.calculateLogP();
    assertEquals(logLBeagle, logLnoBeagle, BEASTTestCase.PRECISION);
}
Also used : Alignment(beast.evolution.alignment.Alignment) BeagleTreeLikelihood(beast.evolution.likelihood.BeagleTreeLikelihood) TreeLikelihood(beast.evolution.likelihood.TreeLikelihood) Tree(beast.evolution.tree.Tree) SiteModel(beast.evolution.sitemodel.SiteModel) Sequence(beast.evolution.alignment.Sequence) JukesCantor(beast.evolution.substitutionmodel.JukesCantor) UncertainAlignmentTest(test.beast.evolution.alignment.UncertainAlignmentTest) Test(org.junit.Test)

Example 20 with TreeLikelihood

use of beast.evolution.likelihood.TreeLikelihood in project beast2 by CompEvol.

the class TreeLikelihoodTest method testSDolloLikelihood.

@Test
public void testSDolloLikelihood() throws Exception {
    UserDataType dataType = new UserDataType();
    dataType.initByName("states", 2, "codeMap", "0=1, 1=0, ?=0 1, -=0 1");
    Alignment data = new Alignment();
    Sequence German_ST = new Sequence("German_ST", BEASTTestCase.German_ST.dataInput.get());
    Sequence Dutch_List = new Sequence("Dutch_List", BEASTTestCase.Dutch_List.dataInput.get());
    ;
    Sequence English_ST = new Sequence("English_ST", BEASTTestCase.English_ST.dataInput.get());
    ;
    Sequence French = new Sequence("French", BEASTTestCase.French.dataInput.get());
    ;
    Sequence Italian = new Sequence("Italian", BEASTTestCase.Italian.dataInput.get());
    ;
    Sequence Spanish = new Sequence("Spanish", BEASTTestCase.Spanish.dataInput.get());
    ;
    data.initByName("sequence", German_ST, "sequence", Dutch_List, "sequence", English_ST, "sequence", French, "sequence", Italian, "sequence", Spanish, "userDataType", dataType);
    Tree tree = BEASTTestCase.getTree(data, "((English_ST:0.22743347188019544,(German_ST:0.10557648379843088,Dutch_List:0.10557648379843088):0.12185698808176457):1.5793160946109988,(Spanish:0.11078392189606047,(Italian:0.10119772534558173,French:0.10119772534558173):0.009586196550478737):1.6959656445951337)");
    RealParameter frequencies = new RealParameter("1 0");
    Frequencies freqs = new Frequencies();
    freqs.initByName("frequencies", frequencies);
    RealParameter deathprob = new RealParameter("1.7");
    MutationDeathModel SDollo = new MutationDeathModel();
    SDollo.initByName("deathprob", deathprob, "frequencies", freqs);
    SiteModel siteModel = new SiteModel();
    siteModel.initByName("mutationRate", "1.0", "gammaCategoryCount", 1, "substModel", SDollo);
    TreeLikelihood likelihood = newTreeLikelihood();
    likelihood.initByName("data", data, "tree", tree, "siteModel", siteModel);
    double logP = 0;
    likelihood.initByName("data", data, "tree", tree, "siteModel", siteModel, "useAmbiguities", true);
    logP = likelihood.calculateLogP();
    // beast1 xml gives -3551.6436
    assertEquals(logP, -3551.6436270344648, BEASTTestCase.PRECISION);
}
Also used : Alignment(beast.evolution.alignment.Alignment) UserDataType(beast.evolution.datatype.UserDataType) BeagleTreeLikelihood(beast.evolution.likelihood.BeagleTreeLikelihood) TreeLikelihood(beast.evolution.likelihood.TreeLikelihood) Tree(beast.evolution.tree.Tree) RealParameter(beast.core.parameter.RealParameter) MutationDeathModel(beast.evolution.substitutionmodel.MutationDeathModel) SiteModel(beast.evolution.sitemodel.SiteModel) Sequence(beast.evolution.alignment.Sequence) Frequencies(beast.evolution.substitutionmodel.Frequencies) UncertainAlignmentTest(test.beast.evolution.alignment.UncertainAlignmentTest) Test(org.junit.Test)

Aggregations

TreeLikelihood (beast.evolution.likelihood.TreeLikelihood)24 Alignment (beast.evolution.alignment.Alignment)22 SiteModel (beast.evolution.sitemodel.SiteModel)22 Tree (beast.evolution.tree.Tree)22 BeagleTreeLikelihood (beast.evolution.likelihood.BeagleTreeLikelihood)20 Test (org.junit.Test)17 UncertainAlignmentTest (test.beast.evolution.alignment.UncertainAlignmentTest)16 Frequencies (beast.evolution.substitutionmodel.Frequencies)13 HKY (beast.evolution.substitutionmodel.HKY)6 RealParameter (beast.core.parameter.RealParameter)5 GeneralSubstitutionModel (beast.evolution.substitutionmodel.GeneralSubstitutionModel)5 Sequence (beast.evolution.alignment.Sequence)4 JukesCantor (beast.evolution.substitutionmodel.JukesCantor)3 BEASTInterface (beast.core.BEASTInterface)2 MarginalTree (bacter.MarginalTree)1 Region (bacter.Region)1 SimulatedAlignment (beast.app.seqgen.SimulatedAlignment)1 Input (beast.core.Input)1 TaxonSet (beast.evolution.alignment.TaxonSet)1 BranchRateModel (beast.evolution.branchratemodel.BranchRateModel)1