Examples with MultivariateElasticModel dr.evomodel.continuous.MultivariateElasticModel used on opensource projects

Example 1 with MultivariateElasticModel

use of dr.evomodel.continuous.MultivariateElasticModel in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    Tree treeModel = (Tree) xo.getChild(Tree.class);
    MultivariateDiffusionModel diffusionModel = (MultivariateDiffusionModel) xo.getChild(MultivariateDiffusionModel.class);
    BranchRateModel rateModel = (BranchRateModel) xo.getChild(BranchRateModel.class);
    boolean useTreeLength = xo.getAttribute(USE_TREE_LENGTH, false);
    boolean scaleByTime = xo.getAttribute(SCALE_BY_TIME, false);
    boolean reciprocalRates = xo.getAttribute(RECIPROCAL_RATES, false);
    if (reciprocalRates) {
        throw new XMLParseException("Reciprocal rates are not yet implemented.");
    }
    if (rateModel == null) {
        rateModel = new DefaultBranchRateModel();
    }
    ContinuousRateTransformation rateTransformation = new ContinuousRateTransformation.Default(treeModel, scaleByTime, useTreeLength);
    final int dim = diffusionModel.getPrecisionmatrix().length;
    String traitName = TreeTraitParserUtilities.DEFAULT_TRAIT_NAME;
    List<Integer> missingIndices;
    // Parameter sampleMissingParameter = null;
    ContinuousTraitPartialsProvider dataModel;
    boolean useMissingIndices = true;
    boolean integratedProcess = xo.getAttribute(INTEGRATED_PROCESS, false);
    if (xo.hasChildNamed(TreeTraitParserUtilities.TRAIT_PARAMETER)) {
        TreeTraitParserUtilities utilities = new TreeTraitParserUtilities();
        TreeTraitParserUtilities.TraitsAndMissingIndices returnValue = utilities.parseTraitsFromTaxonAttributes(xo, traitName, treeModel, true);
        CompoundParameter traitParameter = returnValue.traitParameter;
        missingIndices = returnValue.missingIndices;
        // sampleMissingParameter = returnValue.sampleMissingParameter;
        traitName = returnValue.traitName;
        useMissingIndices = returnValue.useMissingIndices;
        PrecisionType precisionType = PrecisionType.SCALAR;
        if (xo.getAttribute(FORCE_FULL_PRECISION, false) || (useMissingIndices && !xo.getAttribute(FORCE_COMPLETELY_MISSING, false))) {
            precisionType = PrecisionType.FULL;
        }
        if (xo.hasChildNamed(TreeTraitParserUtilities.JITTER)) {
            utilities.jitter(xo, diffusionModel.getPrecisionmatrix().length, missingIndices);
        }
        if (!integratedProcess) {
            dataModel = new ContinuousTraitDataModel(traitName, traitParameter, missingIndices, useMissingIndices, dim, precisionType);
        } else {
            dataModel = new IntegratedProcessTraitDataModel(traitName, traitParameter, missingIndices, useMissingIndices, dim, precisionType);
        }
    } else {
        // Has ContinuousTraitPartialsProvider
        dataModel = (ContinuousTraitPartialsProvider) xo.getChild(ContinuousTraitPartialsProvider.class);
    }
    ConjugateRootTraitPrior rootPrior = ConjugateRootTraitPrior.parseConjugateRootTraitPrior(xo, dataModel.getTraitDimension());
    final boolean allowSingular;
    if (dataModel instanceof IntegratedFactorAnalysisLikelihood) {
        if (traitName == TreeTraitParserUtilities.DEFAULT_TRAIT_NAME) {
            traitName = FACTOR_NAME;
        }
        if (xo.hasAttribute(ALLOW_SINGULAR)) {
            allowSingular = xo.getAttribute(ALLOW_SINGULAR, false);
        } else {
            allowSingular = true;
        }
    } else if (dataModel instanceof RepeatedMeasuresTraitDataModel) {
        traitName = ((RepeatedMeasuresTraitDataModel) dataModel).getTraitName();
        allowSingular = xo.getAttribute(ALLOW_SINGULAR, false);
    } else {
        allowSingular = xo.getAttribute(ALLOW_SINGULAR, false);
    }
    List<BranchRateModel> driftModels = AbstractMultivariateTraitLikelihood.parseDriftModels(xo, diffusionModel);
    List<BranchRateModel> optimalTraitsModels = AbstractMultivariateTraitLikelihood.parseOptimalValuesModels(xo, diffusionModel);
    MultivariateElasticModel elasticModel = null;
    if (xo.hasChildNamed(STRENGTH_OF_SELECTION_MATRIX)) {
        XMLObject cxo = xo.getChild(STRENGTH_OF_SELECTION_MATRIX);
        MatrixParameterInterface strengthOfSelectionMatrixParam;
        strengthOfSelectionMatrixParam = (MatrixParameterInterface) cxo.getChild(MatrixParameterInterface.class);
        if (strengthOfSelectionMatrixParam != null) {
            elasticModel = new MultivariateElasticModel(strengthOfSelectionMatrixParam);
        }
    }
    DiffusionProcessDelegate diffusionProcessDelegate;
    if ((optimalTraitsModels != null && elasticModel != null) || xo.getAttribute(FORCE_OU, false)) {
        if (!integratedProcess) {
            diffusionProcessDelegate = new OUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, elasticModel);
        } else {
            diffusionProcessDelegate = new IntegratedOUDiffusionModelDelegate(treeModel, diffusionModel, optimalTraitsModels, elasticModel);
        }
    } else {
        if (driftModels != null || xo.getAttribute(FORCE_DRIFT, false)) {
            diffusionProcessDelegate = new DriftDiffusionModelDelegate(treeModel, diffusionModel, driftModels);
        } else {
            diffusionProcessDelegate = new HomogeneousDiffusionModelDelegate(treeModel, diffusionModel);
        }
    }
    ContinuousDataLikelihoodDelegate delegate = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModel, rootPrior, rateTransformation, rateModel, allowSingular);
    if (dataModel instanceof IntegratedFactorAnalysisLikelihood) {
        ((IntegratedFactorAnalysisLikelihood) dataModel).setLikelihoodDelegate(delegate);
    }
    TreeDataLikelihood treeDataLikelihood = new TreeDataLikelihood(delegate, treeModel, rateModel);
    boolean reconstructTraits = xo.getAttribute(RECONSTRUCT_TRAITS, true);
    if (reconstructTraits) {
        // if (missingIndices != null && missingIndices.size() == 0) {
        if (!useMissingIndices) {
            ProcessSimulationDelegate simulationDelegate = delegate.getPrecisionType() == PrecisionType.SCALAR ? new ConditionalOnTipsRealizedDelegate(traitName, treeModel, diffusionModel, dataModel, rootPrior, rateTransformation, delegate) : new MultivariateConditionalOnTipsRealizedDelegate(traitName, treeModel, diffusionModel, dataModel, rootPrior, rateTransformation, delegate);
            TreeTraitProvider traitProvider = new ProcessSimulation(treeDataLikelihood, simulationDelegate);
            treeDataLikelihood.addTraits(traitProvider.getTreeTraits());
        } else {
            ProcessSimulationDelegate simulationDelegate = delegate.getPrecisionType() == PrecisionType.SCALAR ? new ConditionalOnTipsRealizedDelegate(traitName, treeModel, diffusionModel, dataModel, rootPrior, rateTransformation, delegate) : new MultivariateConditionalOnTipsRealizedDelegate(traitName, treeModel, diffusionModel, dataModel, rootPrior, rateTransformation, delegate);
            TreeTraitProvider traitProvider = new ProcessSimulation(treeDataLikelihood, simulationDelegate);
            treeDataLikelihood.addTraits(traitProvider.getTreeTraits());
            ProcessSimulationDelegate fullConditionalDelegate = new TipRealizedValuesViaFullConditionalDelegate(traitName, treeModel, diffusionModel, dataModel, rootPrior, rateTransformation, delegate);
            treeDataLikelihood.addTraits(new ProcessSimulation(treeDataLikelihood, fullConditionalDelegate).getTreeTraits());
        // String partialTraitName = getPartiallyMissingTraitName(traitName);
        // 
        // ProcessSimulationDelegate partialSimulationDelegate = new ProcessSimulationDelegate.ConditionalOnPartiallyMissingTipsDelegate(partialTraitName,
        // treeModel, diffusionModel, dataModel, rootPrior, rateTransformation, rateModel, delegate);
        // 
        // TreeTraitProvider partialTraitProvider = new ProcessSimulation(partialTraitName,
        // treeDataLikelihood, partialSimulationDelegate);
        // 
        // treeDataLikelihood.addTraits(partialTraitProvider.getTreeTraits());
        }
    }
    return treeDataLikelihood;
}

