Search in sources :

Example 26 with Sequence

use of dr.evolution.sequence.Sequence in project beast-mcmc by beast-dev.

the class AncestralStateBeagleTreeLikelihoodTest method testJointLikelihood.

public void testJointLikelihood() {
    TreeModel treeModel = new TreeModel("treeModel", tree);
    Sequence[] sequence = new Sequence[3];
    sequence[0] = new Sequence(new Taxon("0"), "A");
    sequence[1] = new Sequence(new Taxon("1"), "C");
    sequence[2] = new Sequence(new Taxon("2"), "C");
    Taxa taxa = new Taxa();
    for (Sequence s : sequence) {
        taxa.addTaxon(s.getTaxon());
    }
    SimpleAlignment alignment = new SimpleAlignment();
    for (Sequence s : sequence) {
        alignment.addSequence(s);
    }
    Parameter mu = new Parameter.Default(1, 1.0);
    Parameter kappa = new Parameter.Default(1, 1.0);
    double[] pi = { 0.25, 0.25, 0.25, 0.25 };
    Parameter freqs = new Parameter.Default(pi);
    FrequencyModel f = new FrequencyModel(Nucleotides.INSTANCE, freqs);
    HKY hky = new HKY(kappa, f);
    GammaSiteRateModel siteRateModel = new GammaSiteRateModel("gammaModel", mu, null, -1, null);
    siteRateModel.setSubstitutionModel(hky);
    BranchModel branchModel = new HomogeneousBranchModel(siteRateModel.getSubstitutionModel());
    BranchRateModel branchRateModel = null;
    AncestralStateBeagleTreeLikelihood treeLikelihood = new AncestralStateBeagleTreeLikelihood(alignment, treeModel, branchModel, siteRateModel, branchRateModel, null, false, PartialsRescalingScheme.DEFAULT, true, null, hky.getDataType(), "stateTag", // useMap = true
    true, false);
    double logLike = treeLikelihood.getLogLikelihood();
    StringBuffer buffer = new StringBuffer();
    //        Tree.Utils.newick(treeModel, treeModel.getRoot(), false, Tree.BranchLengthType.LENGTHS_AS_TIME,
    //                null, null, new NodeAttributeProvider[]{treeLikelihood}, null, null, buffer);
    TreeUtils.newick(treeModel, treeModel.getRoot(), false, TreeUtils.BranchLengthType.LENGTHS_AS_TIME, null, null, new TreeTraitProvider[] { treeLikelihood }, null, buffer);
    System.out.println(buffer);
    System.out.println("t_CA(2) = " + t(false, 2.0));
    System.out.println("t_CC(1) = " + t(true, 1.0));
    double trueValue = 0.25 * t(false, 2.0) * Math.pow(t(true, 1.0), 3.0);
    assertEquals(logLike, Math.log(trueValue), 1e-6);
}
Also used : FrequencyModel(dr.evomodel.substmodel.FrequencyModel) Taxon(dr.evolution.util.Taxon) HomogeneousBranchModel(dr.evomodel.branchmodel.HomogeneousBranchModel) GammaSiteRateModel(dr.evomodel.siteratemodel.GammaSiteRateModel) Sequence(dr.evolution.sequence.Sequence) BranchModel(dr.evomodel.branchmodel.BranchModel) HomogeneousBranchModel(dr.evomodel.branchmodel.HomogeneousBranchModel) TreeModel(dr.evomodel.tree.TreeModel) Taxa(dr.evolution.util.Taxa) SimpleAlignment(dr.evolution.alignment.SimpleAlignment) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) HKY(dr.evomodel.substmodel.nucleotide.HKY) Parameter(dr.inference.model.Parameter) AncestralStateBeagleTreeLikelihood(dr.evomodel.treelikelihood.AncestralStateBeagleTreeLikelihood)

Example 27 with Sequence

use of dr.evolution.sequence.Sequence in project beast-mcmc by beast-dev.

the class PhylipSequentialImporter method importAlignment.

/**
	 * importAlignment. 
	 */
