Search in sources :

Example 16 with StrictClockBranchRates

use of dr.evomodel.branchratemodel.StrictClockBranchRates in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodDiagonalOUBMInd.

public void testLikelihoodDiagonalOUBMInd() {
    System.out.println("\nTest Likelihood using Diagonal OU / BM:");
    // Diffusion
    List<BranchRateModel> optimalTraitsModels = new ArrayList<BranchRateModel>();
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.1", new double[] { 1.0 })));
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.2", new double[] { -3.0 })));
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.3", new double[] { -2.0 })));
    DiagonalMatrix strengthOfSelectionMatrixParamOUBM = new DiagonalMatrix(new Parameter.Default(new double[] { 0.0, 0.0, 50.0 }));
    DiagonalMatrix strengthOfSelectionMatrixParamOU = new DiagonalMatrix(new Parameter.Default(new double[] { 10.0, 20.0, 50.0 }));
    DiagonalMatrix diffusionPrecisionMatrixParameter = new DiagonalMatrix(new Parameter.Default(new double[] { 1.0, 2.0, 3.0 }));
    MultivariateDiffusionModel diffusionModel = new MultivariateDiffusionModel(diffusionPrecisionMatrixParameter);
    DiffusionProcessDelegate diffusionProcessDelegateOUBM = new OUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParamOUBM));
    DiffusionProcessDelegate diffusionProcessDelegateOU = new OUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParamOU));
    DiffusionProcessDelegate diffusionProcessDelegateBM = new HomogeneousDiffusionModelDelegate(treeModel, diffusionModel);
    // CDL
    ContinuousDataLikelihoodDelegate likelihoodDelegateOUBM = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegateOUBM, dataModel, rootPriorInf, rateTransformation, rateModel, false);
    ContinuousDataLikelihoodDelegate likelihoodDelegateOU = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegateOU, dataModel, rootPriorInf, rateTransformation, rateModel, false);
    ContinuousDataLikelihoodDelegate likelihoodDelegateBM = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegateBM, dataModel, rootPriorInf, rateTransformation, rateModel, false);
    // Likelihood Computation
    TreeDataLikelihood dataLikelihoodOUBM = new TreeDataLikelihood(likelihoodDelegateOUBM, treeModel, rateModel);
    TreeDataLikelihood dataLikelihoodOU = new TreeDataLikelihood(likelihoodDelegateOU, treeModel, rateModel);
    TreeDataLikelihood dataLikelihoodBM = new TreeDataLikelihood(likelihoodDelegateBM, treeModel, rateModel);
    // Conditional simulations
    MathUtils.setSeed(17890826);
    double[] traitsOUBM = getConditionalSimulations(dataLikelihoodOUBM, likelihoodDelegateOUBM, diffusionModel, dataModel, rootPriorInf, treeModel, rateTransformation);
    System.err.println(new Vector(traitsOUBM));
    MathUtils.setSeed(17890826);
    double[] traitsOU = getConditionalSimulations(dataLikelihoodOU, likelihoodDelegateOU, diffusionModel, dataModel, rootPriorInf, treeModel, rateTransformation);
    System.err.println(new Vector(traitsOU));
    MathUtils.setSeed(17890826);
    double[] traitsBM = getConditionalSimulations(dataLikelihoodBM, likelihoodDelegateBM, diffusionModel, dataModel, rootPriorInf, treeModel, rateTransformation);
    System.err.println(new Vector(traitsBM));
    // Check that missing dimensions with the same process have the same values
    assertEquals(format.format(traitsBM[3]), format.format(traitsOUBM[3]));
    assertEquals(format.format(traitsBM[4]), format.format(traitsOUBM[4]));
    assertEquals(format.format(traitsBM[7]), format.format(traitsOUBM[7]));
    assertEquals(format.format(traitsOU[5]), format.format(traitsOUBM[5]));
}
Also used : ArrayList(java.util.ArrayList) MultivariateElasticModel(dr.evomodel.continuous.MultivariateElasticModel) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) TreeDataLikelihood(dr.evomodel.treedatalikelihood.TreeDataLikelihood) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) DiagonalMatrix(dr.inference.model.DiagonalMatrix) MultivariateDiffusionModel(dr.evomodel.continuous.MultivariateDiffusionModel) MatrixParameter(dr.inference.model.MatrixParameter) Parameter(dr.inference.model.Parameter) Vector(dr.math.matrixAlgebra.Vector)