Example 2 with MultivariateElasticModel

use of dr.evomodel.continuous.MultivariateElasticModel in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodFullOUNonSymmetricRelaxed.

public void testLikelihoodFullOUNonSymmetricRelaxed() {
    System.out.println("\nTest Likelihood using Full Non symmetric OU Relaxed:");
    // 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.2, 100.0, 0.1 });
    strengthOfSelectionParameters[2] = new Parameter.Default(new double[] { 10.0, 0.1, 50.5 });
    MatrixParameter strengthOfSelectionMatrixParam = new MatrixParameter("strengthOfSelectionMatrix", strengthOfSelectionParameters);
    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("likelihoodFullNonSymmetricOURelaxed", 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("likelihoodFullNonSymmetricOURelaxedInf", dataLikelihoodInf);
    testConditionalMoments(dataLikelihoodInf, likelihoodDelegateInf);
}

Example 3 with MultivariateElasticModel

use of dr.evomodel.continuous.MultivariateElasticModel in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodDiagonalOU.

public void testLikelihoodDiagonalOU() {
    System.out.println("\nTest Likelihood using Diagonal OU:");
    // 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.1, 100.0, 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("likelihoodDiagonalOU", dataLikelihood);
    // Conditional moments (preorder)
    testConditionalMoments(dataLikelihood, likelihoodDelegate);
    // Conditional simulations
    MathUtils.setSeed(17890826);
    double[] expectedTraits = new double[] { -1.0, 2.0, 0.0, 1.0369622398437415, 2.065450266793184, 0.6174755164694558, 0.5, 2.0829935706195615, 5.5, 2.0, 5.0, -8.0, 11.0, 1.0, -1.5, 1.0, 2.5, 4.0 };
    testConditionalSimulations(dataLikelihood, likelihoodDelegate, diffusionModel, dataModel, rootPrior, expectedTraits);
    // Fixed Root
    ContinuousDataLikelihoodDelegate likelihoodDelegateInf = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModel, rootPriorInf, rateTransformation, rateModel, true);
    TreeDataLikelihood dataLikelihoodInf = new TreeDataLikelihood(likelihoodDelegateInf, treeModel, rateModel);
    testLikelihood("likelihoodDiagonalOUInf", dataLikelihoodInf);
    testConditionalMoments(dataLikelihoodInf, likelihoodDelegateInf);
}

