use of dr.evomodel.substmodel.codon.GY94CodonModel in project beast-mcmc by beast-dev.
the class RandomBranchModelParser method parseXMLObject.
@Override
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
Logger.getLogger("dr.evomodel").info("\nUsing random assignment branch model.");
TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
GY94CodonModel baseSubstitutionModel = (GY94CodonModel) xo.getElementFirstChild(BASE_MODEL);
long seed = -1;
boolean hasSeed = false;
if (xo.hasAttribute(SEED)) {
seed = xo.getLongIntegerAttribute(SEED);
hasSeed = true;
}
double rate = 1;
if (xo.hasAttribute(RATE)) {
rate = xo.getDoubleAttribute(RATE);
}
return new RandomBranchModel(treeModel, baseSubstitutionModel, rate, hasSeed, seed);
}
use of dr.evomodel.substmodel.codon.GY94CodonModel in project beast-mcmc by beast-dev.
the class CompleteHistorySimulatorTest method testCodonSimulation.
public void testCodonSimulation() {
Parameter kappa = new Parameter.Default(1, 2.0);
// Expect many more non-syn changes
Parameter omega = new Parameter.Default(1, 5.0);
Codons codons = Codons.UNIVERSAL;
int stateCount = codons.getStateCount();
double[] p = new double[stateCount];
for (int i = 0; i < stateCount; i++) {
p[i] = 1.0 / (double) stateCount;
}
Parameter pi = new Parameter.Default(p);
FrequencyModel f = new FrequencyModel(codons, pi);
GY94CodonModel codonModel = new GY94CodonModel(codons, omega, kappa, f);
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(codonModel);
BranchRateModel branchRateModel = new DefaultBranchRateModel();
double analyticResult = TreeUtils.getTreeLength(tree, tree.getRoot()) * mu.getParameterValue(0);
int nSites = 100;
// use base 61
double[] synRegMatrix = CodonLabeling.getRegisterMatrix(CodonLabeling.SYN, codons, false);
// use base 61
double[] nonSynRegMatrix = CodonLabeling.getRegisterMatrix(CodonLabeling.NON_SYN, codons, false);
runSimulation(tree, siteModel, branchRateModel, nSites, new double[][] { synRegMatrix, nonSynRegMatrix }, analyticResult);
}
use of dr.evomodel.substmodel.codon.GY94CodonModel in project beast-mcmc by beast-dev.
the class BeagleSeqSimTest method simulateRandomBranchAssignment.
// END: main
static void simulateRandomBranchAssignment() {
try {
System.out.println("Test case I dunno which: simulate random branch assignments");
MathUtils.setSeed(666);
int sequenceLength = 10;
ArrayList<Partition> partitionsList = new ArrayList<Partition>();
File treeFile = new File("/home/filip/Dropbox/BeagleSequenceSimulator/SimTree/SimTree.figtree");
Tree tree = Utils.importTreeFromFile(treeFile);
TreeModel treeModel = new TreeModel(tree);
// create Frequency Model
Parameter freqs = new Parameter.Default(Utils.UNIFORM_CODON_FREQUENCIES);
FrequencyModel freqModel = new FrequencyModel(Codons.UNIVERSAL, freqs);
// create base subst model
Parameter omegaParameter = new Parameter.Default("omega", 1, 1.0);
Parameter kappaParameter = new Parameter.Default("kappa", 1, 1.0);
GY94CodonModel baseSubModel = new GY94CodonModel(Codons.UNIVERSAL, omegaParameter, kappaParameter, freqModel);
RandomBranchModel substitutionModel = new RandomBranchModel(treeModel, baseSubModel, 0.25, false, -1);
// create site model
GammaSiteRateModel siteRateModel = new GammaSiteRateModel("siteModel");
// create branch rate model
BranchRateModel branchRateModel = new DefaultBranchRateModel();
// create partition
Partition partition1 = new //
Partition(//
treeModel, //
substitutionModel, //
siteRateModel, //
branchRateModel, //
freqModel, // from
0, // to
sequenceLength - 1, // every
1);
// Sequence ancestralSequence = new Sequence();
// ancestralSequence.appendSequenceString("TCAAGTGAGG");
// partition1.setRootSequence(ancestralSequence);
partitionsList.add(partition1);
// feed to sequence simulator and generate data
BeagleSequenceSimulator simulator = new BeagleSequenceSimulator(partitionsList);
SimpleAlignment alignment = simulator.simulate(simulateInPar, true);
// alignment.setOutputType(SimpleAlignment.OutputType.NEXUS);
alignment.setOutputType(SimpleAlignment.OutputType.XML);
System.out.println(alignment.toString());
} catch (Exception e) {
e.printStackTrace();
}
// END: try-catch
}
use of dr.evomodel.substmodel.codon.GY94CodonModel 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.codon.GY94CodonModel in project beast-mcmc by beast-dev.
the class GY94CodonModelParser method parseXMLObject.
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
Codons codons = Codons.UNIVERSAL;
if (xo.hasAttribute(GeneticCode.GENETIC_CODE)) {
String codeStr = xo.getStringAttribute(GeneticCode.GENETIC_CODE);
codons = Codons.findByName(codeStr);
}
Parameter omegaParameter = (Parameter) xo.getElementFirstChild(OMEGA);
int dim = omegaParameter.getDimension();
double value = omegaParameter.getParameterValue(dim - 1);
if (value < 0) {
throw new RuntimeException("Negative Omega parameter value " + value);
}
//END: negative check
Parameter kappaParameter = (Parameter) xo.getElementFirstChild(KAPPA);
dim = kappaParameter.getDimension();
value = kappaParameter.getParameterValue(dim - 1);
if (value < 0) {
throw new RuntimeException("Negative kappa parameter value value " + value);
}
//END: negative check
FrequencyModel freqModel = (FrequencyModel) xo.getChild(FrequencyModel.class);
GY94CodonModel codonModel = new GY94CodonModel(codons, omegaParameter, kappaParameter, freqModel);
return codonModel;
}
Aggregations