use of dr.evomodel.substmodel.EmpiricalRateMatrix in project beast-mcmc by beast-dev.
the class BeagleSeqSimTest method simulateAminoAcid.
// END: simulateThreePartitions
static void simulateAminoAcid() {
try {
System.out.println("Test case 4: simulateAminoAcid");
MathUtils.setSeed(666);
int sequenceLength = 10;
ArrayList<Partition> partitionsList = new ArrayList<Partition>();
// create tree
NewickImporter importer = new NewickImporter("(SimSeq1:73.7468,(SimSeq2:25.256989999999995,SimSeq3:45.256989999999995):18.48981);");
Tree tree = importer.importTree(null);
TreeModel treeModel = new TreeModel(tree);
// create site model
GammaSiteRateModel siteRateModel = new GammaSiteRateModel("siteModel");
// create branch rate model
BranchRateModel branchRateModel = new DefaultBranchRateModel();
// create Frequency Model
Parameter freqs = new Parameter.Default(new double[] { 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05 });
FrequencyModel freqModel = new FrequencyModel(AminoAcids.INSTANCE, freqs);
// create substitution model
EmpiricalRateMatrix rateMatrix = Blosum62.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, freqModel);
HomogeneousBranchModel substitutionModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
// create partition
Partition partition1 = new //
Partition(//
treeModel, //
substitutionModel, //
siteRateModel, //
branchRateModel, //
freqModel, // from
0, // to
sequenceLength - 1, // every
1);
partitionsList.add(partition1);
// feed to sequence simulator and generate data
BeagleSequenceSimulator simulator = new BeagleSequenceSimulator(partitionsList);
System.out.println(simulator.simulate(simulateInPar, false).toString());
} catch (Exception e) {
e.printStackTrace();
System.exit(-1);
}
// END: try-catch
}
use of dr.evomodel.substmodel.EmpiricalRateMatrix in project beast-mcmc by beast-dev.
the class PartitionData method createBranchModel.
public BranchModel createBranchModel() {
BranchModel branchModel = null;
if (this.substitutionModelIndex == 0) {
// HKY
Parameter kappa = new Parameter.Default(1, substitutionParameterValues[0]);
FrequencyModel frequencyModel = this.createFrequencyModel();
HKY hky = new HKY(kappa, frequencyModel);
branchModel = new HomogeneousBranchModel(hky);
} else if (this.substitutionModelIndex == 1) {
// GTR
Parameter ac = new Parameter.Default(1, substitutionParameterValues[1]);
Parameter ag = new Parameter.Default(1, substitutionParameterValues[2]);
Parameter at = new Parameter.Default(1, substitutionParameterValues[3]);
Parameter cg = new Parameter.Default(1, substitutionParameterValues[4]);
Parameter ct = new Parameter.Default(1, substitutionParameterValues[5]);
Parameter gt = new Parameter.Default(1, substitutionParameterValues[6]);
FrequencyModel frequencyModel = this.createFrequencyModel();
GTR gtr = new GTR(ac, ag, at, cg, ct, gt, frequencyModel);
branchModel = new HomogeneousBranchModel(gtr);
} else if (this.substitutionModelIndex == 2) {
// TN93
Parameter kappa1 = new Parameter.Default(1, substitutionParameterValues[7]);
Parameter kappa2 = new Parameter.Default(1, substitutionParameterValues[8]);
FrequencyModel frequencyModel = this.createFrequencyModel();
TN93 tn93 = new TN93(kappa1, kappa2, frequencyModel);
branchModel = new HomogeneousBranchModel(tn93);
} else if (this.substitutionModelIndex == 3) {
// Yang Codon Model
FrequencyModel frequencyModel = this.createFrequencyModel();
Parameter kappa = new Parameter.Default(1, substitutionParameterValues[9]);
Parameter omega = new Parameter.Default(1, substitutionParameterValues[10]);
GY94CodonModel yangCodonModel = new GY94CodonModel(Codons.UNIVERSAL, omega, kappa, frequencyModel);
branchModel = new HomogeneousBranchModel(yangCodonModel);
} else if (this.substitutionModelIndex == 4) {
// MG94CodonModel
FrequencyModel frequencyModel = this.createFrequencyModel();
Parameter alpha = new Parameter.Default(1, substitutionParameterValues[11]);
Parameter beta = new Parameter.Default(1, substitutionParameterValues[12]);
Parameter kappa = new Parameter.Default(1, substitutionParameterValues[13]);
MG94HKYCodonModel mg94 = new MG94HKYCodonModel(Codons.UNIVERSAL, alpha, beta, kappa, frequencyModel);
branchModel = new HomogeneousBranchModel(mg94);
} else if (this.substitutionModelIndex == 5) {
// Blosum62
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = Blosum62.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else if (this.substitutionModelIndex == 6) {
// CPREV
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = CPREV.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else if (this.substitutionModelIndex == 7) {
// Dayhoff
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = Dayhoff.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else if (this.substitutionModelIndex == 8) {
// JTT
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = JTT.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else if (this.substitutionModelIndex == 9) {
// LG
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = LG.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else if (this.substitutionModelIndex == 10) {
// MTREV
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = MTREV.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else if (this.substitutionModelIndex == 11) {
// WAG
FrequencyModel frequencyModel = this.createFrequencyModel();
EmpiricalRateMatrix rateMatrix = WAG.INSTANCE;
EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, frequencyModel);
branchModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
} else {
System.out.println("Not yet implemented");
}
return branchModel;
}
use of dr.evomodel.substmodel.EmpiricalRateMatrix in project beast-mcmc by beast-dev.
the class EmpiricalAminoAcidModelParser method parseXMLObject.
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
FrequencyModel freqModel = null;
if (xo.hasAttribute(FREQUENCIES)) {
XMLObject cxo = xo.getChild(FREQUENCIES);
freqModel = (FrequencyModel) cxo.getChild(FrequencyModel.class);
}
EmpiricalRateMatrix rateMatrix = null;
String type = xo.getStringAttribute(TYPE);
if (type.equals(AminoAcidModelType.BLOSUM_62.getXMLName())) {
rateMatrix = Blosum62.INSTANCE;
} else if (type.equals(AminoAcidModelType.DAYHOFF.getXMLName())) {
rateMatrix = Dayhoff.INSTANCE;
} else if (type.equals(AminoAcidModelType.JTT.getXMLName())) {
rateMatrix = JTT.INSTANCE;
} else if (type.equals(AminoAcidModelType.MT_REV_24.getXMLName())) {
rateMatrix = MTREV.INSTANCE;
} else if (type.equals(AminoAcidModelType.CP_REV_45.getXMLName())) {
rateMatrix = CPREV.INSTANCE;
} else if (type.equals(AminoAcidModelType.WAG.getXMLName())) {
rateMatrix = WAG.INSTANCE;
} else if (type.equals(AminoAcidModelType.LG.getXMLName())) {
rateMatrix = LG.INSTANCE;
} else if (type.equals(AminoAcidModelType.FLU.getXMLName())) {
rateMatrix = FLU.INSTANCE;
} else {
throw new XMLParseException("Unrecognized empirical amino acid model: " + type);
}
return new EmpiricalAminoAcidModel(rateMatrix, freqModel);
}
Aggregations