public Alignment importAlignment() throws IOException, Importer.ImportException {
    SimpleAlignment alignment = null;
    try {
        int taxonCount = readInteger();
        int siteCount = readInteger();
        String firstSeq = null;
        for (int i = 0; i < taxonCount; i++) {
            StringBuffer name = new StringBuffer();
            char ch = read();
            int n = 0;
            while (!Character.isWhitespace(ch) && (maxNameLength < 1 || n < maxNameLength)) {
                name.append(ch);
                ch = read();
                n++;
            }
            StringBuffer seq = new StringBuffer(siteCount);
            readSequence(seq, dataType, "", siteCount, "-", "?", ".", firstSeq);
            if (firstSeq == null) {
                firstSeq = seq.toString();
            }
            if (alignment == null) {
                alignment = new SimpleAlignment();
            }
            alignment.addSequence(new Sequence(new Taxon(name.toString()), seq.toString()));
        }
    } catch (EOFException e) {
    }
    return alignment;
}
Also used : SimpleAlignment(dr.evolution.alignment.SimpleAlignment) Taxon(dr.evolution.util.Taxon) EOFException(java.io.EOFException) Sequence(dr.evolution.sequence.Sequence)

Example 28 with Sequence

use of dr.evolution.sequence.Sequence in project beast-mcmc by beast-dev.

the class DataLikelihoodTester2 method createAlignment.

private static SimpleAlignment createAlignment(Object[][] taxa_sequence, DataType dataType) {
    SimpleAlignment alignment = new SimpleAlignment();
    alignment.setDataType(dataType);
    //        alignment.setDataType(Nucleotides.INSTANCE);
    // 6, 17
    Taxon[] taxa = new Taxon[taxa_sequence[0].length];
    System.out.println("Taxon len = " + taxa_sequence[0].length);
    System.out.println("Alignment len = " + taxa_sequence[1].length);
    if (taxa_sequence.length > 2)
        System.out.println("Date len = " + taxa_sequence[2].length);
    for (int i = 0; i < taxa_sequence[0].length; i++) {
        taxa[i] = new Taxon(taxa_sequence[0][i].toString());
        if (taxa_sequence.length > 2) {
            Date date = new Date((Double) taxa_sequence[2][i], Units.Type.YEARS, (Boolean) taxa_sequence[3][0]);
            taxa[i].setDate(date);
        }
        //taxonList.addTaxon(taxon);
        Sequence sequence = new Sequence(taxa_sequence[1][i].toString());
        sequence.setTaxon(taxa[i]);
        sequence.setDataType(dataType);
        alignment.addSequence(sequence);
    }
    return alignment;
}
Also used : SimpleAlignment(dr.evolution.alignment.SimpleAlignment) Taxon(dr.evolution.util.Taxon) Sequence(dr.evolution.sequence.Sequence) Date(dr.evolution.util.Date)

Example 29 with Sequence

use of dr.evolution.sequence.Sequence in project beast-mcmc by beast-dev.

the class SequenceLikelihoodTest method createAlignmentWithAllUniquePatterns.

protected void createAlignmentWithAllUniquePatterns(Object[][] taxa_sequence, DataType dataType) {
    alignment = new SimpleAlignment();
    alignment.setDataType(dataType);
    int nTaxa = taxa_sequence[0].length;
    String[] allUniquePatterns = createAllUniquePatterns(nTaxa, dataType);
    taxa_sequence[1] = allUniquePatterns;
    // 6, 17
    taxa = new Taxon[nTaxa];
    System.out.println("Taxon len = " + taxa_sequence[0].length);
    System.out.println("Alignment len = " + taxa_sequence[1].length);
    if (taxa_sequence.length > 2)
        System.out.println("Date len = " + taxa_sequence[2].length);
    for (int i = 0; i < taxa_sequence[0].length; i++) {
        taxa[i] = new Taxon(taxa_sequence[0][i].toString());
        if (taxa_sequence.length > 2) {
            Date date = new Date((Double) taxa_sequence[2][i], Units.Type.YEARS, (Boolean) taxa_sequence[3][0]);
            taxa[i].setDate(date);
        }
        //taxonList.addTaxon(taxon);
        Sequence sequence = new Sequence(taxa_sequence[1][i].toString());
        sequence.setTaxon(taxa[i]);
        sequence.setDataType(dataType);
        alignment.addSequence(sequence);
    }
    System.out.println("Sequence pattern count = " + alignment.getPatternCount());
}
Also used : SimpleAlignment(dr.evolution.alignment.SimpleAlignment) Taxon(dr.evolution.util.Taxon) Sequence(dr.evolution.sequence.Sequence) Date(dr.evolution.util.Date)

