Examples with TwoPhaseModel dr.oldevomodel.substmodel.TwoPhaseModel used on opensource projects

Example 1 with TwoPhaseModel

use of dr.oldevomodel.substmodel.TwoPhaseModel 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 TreeModel(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);
}

Example 2 with TwoPhaseModel

use of dr.oldevomodel.substmodel.TwoPhaseModel in project beast-mcmc by beast-dev.

the class TwoPhaseModelParser method parseXMLObject.

//AbstractXMLObjectParser implementation
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    OnePhaseModel subModel = (OnePhaseModel) xo.getElementFirstChild(SUBMODEL);
    Microsatellite dataType = (Microsatellite) xo.getChild(Microsatellite.class);
    Parameter.Default geoParam = (Parameter.Default) xo.getElementFirstChild(GEO_PARAM);
    Parameter paramP = (Parameter) xo.getElementFirstChild(ONEPHASEPR_PARAM);
    Parameter limitE = null;
    if (xo.hasChildNamed(TRANS_PARAM)) {
        limitE = (Parameter) xo.getElementFirstChild(TRANS_PARAM);
    }
    boolean estimateSubmodelParams = xo.getAttribute(ESTIMATE_SUBMODEL_PARAMS, false);
    FrequencyModel freqModel = null;
    if (xo.hasChildNamed(FrequencyModelParser.FREQUENCIES)) {
        freqModel = (FrequencyModel) xo.getElementFirstChild(FrequencyModelParser.FREQUENCIES);
    }
    return new TwoPhaseModel(dataType, freqModel, subModel, paramP, geoParam, limitE, estimateSubmodelParams);
}

Example 3 with TwoPhaseModel

use of dr.oldevomodel.substmodel.TwoPhaseModel in project beast-mcmc by beast-dev.

the class TwoPhaseModelTest method testTwoPhaseModel.

public void testTwoPhaseModel() {
    for (Instance test : all) {
        OnePhaseModel subModel = test.getSubModel();
        Microsatellite microsat = (Microsatellite) subModel.getDataType();
        Parameter pParam = new Parameter.Default(test.getPParam());
        Parameter mParam = new Parameter.Default(test.getMParam());
        TwoPhaseModel tpm = new TwoPhaseModel(microsat, null, subModel, pParam, mParam, null, false);
        int k;
        tpm.computeStationaryDistribution();
        double[] statDist = tpm.getStationaryDistribution();
        final double[] expectedStatDist = test.getPi();
        for (k = 0; k < statDist.length; ++k) {
            assertEquals(statDist[k], expectedStatDist[k], 1e-10);
        }
        int stateCount = microsat.getStateCount();
        double[] mat = new double[stateCount * stateCount];
        tpm.getTransitionProbabilities(test.getDistance(), mat);
        final double[] result = test.getExpectedResult();
        for (k = 0; k < mat.length; ++k) {
            assertEquals(result[k], mat[k], 5e-9);
        //System.out.print(" " + (mat[k]));// - result[k]));
        }
        k = 0;
        for (int i = 0; i < microsat.getStateCount(); i++) {
            for (int j = 0; j < microsat.getStateCount(); j++) {
                assertEquals(result[k++], tpm.getOneTransitionProbabilityEntry(test.getDistance(), i, j), 1e-10);
            }
        }
        for (int j = 0; j < microsat.getStateCount(); j++) {
            double[] colTransitionProb = tpm.getColTransitionProbabilities(test.getDistance(), j);
            for (int i = 0; i < microsat.getStateCount(); i++) {
                assertEquals(result[i * microsat.getStateCount() + j], colTransitionProb[i], 1e-10);
            }
        }
        for (int i = 0; i < microsat.getStateCount(); i++) {
            double[] rowTransitionProb = tpm.getRowTransitionProbabilities(test.getDistance(), i);
            for (int j = 0; j < microsat.getStateCount(); j++) {
                assertEquals(result[i * microsat.getStateCount() + j], rowTransitionProb[j], 1e-10);
            }
        }
    }
}