Example 17 with StrictClockBranchRates

use of dr.evomodel.branchratemodel.StrictClockBranchRates in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodFullDiagonalOUFactor.

public void testLikelihoodFullDiagonalOUFactor() {
    System.out.println("\nTest Likelihood comparing full and diagonal OU and factor:");
    // Diffusion
    List<BranchRateModel> optimalTraitsModels = new ArrayList<BranchRateModel>();
    ArbitraryBranchRates.BranchRateTransform transform = make(false, false, false);
    optimalTraitsModels.add(new ArbitraryBranchRates(treeModel, new Parameter.Default("rate.1", new double[] { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }), transform, false));
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.2", new double[] { 1.5 })));
    Parameter[] strengthOfSelectionParameters = new Parameter[2];
    strengthOfSelectionParameters[0] = new Parameter.Default(new double[] { 0.5, 0.0 });
    strengthOfSelectionParameters[1] = new Parameter.Default(new double[] { 0.0, 1.5 });
    MatrixParameter strengthOfSelectionMatrixParam = new MatrixParameter("strengthOfSelectionMatrix", strengthOfSelectionParameters);
    DiagonalMatrix strengthOfSelectionMatrixParamDiagonal = new DiagonalMatrix(new Parameter.Default(new double[] { 0.5, 1.5 }));
    DiffusionProcessDelegate diffusionProcessDelegate = new OUDiffusionModelDelegate(treeModel, diffusionModelFactor, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParam));
    DiffusionProcessDelegate diffusionProcessDelegateDiagonal = new OUDiffusionModelDelegate(treeModel, diffusionModelFactor, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParamDiagonal));
    // CDL
    ContinuousDataLikelihoodDelegate likelihoodDelegateFactors = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModelFactor, rootPriorFactor, rateTransformation, rateModel, false);
    ContinuousDataLikelihoodDelegate likelihoodDelegateFactorsDiagonal = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegateDiagonal, dataModelFactor, rootPriorFactor, rateTransformation, rateModel, false);
    dataModelFactor.setLikelihoodDelegate(likelihoodDelegateFactors);
    // Likelihood Computation
    TreeDataLikelihood dataLikelihoodFactors = new TreeDataLikelihood(likelihoodDelegateFactors, treeModel, rateModel);
    TreeDataLikelihood dataLikelihoodFactorsDiagonal = new TreeDataLikelihood(likelihoodDelegateFactorsDiagonal, treeModel, rateModel);
    double likelihoodFactorData = dataLikelihoodFactors.getLogLikelihood();
    double likelihoodFactorDiffusion = dataModelFactor.getLogLikelihood();
    double likelihoodFactorDataDiagonal = dataLikelihoodFactorsDiagonal.getLogLikelihood();
    double likelihoodFactorDiffusionDiagonal = dataModelFactor.getLogLikelihood();
    assertEquals("likelihoodFullDiagonalOUFactor", format.format(likelihoodFactorData + likelihoodFactorDiffusion), format.format(likelihoodFactorDataDiagonal + likelihoodFactorDiffusionDiagonal));
}
Also used : MatrixParameter(dr.inference.model.MatrixParameter) ArrayList(java.util.ArrayList) MultivariateElasticModel(dr.evomodel.continuous.MultivariateElasticModel) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) ArbitraryBranchRates(dr.evomodel.branchratemodel.ArbitraryBranchRates) TreeDataLikelihood(dr.evomodel.treedatalikelihood.TreeDataLikelihood) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) DiagonalMatrix(dr.inference.model.DiagonalMatrix) MatrixParameter(dr.inference.model.MatrixParameter) Parameter(dr.inference.model.Parameter)

Example 18 with StrictClockBranchRates

use of dr.evomodel.branchratemodel.StrictClockBranchRates in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodDiagonalOUBM.