Example 30 with Sequence

use of dr.evolution.sequence.Sequence in project beast-mcmc by beast-dev.

the class AncestralStateTreeLikelihoodTest method testJointLikelihood.

public void testJointLikelihood() {
    TreeModel treeModel = new TreeModel("treeModel", tree);
    Sequence[] sequence = new Sequence[3];
    sequence[0] = new Sequence(new Taxon("0"), "A");
    sequence[1] = new Sequence(new Taxon("1"), "C");
    sequence[2] = new Sequence(new Taxon("2"), "C");
    Taxa taxa = new Taxa();
    for (Sequence s : sequence) {
        taxa.addTaxon(s.getTaxon());
    }
    SimpleAlignment alignment = new SimpleAlignment();
    for (Sequence s : sequence) {
        alignment.addSequence(s);
    }
    Parameter mu = new Parameter.Default(1, 1.0);
    Parameter kappa = new Parameter.Default(1, 1.0);
    double[] pi = { 0.25, 0.25, 0.25, 0.25 };
    Parameter freqs = new Parameter.Default(pi);
    FrequencyModel f = new FrequencyModel(Nucleotides.INSTANCE, freqs);
    HKY hky = new HKY(kappa, f);
    AncestralStateTreeLikelihood treeLikelihood = new AncestralStateTreeLikelihood(alignment, treeModel, new GammaSiteModel(hky), new StrictClockBranchRates(mu), false, true, Nucleotides.INSTANCE, "state", false, // useMap = true
    true, false);
    double logLike = treeLikelihood.getLogLikelihood();
    StringBuffer buffer = new StringBuffer();
    TreeUtils.newick(treeModel, treeModel.getRoot(), false, TreeUtils.BranchLengthType.LENGTHS_AS_TIME, null, null, new TreeTraitProvider[] { treeLikelihood }, null, buffer);
    System.out.println(buffer);
    System.out.println("t_CA(2) = " + t(false, 2.0));
    System.out.println("t_CC(1) = " + t(true, 1.0));
    double trueValue = 0.25 * t(false, 2.0) * Math.pow(t(true, 1.0), 3.0);
    assertEquals(logLike, Math.log(trueValue), 1e-6);
}
Also used : FrequencyModel(dr.oldevomodel.substmodel.FrequencyModel) Taxon(dr.evolution.util.Taxon) AncestralStateTreeLikelihood(dr.oldevomodel.treelikelihood.AncestralStateTreeLikelihood) Sequence(dr.evolution.sequence.Sequence) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) TreeModel(dr.evomodel.tree.TreeModel) Taxa(dr.evolution.util.Taxa) GammaSiteModel(dr.oldevomodel.sitemodel.GammaSiteModel) SimpleAlignment(dr.evolution.alignment.SimpleAlignment) HKY(dr.oldevomodel.substmodel.HKY) Parameter(dr.inference.model.Parameter)

Aggregations

Sequence (dr.evolution.sequence.Sequence)31 SimpleAlignment (dr.evolution.alignment.SimpleAlignment)15 Taxon (dr.evolution.util.Taxon)12 DataType (dr.evolution.datatype.DataType)6 BranchRateModel (dr.evomodel.branchratemodel.BranchRateModel)6 GammaSiteRateModel (dr.evomodel.siteratemodel.GammaSiteRateModel)5 TreeModel (dr.evomodel.tree.TreeModel)5 Parameter (dr.inference.model.Parameter)5 Tree (dr.evolution.tree.Tree)4 Date (dr.evolution.util.Date)4 HomogeneousBranchModel (dr.evomodel.branchmodel.HomogeneousBranchModel)4 DefaultBranchRateModel (dr.evomodel.branchratemodel.DefaultBranchRateModel)4 FrequencyModel (dr.evomodel.substmodel.FrequencyModel)4 Partition (dr.app.beagle.tools.Partition)3 NewickImporter (dr.evolution.io.NewickImporter)3 HKY (dr.evomodel.substmodel.nucleotide.HKY)3 ArrayList (java.util.ArrayList)3 BeagleSequenceSimulator (dr.app.beagle.tools.BeagleSequenceSimulator)2 ImportException (dr.evolution.io.Importer.ImportException)2 Taxa (dr.evolution.util.Taxa)2