Example 4 with MultivariateElasticModel

use of dr.evomodel.continuous.MultivariateElasticModel in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodDiagonalOUFactor.

public void testLikelihoodDiagonalOUFactor() {
    System.out.println("\nTest Likelihood using diagonal OU and factor:");
    // 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[] { -1.5 })));
    DiagonalMatrix strengthOfSelectionMatrixParam = new DiagonalMatrix(new Parameter.Default(new double[] { 0.5, 50.0 }));
    DiffusionProcessDelegate diffusionProcessDelegate = new OUDiffusionModelDelegate(treeModel, diffusionModelFactor, optimalTraitsModels, new MultivariateElasticModel(strengthOfSelectionMatrixParam));
    // CDL
    ContinuousDataLikelihoodDelegate likelihoodDelegateFactors = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModelFactor, rootPriorFactor, rateTransformation, rateModel, false);
    dataModelFactor.setLikelihoodDelegate(likelihoodDelegateFactors);
    // Likelihood Computation
    TreeDataLikelihood dataLikelihoodFactors = new TreeDataLikelihood(likelihoodDelegateFactors, treeModel, rateModel);
    testLikelihood("likelihoodDiagonalOUFactor", dataModelFactor, dataLikelihoodFactors);
    // Conditional simulations
    MathUtils.setSeed(17890826);
    double[] expectedTraits = new double[] { 1.3270345780274333, -1.5589839744569975, 1.241407854756886, -1.4525648723106128, 1.388017192005544, -1.533399261149814, 1.8040948421311085, -1.4189758121385794, 1.1408165195832969, -1.4607180451268982, 1.6048925583434688, -1.4333922414628846 };
    testConditionalSimulations(dataLikelihoodFactors, likelihoodDelegateFactors, diffusionModelFactor, dataModelFactor, rootPriorFactor, expectedTraits);
}

Example 5 with MultivariateElasticModel

use of dr.evomodel.continuous.MultivariateElasticModel in project beast-mcmc by beast-dev.

the class ContinuousDataLikelihoodDelegateTest method testLikelihoodFullOURelaxed.

public void testLikelihoodFullOURelaxed() {
    System.out.println("\nTest Likelihood using Full OU Relaxed:");
    // 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.2, 0.0 });
    strengthOfSelectionParameters[1] = new Parameter.Default(new double[] { 0.2, 10.5, 0.1 });
    strengthOfSelectionParameters[2] = new Parameter.Default(new double[] { 0.0, 0.1, 100.0 });
    MatrixParameter strengthOfSelectionMatrixParam = new MatrixParameter("strengthOfSelectionMatrix", strengthOfSelectionParameters);
    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("likelihoodFullOURelaxed", dataLikelihood);
    // Conditional moments (preorder)
    testConditionalMoments(dataLikelihood, likelihoodDelegate);
    // Conditional simulations
    MathUtils.setSeed(17890826);
    double[] expectedTraits = new double[] { -1.0, 2.0, 0.0, 1.6349449153945943, 2.8676718538313635, -1.0653412418514505, 0.5, 3.3661883786009166, 5.5, 2.0, 5.0, -8.0, 11.0, 1.0, -1.5, 1.0, 2.5, 4.0 };
    testConditionalSimulations(dataLikelihood, likelihoodDelegate, diffusionModel, dataModel, rootPrior, expectedTraits);
    // Fixed Root
    ContinuousDataLikelihoodDelegate likelihoodDelegateInf = new ContinuousDataLikelihoodDelegate(treeModel, diffusionProcessDelegate, dataModel, rootPriorInf, rateTransformation, rateModel, true);
    TreeDataLikelihood dataLikelihoodInf = new TreeDataLikelihood(likelihoodDelegateInf, treeModel, rateModel);
    testLikelihood("likelihoodFullOURelaxedInf", dataLikelihoodInf);
    testConditionalMoments(dataLikelihoodInf, likelihoodDelegateInf);
}