public void testLikelihoodDiagonalOUBM() {
    System.out.println("\nTest Likelihood using Diagonal OU / BM:");
    // Diffusion
    List<BranchRateModel> optimalTraitsModels = new ArrayList<BranchRateModel>();
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.1", new double[] { 1.0 })));
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.2", new double[] { 2.0 })));
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.3", new double[] { -2.0 })));
    DiagonalMatrix strengthOfSelectionMatrixParam = new DiagonalMatrix(new Parameter.Default(new double[] { 0.0, 0.000001, 50.0 }));
    DiffusionProcessDelegate diffusionProcessDelegate = new OUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParam));
    // CDL
    ContinuousDataLikelihoodDelegate likelihoodDelegate = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModel, rootPrior, rateTransformation, rateModel, false);
    // Likelihood Computation
    TreeDataLikelihood dataLikelihood = new TreeDataLikelihood(likelihoodDelegate, treeModel, rateModel);
    testLikelihood("likelihoodDiagonalOUBM", dataLikelihood);
    // Conditional moments (preorder)
    testConditionalMoments(dataLikelihood, likelihoodDelegate);
    // Fixed Root
    ContinuousDataLikelihoodDelegate likelihoodDelegateInf = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModel, rootPriorInf, rateTransformation, rateModel, true);
    TreeDataLikelihood dataLikelihoodInf = new TreeDataLikelihood(likelihoodDelegateInf, treeModel, rateModel);
    testLikelihood("likelihoodDiagonalOUBMInf", dataLikelihoodInf);
    testConditionalMoments(dataLikelihoodInf, likelihoodDelegateInf);
}
Also used : ArrayList(java.util.ArrayList) MultivariateElasticModel(dr.evomodel.continuous.MultivariateElasticModel) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) TreeDataLikelihood(dr.evomodel.treedatalikelihood.TreeDataLikelihood) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) DiagonalMatrix(dr.inference.model.DiagonalMatrix) MatrixParameter(dr.inference.model.MatrixParameter) Parameter(dr.inference.model.Parameter)

Example 19 with StrictClockBranchRates

use of dr.evomodel.branchratemodel.StrictClockBranchRates in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodFullAndDiagonalOU.

public void testLikelihoodFullAndDiagonalOU() {
    System.out.println("\nTest Likelihood comparing Full and Diagonal OU:");
    // Diffusion
    List<BranchRateModel> optimalTraitsModels = new ArrayList<BranchRateModel>();
    ArbitraryBranchRates.BranchRateTransform transform = make(false, false, false);
    optimalTraitsModels.add(new ArbitraryBranchRates(treeModel, new Parameter.Default("rate.1", new double[] { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }), transform, false));
    optimalTraitsModels.add(new ArbitraryBranchRates(treeModel, new Parameter.Default("rate.2", new double[] { 0, -1, 2, -3, 4, -5, 6, -7, 8, -9 }), transform, false));
    optimalTraitsModels.add(new StrictClockBranchRates(new Parameter.Default("rate.3", new double[] { -2.0 })));
    Parameter[] strengthOfSelectionParameters = new Parameter[3];
    strengthOfSelectionParameters[0] = new Parameter.Default(new double[] { 0.5, 0.0, 0.0 });
    strengthOfSelectionParameters[1] = new Parameter.Default(new double[] { 0.0, 10.5, 0.0 });
    strengthOfSelectionParameters[2] = new Parameter.Default(new double[] { 0.0, 0.0, 100.0 });
    MatrixParameter strengthOfSelectionMatrixParam = new MatrixParameter("strengthOfSelectionMatrix", strengthOfSelectionParameters);
    DiagonalMatrix strengthOfSelectionMatrixParamDiagonal = new DiagonalMatrix(new Parameter.Default(new double[] { 0.5, 10.5, 100.0 }));
    DiffusionProcessDelegate diffusionProcessDelegate = new OUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParam));
    DiffusionProcessDelegate diffusionProcessDelegateDiagonal = new OUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParamDiagonal));
    // CDL
    ContinuousDataLikelihoodDelegate likelihoodDelegate = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModel, rootPrior, rateTransformation, rateModel, false);
    ContinuousDataLikelihoodDelegate likelihoodDelegateDiagonal = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegateDiagonal, dataModel, rootPrior, rateTransformation, rateModel, false);
    // Likelihood Computation
    TreeDataLikelihood dataLikelihood = new TreeDataLikelihood(likelihoodDelegate, treeModel, rateModel);
    TreeDataLikelihood dataLikelihoodDiagonal = new TreeDataLikelihood(likelihoodDelegateDiagonal, treeModel, rateModel);
    assertEquals("likelihoodFullDiagonalOU", format.format(dataLikelihood.getLogLikelihood()), format.format(dataLikelihoodDiagonal.getLogLikelihood()));
}
Also used : MatrixParameter(dr.inference.model.MatrixParameter) ArrayList(java.util.ArrayList) MultivariateElasticModel(dr.evomodel.continuous.MultivariateElasticModel) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) ArbitraryBranchRates(dr.evomodel.branchratemodel.ArbitraryBranchRates) TreeDataLikelihood(dr.evomodel.treedatalikelihood.TreeDataLikelihood) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) DiagonalMatrix(dr.inference.model.DiagonalMatrix) MatrixParameter(dr.inference.model.MatrixParameter) Parameter(dr.inference.model.Parameter)

