use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.
the class ComplexSubstitutionModelParser method parseXMLObject.
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
Parameter ratesParameter;
XMLObject cxo;
if (xo.hasChildNamed(FREQUENCIES)) {
cxo = xo.getChild(FREQUENCIES);
} else {
cxo = xo.getChild(ROOT_FREQUENCIES);
}
FrequencyModel freqModel = (FrequencyModel) cxo.getChild(FrequencyModel.class);
DataType dataType = freqModel.getDataType();
cxo = xo.getChild(RATES);
int states = dataType.getStateCount();
Logger.getLogger("dr.app.beagle.evomodel").info(" Complex Substitution Model (stateCount=" + states + ")");
ratesParameter = (Parameter) cxo.getChild(Parameter.class);
int rateCount = (dataType.getStateCount() - 1) * dataType.getStateCount();
if (ratesParameter == null) {
if (rateCount == 1) {
// simplest model for binary traits...
} else {
throw new XMLParseException("No rates parameter found in " + getParserName());
}
} else if (ratesParameter.getDimension() != rateCount) {
throw new XMLParseException("Rates parameter in " + getParserName() + " element should have " + rateCount + " dimensions.");
}
boolean checkConditioning = xo.getAttribute(CHECK_CONDITIONING, true);
if (!xo.hasChildNamed(INDICATOR)) {
if (!checkConditioning) {
return new ComplexSubstitutionModel(COMPLEX_SUBSTITUTION_MODEL, dataType, freqModel, ratesParameter) {
protected EigenSystem getDefaultEigenSystem(int stateCount) {
return new ComplexColtEigenSystem(stateCount, false, ColtEigenSystem.defaultMaxConditionNumber, ColtEigenSystem.defaultMaxIterations);
}
};
} else {
return new ComplexSubstitutionModel(COMPLEX_SUBSTITUTION_MODEL, dataType, freqModel, ratesParameter);
}
}
cxo = xo.getChild(INDICATOR);
Parameter indicatorParameter = (Parameter) cxo.getChild(Parameter.class);
if (indicatorParameter == null || ratesParameter == null || indicatorParameter.getDimension() != ratesParameter.getDimension())
throw new XMLParseException("Rates and indicator parameters in " + getParserName() + " element must be the same dimension.");
if (xo.hasAttribute(BSSVS_TOLERANCE)) {
double tolerance = xo.getAttribute(BSSVS_TOLERANCE, BayesianStochasticSearchVariableSelection.Utils.getTolerance());
if (tolerance > BayesianStochasticSearchVariableSelection.Utils.getTolerance()) {
// Only increase smallest allowed tolerance
BayesianStochasticSearchVariableSelection.Utils.setTolerance(tolerance);
Logger.getLogger("dr.app.beagle.evomodel").info("\tIncreasing BSSVS tolerance to " + tolerance);
}
}
if (xo.hasAttribute(BSSVS_SCALAR)) {
double scalar = xo.getAttribute(BSSVS_SCALAR, BayesianStochasticSearchVariableSelection.Utils.getScalar());
if (scalar < BayesianStochasticSearchVariableSelection.Utils.getScalar()) {
BayesianStochasticSearchVariableSelection.Utils.setScalar(scalar);
Logger.getLogger("dr.app.beagle.evomodel").info("\tDecreasing BSSVS scalar to " + scalar);
}
}
SVSComplexSubstitutionModel model;
if (!checkConditioning) {
model = new SVSComplexSubstitutionModel(SVS_COMPLEX_SUBSTITUTION_MODEL, dataType, freqModel, ratesParameter, indicatorParameter) {
protected EigenSystem getDefaultEigenSystem(int stateCount) {
return new ComplexColtEigenSystem(stateCount, false, ColtEigenSystem.defaultMaxConditionNumber, ColtEigenSystem.defaultMaxIterations);
}
};
} else {
model = new SVSComplexSubstitutionModel(SVS_COMPLEX_SUBSTITUTION_MODEL, dataType, freqModel, ratesParameter, indicatorParameter);
}
boolean randomize = xo.getAttribute(RANDOMIZE, false);
if (randomize) {
// Randomization may need multiple tries
int tries = 0;
boolean valid = false;
while (!valid && tries < maxRandomizationTries) {
BayesianStochasticSearchVariableSelection.Utils.randomize(indicatorParameter, dataType.getStateCount(), false);
valid = !Double.isInfinite(model.getLogLikelihood());
tries++;
}
Logger.getLogger("dr.app.beagle.evomodel").info("\tRandomization attempts: " + tries);
}
if (!xo.getAttribute(NORMALIZED, true)) {
model.setNormalization(false);
Logger.getLogger("dr.app.beagle.evomodel").info("\tNormalization: false");
}
Logger.getLogger("dr.app.beagle.evomodel").info("\t\tPlease cite: Edwards, Suchard et al. (2011)\n");
return model;
}
use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.
the class MarkovJumpsTreeLikelihoodParser method createTreeLikelihood.
protected BeagleTreeLikelihood createTreeLikelihood(PatternList patternList, MutableTreeModel treeModel, BranchModel branchModel, GammaSiteRateModel siteRateModel, BranchRateModel branchRateModel, TipStatesModel tipStatesModel, boolean useAmbiguities, PartialsRescalingScheme scalingScheme, boolean delayScaling, Map<Set<String>, Parameter> partialsRestrictions, XMLObject xo) throws XMLParseException {
DataType dataType = branchModel.getRootSubstitutionModel().getDataType();
String stateTag = xo.getAttribute(RECONSTRUCTION_TAG_NAME, RECONSTRUCTION_TAG);
String jumpTag = xo.getAttribute(JUMP_TAG_NAME, JUMP_TAG);
boolean scaleRewards = xo.getAttribute(SCALE_REWARDS, true);
boolean useMAP = xo.getAttribute(MAP_RECONSTRUCTION, false);
boolean useMarginalLogLikelihood = xo.getAttribute(MARGINAL_LIKELIHOOD, true);
boolean conditionalProbabilitiesInLogSpace = xo.getAttribute(CONDITIONAL_PROBABILITIES_IN_LOG_SPACE, false);
boolean useUniformization = xo.getAttribute(USE_UNIFORMIZATION, false);
boolean reportUnconditionedColumns = xo.getAttribute(REPORT_UNCONDITIONED_COLUMNS, false);
int nSimulants = xo.getAttribute(NUMBER_OF_SIMULANTS, 1);
if (patternList.areUnique()) {
throw new XMLParseException("Markov Jumps reconstruction cannot be used with compressed (unique) patterns.");
}
MarkovJumpsBeagleTreeLikelihood treeLikelihood = new MarkovJumpsBeagleTreeLikelihood(patternList, treeModel, branchModel, siteRateModel, branchRateModel, tipStatesModel, useAmbiguities, scalingScheme, delayScaling, partialsRestrictions, dataType, stateTag, useMAP, useMarginalLogLikelihood, useUniformization, reportUnconditionedColumns, nSimulants, conditionalProbabilitiesInLogSpace);
int registersFound = parseAllChildren(xo, treeLikelihood, dataType.getStateCount(), jumpTag, MarkovJumpsType.COUNTS, // For backwards compatibility
false);
XMLObject cxo = xo.getChild(COUNTS);
if (cxo != null) {
registersFound += parseAllChildren(cxo, treeLikelihood, dataType.getStateCount(), jumpTag, MarkovJumpsType.COUNTS, false);
}
cxo = xo.getChild(REWARDS);
if (cxo != null) {
registersFound += parseAllChildren(cxo, treeLikelihood, dataType.getStateCount(), jumpTag, MarkovJumpsType.REWARDS, scaleRewards);
}
if (registersFound == 0) {
// Some default values for testing
// double[] registration = new double[dataType.getStateCount()*dataType.getStateCount()];
// MarkovJumpsCore.fillRegistrationMatrix(registration,dataType.getStateCount()); // Count all transitions
// Parameter registerParameter = new Parameter.Default(registration);
// registerParameter.setId(jumpTag);
// treeLikelihood.addRegister(registerParameter,
// MarkovJumpsType.COUNTS,
// false);
// Do nothing, should run the same as AncestralStateBeagleTreeLikelihood
}
boolean saveCompleteHistory = xo.getAttribute(SAVE_HISTORY, false);
if (saveCompleteHistory) {
Parameter allCounts = new Parameter.Default(dataType.getStateCount() * dataType.getStateCount());
for (int i = 0; i < dataType.getStateCount(); ++i) {
for (int j = 0; j < dataType.getStateCount(); ++j) {
if (j == i) {
allCounts.setParameterValue(i * dataType.getStateCount() + j, 0.0);
} else {
allCounts.setParameterValue(i * dataType.getStateCount() + j, 1.0);
}
}
}
allCounts.setId(MarkovJumpsBeagleTreeLikelihood.TOTAL_COUNTS);
treeLikelihood.setLogHistories(xo.getAttribute(LOG_HISTORY, false));
treeLikelihood.setUseCompactHistory(xo.getAttribute(COMPACT_HISTORY, false));
treeLikelihood.addRegister(allCounts, MarkovJumpsType.HISTORY, false);
}
return treeLikelihood;
}
use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.
the class MarkovModulatedSubstitutionModelParser method parseXMLObject.
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
DataType dataType = DataTypeUtils.getDataType(xo);
if (!(dataType instanceof HiddenDataType)) {
throw new XMLParseException("Must construct " + MARKOV_MODULATED_MODEL + " with hidden data types. " + "You may need to provide the `-universal` extension to your hidden code type.");
}
Parameter switchingRates = (Parameter) xo.getElementFirstChild(SWITCHING_RATES);
List<SubstitutionModel> substModels = new ArrayList<SubstitutionModel>();
for (int i = 0; i < xo.getChildCount(); i++) {
Object cxo = xo.getChild(i);
if (cxo instanceof SubstitutionModel) {
substModels.add((SubstitutionModel) cxo);
}
}
boolean geometricRates = xo.getAttribute(GEOMETRIC_RATES, false);
Parameter rateScalar = xo.hasChildNamed(RATE_SCALAR) ? (Parameter) xo.getChild(RATE_SCALAR).getChild(Parameter.class) : null;
SiteRateModel siteRateModel = (SiteRateModel) xo.getChild(SiteRateModel.class);
if (siteRateModel != null) {
if (siteRateModel.getCategoryCount() != substModels.size() && substModels.size() % siteRateModel.getCategoryCount() != 0) {
throw new XMLParseException("Number of gamma categories must equal number of substitution models in " + xo.getId());
}
}
MarkovModulatedSubstitutionModel mmsm = new MarkovModulatedSubstitutionModel(xo.getId(), substModels, switchingRates, dataType, null, rateScalar, geometricRates, siteRateModel);
if (xo.getAttribute(RENORMALIZE, false)) {
mmsm.setNormalization(true);
}
return mmsm;
}
use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.
the class AlignmentParser method parseXMLObject.
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
final SimpleAlignment alignment = new SimpleAlignment();
final DataType dataType = DataTypeUtils.getDataType(xo);
if (dataType == null) {
throw new XMLParseException("dataType attribute expected for alignment element");
}
alignment.setDataType(dataType);
for (int i = 0; i < xo.getChildCount(); i++) {
final Object child = xo.getChild(i);
if (child instanceof UncertainSequence) {
alignment.addSequence((UncertainSequence) child);
} else if (child instanceof Sequence) {
alignment.addSequence((Sequence) child);
} else if (child instanceof DataType) {
// already dealt with
} else {
throw new XMLParseException("Unknown child element found in alignment");
}
}
final Logger logger = Logger.getLogger("dr.evoxml");
logger.info("\nRead alignment" + (xo.hasAttribute(XMLParser.ID) ? ": " + xo.getId() : "") + "\n Sequences = " + alignment.getSequenceCount() + "\n Sites = " + alignment.getSiteCount() + "\n Datatype = " + alignment.getDataType().getDescription());
return alignment;
}
use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.
the class MutationDeathTypeParser method parseXMLObject.
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
Logger.getLogger("dr.evolution").info("\nCreating an extended data type.");
Logger.getLogger("dr.evolution").info("\tIf you publish results using this model, please reference Alekseyenko, Lee and Suchard (in submision).\n");
DataType dataType = DataTypeUtils.getDataType(xo);
char extantChar = '\0';
XMLObject cxo = xo.getChild(EXTANT);
if (cxo != null) {
extantChar = cxo.getStringAttribute(CODE).charAt(0);
}
cxo = xo.getChild(STATE);
char stateChar;
if (cxo != null) {
stateChar = cxo.getStringAttribute(CODE).charAt(0);
} else {
stateChar = dataType.getChar(dataType.getGapState());
}
Logger.getLogger("dr.evolution").info("\tNon-existent code: " + stateChar);
if (dataType == null && extantChar == '\0')
throw new XMLParseException("In " + xo.getName() + " you must either provide a data type or a code for extant state");
MutationDeathType mdt;
if (dataType != null) {
Logger.getLogger("dr.evolution").info("\tBase type: " + dataType.toString());
mdt = new MutationDeathType(dataType, stateChar);
} else {
Logger.getLogger("dr.evolution").info("\tExtant code: " + extantChar);
mdt = new MutationDeathType(stateChar, extantChar);
}
char ambiguityChar = '\0';
String states;
Object dxo;
for (int i = 0; i < xo.getChildCount(); ++i) {
dxo = xo.getChild(i);
if (dxo instanceof XMLObject) {
cxo = (XMLObject) dxo;
if (cxo.getName().equals(AMBIGUITY)) {
ambiguityChar = cxo.getStringAttribute(CODE).charAt(0);
if (cxo.hasAttribute(STATES)) {
states = cxo.getStringAttribute(STATES);
} else {
states = "";
}
mdt.addAmbiguity(ambiguityChar, states);
Logger.getLogger("dr.evolution").info("\tAmbiguity code: " + ambiguityChar);
}
}
}
return mdt;
}
Aggregations