Examples with FrequencyModel dr.evomodel.substmodel.FrequencyModel used on opensource projects

Example 21 with FrequencyModel

use of dr.evomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class ProductChainFrequencyModelTest method testFrequencyModel.

public void testFrequencyModel() {
    FrequencyModel firstPosition = new FrequencyModel(Nucleotides.INSTANCE, freq1);
    FrequencyModel secondPosition = new FrequencyModel(Nucleotides.INSTANCE, freq2);
    FrequencyModel thirdPosition = new FrequencyModel(Nucleotides.INSTANCE, freq3);
    List<FrequencyModel> freqModels = new ArrayList<FrequencyModel>(3);
    freqModels.add(firstPosition);
    freqModels.add(secondPosition);
    freqModels.add(thirdPosition);
    ProductChainFrequencyModel pcFreqModel = new ProductChainFrequencyModel("freq", freqModels);
    double[] freqs = pcFreqModel.getFrequencies();
    System.out.println("Freq.length = " + freqs.length);
    int pos1 = 2;
    int pos2 = 1;
    int pos3 = 3;
    int index = computeIndex(pos1, pos2, pos3);
    System.out.println("Entry: " + new Vector(pcFreqModel.decomposeEntry(index)));
    System.out.println("Freq = " + freqs[index]);
    System.out.println("Freq = " + computeFreq(pos1, pos2, pos3));
    assertEquals(computeFreq(pos1, pos2, pos3), freqs[index]);
}

Example 22 with FrequencyModel

use of dr.evomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class MarkovJumpsSubstitutionModelTest method testMarginalRates.

public void testMarginalRates() {
    HKY substModel = new HKY(2.0, new FrequencyModel(Nucleotides.INSTANCE, // A,C,G,T
    new double[] { 0.3, 0.2, 0.25, 0.25 }));
    int states = substModel.getDataType().getStateCount();
    MarkovJumpsSubstitutionModel markovjumps = new MarkovJumpsSubstitutionModel(substModel, MarkovJumpsType.COUNTS);
    double[] r = new double[states * states];
    // A
    int from = 0;
    // C
    int to = 1;
    MarkovJumpsCore.fillRegistrationMatrix(r, from, to, states, 1.0);
    markovjumps.setRegistration(r);
    double marginalRate = markovjumps.getMarginalRate();
    System.out.println("Marginal rate = " + marginalRate);
    assertEquals(rMarkovMarginalRate, marginalRate, tolerance);
    MarkovJumpsCore.fillRegistrationMatrix(r, states);
    markovjumps.setRegistration(r);
    marginalRate = markovjumps.getMarginalRate();
    System.out.println("Marginal rate = " + marginalRate);
    assertEquals(1.0, marginalRate, tolerance);
}

Example 23 with FrequencyModel

use of dr.evomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class CompleteHistorySimulatorTest method testHKYSimulation.

public void testHKYSimulation() {
    Parameter kappa = new Parameter.Default(1, 2.0);
    double[] pi = { 0.45, 0.05, 0.25, 0.25 };
    Parameter freqs = new Parameter.Default(pi);
    FrequencyModel f = new FrequencyModel(Nucleotides.INSTANCE, freqs);
    HKY hky = new HKY(kappa, f);
    int stateCount = hky.getDataType().getStateCount();
    Parameter mu = new Parameter.Default(1, 0.5);
    Parameter alpha = new Parameter.Default(1, 0.5);
    GammaSiteRateModel siteModel = new GammaSiteRateModel("gammaModel", mu, alpha, 4, null);
    siteModel.setSubstitutionModel(hky);
    BranchRateModel branchRateModel = new DefaultBranchRateModel();
    double analyticResult = TreeUtils.getTreeLength(tree, tree.getRoot()) * mu.getParameterValue(0);
    int nSites = 200;
    double[] register1 = new double[stateCount * stateCount];
    double[] register2 = new double[stateCount * stateCount];
    // Count all jumps
    MarkovJumpsCore.fillRegistrationMatrix(register1, stateCount);
    // Move some jumps from 1 to 2
    register1[1 * stateCount + 2] = 0;
    register2[1 * stateCount + 2] = 1;
    register1[1 * stateCount + 3] = 0;
    register2[1 * stateCount + 3] = 1;
    register1[2 * stateCount + 3] = 0;
    register2[2 * stateCount + 3] = 1;
    runSimulation(tree, siteModel, branchRateModel, nSites, new double[][] { register1, register2 }, analyticResult);
}

Example 24 with FrequencyModel

use of dr.evomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class TwoStateMarkovRewardsTest method testTwoStateRewards.

public void testTwoStateRewards() {
    DataType dataType = TwoStates.INSTANCE;
    FrequencyModel freqModel = new FrequencyModel(TwoStates.INSTANCE, new double[] { 0.5, 0.5 });
    Parameter rates = new Parameter.Default(new double[] { 4.0, 6.0 });
    ComplexSubstitutionModel twoStateModel = new ComplexSubstitutionModel("two", dataType, freqModel, rates) {
    };
    twoStateModel.setNormalization(false);
    MarkovJumpsSubstitutionModel markovRewards = new MarkovJumpsSubstitutionModel(twoStateModel, MarkovJumpsType.REWARDS);
    double[] r = new double[2];
    double[] q = new double[4];
    double[] c = new double[4];
    int mark = 0;
    double weight = 1.0;
    r[mark] = weight;
    markovRewards.setRegistration(r);
    twoStateModel.getInfinitesimalMatrix(q);
    System.out.println("Q = " + new Vector(q));
    System.out.println("Reward for state 0");
    double time = 1.0;
    markovRewards.computeCondStatMarkovJumps(time, c);
    System.out.println("Reward conditional on X(0) = i, X(t) = j: " + new Vector(c));
    double endTime = 10.0;
    int steps = 10;
    for (time = 0.0; time < endTime; time += (endTime / steps)) {
        markovRewards.computeCondStatMarkovJumps(time, c);
        // start = 0, end = 0
        System.out.println(time + "," + c[0]);
    }
}

Example 25 with FrequencyModel

use of dr.evomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class LineageSpecificBranchModelParser method parseXMLObject.

@Override
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    FrequencyModel rootFrequencyModel = (FrequencyModel) xo.getChild(FrequencyModel.class);
    List<SubstitutionModel> substitutionModels = new ArrayList<SubstitutionModel>();
    XMLObject cxo = xo.getChild(MODELS);
    for (int i = 0; i < cxo.getChildCount(); i++) {
        SubstitutionModel substModel = (SubstitutionModel) cxo.getChild(i);
        substitutionModels.add(substModel);
    }
    //END: models loop
    // TODO: check if categories numbering starts from zero
    Parameter categories = (Parameter) xo.getElementFirstChild(CATEGORIES);
    return new //provider, 
    LineageSpecificBranchModel(//provider, 
    treeModel, //provider, 
    rootFrequencyModel, //provider, 
    substitutionModels, categories);
}