Example 20 with StrictClockBranchRates

use of dr.evomodel.branchratemodel.StrictClockBranchRates in project beast-mcmc by beast-dev.

the class MsatSamplingTreeLikelihoodTest method setUp.

public void setUp() throws Exception {
    super.setUp();
    // taxa
    ArrayList<Taxon> taxonList3 = new ArrayList<Taxon>();
    Collections.addAll(taxonList3, new Taxon("Taxon1"), new Taxon("Taxon2"), new Taxon("Taxon3"), new Taxon("Taxon4"), new Taxon("Taxon5"), new Taxon("Taxon6"), new Taxon("Taxon7"));
    Taxa taxa3 = new Taxa(taxonList3);
    // msat datatype
    Microsatellite msat = new Microsatellite(1, 6);
    Patterns msatPatterns = new Patterns(msat, taxa3);
    // pattern in the correct code form.
    msatPatterns.addPattern(new int[] { 0, 1, 3, 2, 4, 5, 1 });
    // create tree
    NewickImporter importer = new NewickImporter("(((Taxon1:0.3,Taxon2:0.3):0.6,Taxon3:0.9):0.9,((Taxon4:0.5,Taxon5:0.5):0.3,(Taxon6:0.7,Taxon7:0.7):0.1):1.0);");
    Tree tree = importer.importTree(null);
    // treeModel
    TreeModel treeModel = new DefaultTreeModel(tree);
    // msatsubstModel
    AsymmetricQuadraticModel eu1 = new AsymmetricQuadraticModel(msat, null);
    // create msatSamplerTreeModel
    Parameter internalVal = new Parameter.Default(new double[] { 2, 3, 4, 2, 1, 5 });
    int[] externalValues = msatPatterns.getPattern(0);
    HashMap<String, Integer> taxaMap = new HashMap<String, Integer>(externalValues.length);
    boolean internalValuesProvided = true;
    for (int i = 0; i < externalValues.length; i++) {
        taxaMap.put(msatPatterns.getTaxonId(i), i);
    }
    MicrosatelliteSamplerTreeModel msatTreeModel = new MicrosatelliteSamplerTreeModel("JUnitTestEx", treeModel, internalVal, msatPatterns, externalValues, taxaMap, internalValuesProvided);
    // create msatSamplerTreeLikelihood
    BranchRateModel branchRateModel = new StrictClockBranchRates(new Parameter.Default(1.0));
    eu1Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, eu1, branchRateModel);
    // eu2
    TwoPhaseModel eu2 = new TwoPhaseModel(msat, null, eu1, new Parameter.Default(0.0), new Parameter.Default(0.4), null, false);
    eu2Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, eu2, branchRateModel);
    // ec1
    LinearBiasModel ec1 = new LinearBiasModel(msat, null, eu1, new Parameter.Default(0.48), new Parameter.Default(0.0), false, false, false);
    ec1Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, ec1, branchRateModel);
    // ec2
    TwoPhaseModel ec2 = new TwoPhaseModel(msat, null, ec1, new Parameter.Default(0.0), new Parameter.Default(0.4), null, false);
    ec2Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, ec2, branchRateModel);
    // el1
    LinearBiasModel el1 = new LinearBiasModel(msat, null, eu1, new Parameter.Default(0.2), new Parameter.Default(-0.018), true, false, false);
    el1Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, el1, branchRateModel);
    AsymmetricQuadraticModel pu1 = new AsymmetricQuadraticModel(msat, null, new Parameter.Default(1.0), new Parameter.Default(0.015), new Parameter.Default(0.0), new Parameter.Default(1.0), new Parameter.Default(0.015), new Parameter.Default(0.0), false);
    pu1Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, pu1, branchRateModel);
    // ec2
    TwoPhaseModel pu2 = new TwoPhaseModel(msat, null, pu1, new Parameter.Default(0.0), new Parameter.Default(0.4), null, false);
    pu2Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, pu2, branchRateModel);
    // ec1
    LinearBiasModel pc1 = new LinearBiasModel(msat, null, pu1, new Parameter.Default(0.48), new Parameter.Default(0.0), false, false, false);
    pc1Likelihood = new MicrosatelliteSamplerTreeLikelihood(msatTreeModel, pc1, branchRateModel);
}
Also used : HashMap(java.util.HashMap) ArrayList(java.util.ArrayList) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) Taxa(dr.evolution.util.Taxa) MicrosatelliteSamplerTreeModel(dr.evomodel.tree.MicrosatelliteSamplerTreeModel) TreeModel(dr.evomodel.tree.TreeModel) DefaultTreeModel(dr.evomodel.tree.DefaultTreeModel) NewickImporter(dr.evolution.io.NewickImporter) AsymmetricQuadraticModel(dr.oldevomodel.substmodel.AsymmetricQuadraticModel) Tree(dr.evolution.tree.Tree) Patterns(dr.evolution.alignment.Patterns) Microsatellite(dr.evolution.datatype.Microsatellite) MicrosatelliteSamplerTreeModel(dr.evomodel.tree.MicrosatelliteSamplerTreeModel) Taxon(dr.evolution.util.Taxon) DefaultTreeModel(dr.evomodel.tree.DefaultTreeModel) LinearBiasModel(dr.oldevomodel.substmodel.LinearBiasModel) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) Parameter(dr.inference.model.Parameter) MicrosatelliteSamplerTreeLikelihood(dr.oldevomodel.treelikelihood.MicrosatelliteSamplerTreeLikelihood) TwoPhaseModel(dr.oldevomodel.substmodel.TwoPhaseModel)

