Search in sources :

Example 1 with BasicAlignment

use of jebl.evolution.alignments.BasicAlignment in project beast-mcmc by beast-dev.

the class SeqGen method simulate.

public Alignment simulate(Tree tree) {
    int[] initialSequence = new int[length];
    drawSequence(initialSequence, freqModel);
    int[] siteCategories = new int[length];
    drawSiteCategories(siteModel, siteCategories);
    double[] rates = new double[siteModel.getCategoryCount()];
    for (int i = 0; i < rates.length; i++) {
        rates[i] = siteModel.getRateForCategory(i) * substitutionRate;
    }
    for (int i = 0; i < tree.getChildCount(tree.getRoot()); i++) {
        NodeRef child = tree.getChild(tree.getRoot(), i);
        evolveSequences(initialSequence, tree, child, substModel, siteCategories, rates);
    }
    Map<State, State[]> damageMap = new HashMap<State, State[]>();
    damageMap.put(Nucleotides.A_STATE, new State[] { Nucleotides.G_STATE });
    damageMap.put(Nucleotides.C_STATE, new State[] { Nucleotides.T_STATE });
    damageMap.put(Nucleotides.G_STATE, new State[] { Nucleotides.A_STATE });
    damageMap.put(Nucleotides.T_STATE, new State[] { Nucleotides.C_STATE });
    Sequence[] sequences = new Sequence[tree.getExternalNodeCount()];
    List<NucleotideState> nucs = jebl.evolution.sequences.Nucleotides.getCanonicalStates();
    for (int i = 0; i < tree.getExternalNodeCount(); i++) {
        NodeRef node = tree.getExternalNode(i);
        int[] seq = (int[]) tree.getNodeTaxon(node).getAttribute("seq");
        State[] states = new State[seq.length];
        for (int j = 0; j < states.length; j++) {
            states[j] = nucs.get(seq[j]);
        }
        if (damageRate > 0) {
            damageSequence(states, damageRate, tree.getNodeHeight(node), damageMap);
        }
        sequences[i] = new BasicSequence(SequenceType.NUCLEOTIDE, Taxon.getTaxon(tree.getNodeTaxon(node).getId()), states);
    }
    BasicAlignment alignment = new BasicAlignment(sequences);
    return alignment;
}
Also used : HashMap(java.util.HashMap) NodeRef(dr.evolution.tree.NodeRef) BasicAlignment(jebl.evolution.alignments.BasicAlignment)

Aggregations

NodeRef (dr.evolution.tree.NodeRef)1 HashMap (java.util.HashMap)1 BasicAlignment (jebl.evolution.alignments.BasicAlignment)1