Aggregations

StrictClockBranchRates (dr.evomodel.branchratemodel.StrictClockBranchRates)41 BranchRateModel (dr.evomodel.branchratemodel.BranchRateModel)32 ArrayList (java.util.ArrayList)29 DefaultBranchRateModel (dr.evomodel.branchratemodel.DefaultBranchRateModel)26 Parameter (dr.inference.model.Parameter)26 TreeDataLikelihood (dr.evomodel.treedatalikelihood.TreeDataLikelihood)21 MultivariateElasticModel (dr.evomodel.continuous.MultivariateElasticModel)18 MatrixParameter (dr.inference.model.MatrixParameter)15 ArbitraryBranchRates (dr.evomodel.branchratemodel.ArbitraryBranchRates)10 NewickImporter (dr.evolution.io.NewickImporter)8 DiagonalMatrix (dr.inference.model.DiagonalMatrix)7 CladeNodeModel (dr.evomodel.bigfasttree.constrainedtree.CladeNodeModel)5 ConstrainedTreeBranchLengthProvider (dr.evomodel.bigfasttree.constrainedtree.ConstrainedTreeBranchLengthProvider)5 NodeRef (dr.evolution.tree.NodeRef)4 DefaultTreeModel (dr.evomodel.tree.DefaultTreeModel)4 TreeModel (dr.evomodel.tree.TreeModel)4 Tree (dr.evolution.tree.Tree)3 Taxon (dr.evolution.util.Taxon)3 CladeRef (dr.evomodel.bigfasttree.constrainedtree.CladeRef)3 GammaSiteRateModel (dr.evomodel.siteratemodel.GammaSiteRateModel)3