Examples with DataType dr.evolution.datatype.DataType used on opensource projects

Example 1 with DataType

use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.

the class AncestralStateTreeLikelihoodParser method createTreeLikelihood.

protected BeagleTreeLikelihood createTreeLikelihood(//
PatternList patternList, //
TreeModel treeModel, //
BranchModel branchModel, //
GammaSiteRateModel siteRateModel, //
BranchRateModel branchRateModel, //
TipStatesModel tipStatesModel, //
boolean useAmbiguities, //
PartialsRescalingScheme scalingScheme, boolean delayScaling, Map<//
Set<String>, Parameter> partialsRestrictions, //
XMLObject xo) throws XMLParseException {
    //		System.err.println("XML object: " + xo.toString());
    DataType dataType = branchModel.getRootSubstitutionModel().getDataType();
    // default tag is RECONSTRUCTION_TAG
    String tag = xo.getAttribute(RECONSTRUCTION_TAG_NAME, RECONSTRUCTION_TAG);
    boolean useMAP = xo.getAttribute(MAP_RECONSTRUCTION, false);
    boolean useMarginalLogLikelihood = xo.getAttribute(MARGINAL_LIKELIHOOD, true);
    if (patternList.areUnique()) {
        throw new XMLParseException("Ancestral state reconstruction cannot be used with compressed (unique) patterns.");
    }
    return new // Current just returns a OldBeagleTreeLikelihood
    AncestralStateBeagleTreeLikelihood(patternList, treeModel, branchModel, siteRateModel, branchRateModel, tipStatesModel, useAmbiguities, scalingScheme, delayScaling, partialsRestrictions, dataType, tag, useMAP, useMarginalLogLikelihood);
}

Example 2 with DataType

use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.

the class DiscreteTraitBranchRateModelParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    PatternList patternList = (PatternList) xo.getChild(PatternList.class);
    TreeTraitProvider traitProvider = (TreeTraitProvider) xo.getChild(TreeTraitProvider.class);
    DataType dataType = DataTypeUtils.getDataType(xo);
    Parameter rateParameter = null;
    Parameter relativeRatesParameter = null;
    Parameter indicatorsParameter = null;
    if (xo.getChild(RATE) != null) {
        rateParameter = (Parameter) xo.getElementFirstChild(RATE);
    }
    if (xo.getChild(RATES) != null) {
        rateParameter = (Parameter) xo.getElementFirstChild(RATES);
    }
    if (xo.getChild(RELATIVE_RATES) != null) {
        relativeRatesParameter = (Parameter) xo.getElementFirstChild(RELATIVE_RATES);
    }
    if (xo.getChild(INDICATORS) != null) {
        indicatorsParameter = (Parameter) xo.getElementFirstChild(INDICATORS);
    }
    int traitIndex = xo.getAttribute(TRAIT_INDEX, 1) - 1;
    String traitName = "states";
    Logger.getLogger("dr.evomodel").info("Using discrete trait branch rate model.\n" + "\tIf you use this model, please cite:\n" + "\t\tDrummond and Suchard (in preparation)");
    if (traitProvider == null) {
        // Use the version that reconstructs the trait using parsimony:
        return new DiscreteTraitBranchRateModel(treeModel, patternList, traitIndex, rateParameter);
    } else {
        if (traitName != null) {
            TreeTrait trait = traitProvider.getTreeTrait(traitName);
            if (trait == null) {
                throw new XMLParseException("A trait called, " + traitName + ", was not available from the TreeTraitProvider supplied to " + getParserName() + ", with ID " + xo.getId());
            }
            if (relativeRatesParameter != null) {
                return new DiscreteTraitBranchRateModel(traitProvider, dataType, treeModel, trait, traitIndex, rateParameter, relativeRatesParameter, indicatorsParameter);
            } else {
                return new DiscreteTraitBranchRateModel(traitProvider, dataType, treeModel, trait, traitIndex, rateParameter);
            }
        } else {
            TreeTrait[] traits = new TreeTrait[dataType.getStateCount()];
            for (int i = 0; i < dataType.getStateCount(); i++) {
                traits[i] = traitProvider.getTreeTrait(dataType.getCode(i));
                if (traits[i] == null) {
                    throw new XMLParseException("A trait called, " + dataType.getCode(i) + ", was not available from the TreeTraitProvider supplied to " + getParserName() + ", with ID " + xo.getId());
                }
            }
            return new DiscreteTraitBranchRateModel(traitProvider, traits, treeModel, rateParameter);
        }
    }
}

Example 3 with DataType

use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.

the class GLMSubstitutionModelParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    DataType dataType = DataTypeUtils.getDataType(xo);
    if (dataType == null)
        dataType = (DataType) xo.getChild(DataType.class);
    int rateCount = (dataType.getStateCount() - 1) * dataType.getStateCount();
    // Should be constructed as a log-linear model
    GeneralizedLinearModel glm = (GeneralizedLinearModel) xo.getChild(GeneralizedLinearModel.class);
    // LogLinearModel glm = (LogLinearModel) xo.getChild(GeneralizedLinearModel.class);
    int length = glm.getXBeta().length;
    if (length != rateCount) {
        throw new XMLParseException("Rates parameter in " + getParserName() + " element should have " + (rateCount) + " dimensions.  However GLM dimension is " + length);
    }
    XMLObject cxo = xo.getChild(ComplexSubstitutionModelParser.ROOT_FREQUENCIES);
    FrequencyModel rootFreq = (FrequencyModel) cxo.getChild(FrequencyModel.class);
    if (dataType != rootFreq.getDataType()) {
        throw new XMLParseException("Data type of " + getParserName() + " element does not match that of its rootFrequencyModel.");
    }
    return new GLMSubstitutionModel(xo.getId(), dataType, rootFreq, glm);
}

Example 4 with DataType

use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.

the class GeneralSubstitutionModelParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    Parameter ratesParameter = null;
    FrequencyModel freqModel = null;
    if (xo.hasChildNamed(FREQUENCIES)) {
        XMLObject cxo = xo.getChild(FREQUENCIES);
        freqModel = (FrequencyModel) cxo.getChild(FrequencyModel.class);
    }
    DataType dataType = DataTypeUtils.getDataType(xo);
    if (dataType == null)
        dataType = (DataType) xo.getChild(DataType.class);
    if (dataType == null)
        dataType = freqModel.getDataType();
    if (dataType != freqModel.getDataType()) {
        throw new XMLParseException("Data type of " + getParserName() + " element does not match that of its frequencyModel.");
    }
    XMLObject cxo = xo.getChild(RATES);
    ratesParameter = (Parameter) cxo.getChild(Parameter.class);
    int states = dataType.getStateCount();
    Logger.getLogger("dr.evomodel").info("  General Substitution Model (stateCount=" + states + ")");
    boolean hasRelativeRates = cxo.hasChildNamed(RELATIVE_TO) || (cxo.hasAttribute(RELATIVE_TO) && cxo.getIntegerAttribute(RELATIVE_TO) > 0);
    int nonReversibleRateCount = ((dataType.getStateCount() - 1) * dataType.getStateCount());
    int reversibleRateCount = (nonReversibleRateCount / 2);
    boolean isNonReversible = ratesParameter.getDimension() == nonReversibleRateCount;
    boolean hasIndicator = xo.hasChildNamed(INDICATOR);
    if (!hasRelativeRates) {
        Parameter indicatorParameter = null;
        if (ratesParameter.getDimension() != reversibleRateCount && ratesParameter.getDimension() != nonReversibleRateCount) {
            throw new XMLParseException("Rates parameter in " + getParserName() + " element should have " + (reversibleRateCount) + " dimensions for reversible model or " + nonReversibleRateCount + " dimensions for non-reversible. " + "However parameter dimension is " + ratesParameter.getDimension());
        }
        if (hasIndicator) {
            // this is using BSSVS
            cxo = xo.getChild(INDICATOR);
            indicatorParameter = (Parameter) cxo.getChild(Parameter.class);
            if (indicatorParameter.getDimension() != ratesParameter.getDimension()) {
                throw new XMLParseException("Rates and indicator parameters in " + getParserName() + " element must be the same dimension.");
            }
            boolean randomize = xo.getAttribute(ComplexSubstitutionModelParser.RANDOMIZE, false);
            if (randomize) {
                BayesianStochasticSearchVariableSelection.Utils.randomize(indicatorParameter, dataType.getStateCount(), !isNonReversible);
            }
        }
        if (isNonReversible) {
            // if (xo.hasChildNamed(ROOT_FREQ)) {
            // cxo = xo.getChild(ROOT_FREQ);
            // FrequencyModel rootFreq = (FrequencyModel) cxo.getChild(FrequencyModel.class);
            // 
            // if (dataType != rootFreq.getDataType()) {
            // throw new XMLParseException("Data type of " + getParserName() + " element does not match that of its rootFrequencyModel.");
            // }
            // 
            // Logger.getLogger("dr.evomodel").info("  Using BSSVS Complex Substitution Model");
            // return new SVSComplexSubstitutionModel(getParserName(), dataType, freqModel, ratesParameter, indicatorParameter);
            // 
            // } else {
            // throw new XMLParseException("Non-reversible model missing " + ROOT_FREQ + " element");
            // }
            Logger.getLogger("dr.evomodel").info("  Using BSSVS Complex Substitution Model");
            return new SVSComplexSubstitutionModel(getParserName(), dataType, freqModel, ratesParameter, indicatorParameter);
        } else {
            Logger.getLogger("dr.evomodel").info("  Using BSSVS General Substitution Model");
            return new SVSGeneralSubstitutionModel(getParserName(), dataType, freqModel, ratesParameter, indicatorParameter);
        }
    } else {
        if (ratesParameter.getDimension() != reversibleRateCount - 1) {
            throw new XMLParseException("Rates parameter in " + getParserName() + " element should have " + (reversibleRateCount - 1) + " dimensions. However parameter dimension is " + ratesParameter.getDimension());
        }
        int relativeTo = 0;
        if (hasRelativeRates) {
            relativeTo = cxo.getIntegerAttribute(RELATIVE_TO) - 1;
        }
        if (relativeTo < 0 || relativeTo >= reversibleRateCount) {
            throw new XMLParseException(RELATIVE_TO + " must be 1 or greater");
        } else {
            int t = relativeTo;
            int s = states - 1;
            int row = 0;
            while (t >= s) {
                t -= s;
                s -= 1;
                row += 1;
            }
            int col = t + row + 1;
            Logger.getLogger("dr.evomodel").info("  Rates relative to " + dataType.getCode(row) + "<->" + dataType.getCode(col));
        }
        if (ratesParameter == null) {
            if (reversibleRateCount == 1) {
            // simplest model for binary traits...
            } else {
                throw new XMLParseException("No rates parameter found in " + getParserName());
            }
        }
        return new GeneralSubstitutionModel(getParserName(), dataType, freqModel, ratesParameter, relativeTo);
    }
}

Example 5 with DataType

use of dr.evolution.datatype.DataType in project beast-mcmc by beast-dev.

the class FrequencyModelParser method getEmpirical3x4Freqs.

// END: parseXMLObject
private double[] getEmpirical3x4Freqs(PatternList patternList) {
    DataType nucleotideDataType = Nucleotides.INSTANCE;
    Codons codonDataType = Codons.UNIVERSAL;
    List<String> stopCodonsList = Arrays.asList(STOP_CODONS);
    int cStateCount = codonDataType.getStateCount();
    int nStateCount = nucleotideDataType.getStateCount();
    int nPosition = 3;
    double[] stopCodonFreqs = new double[STOP_CODONS.length];
    double[][] counts = new double[nStateCount][nPosition];
    int[] countsPos = new int[nPosition];
    int patternCount = patternList.getPatternCount();
    for (int i = 0; i < patternCount; i++) {
        int[] sitePatterns = patternList.getPattern(i);
        for (int codonState : sitePatterns) {
            int[] nucleotideStates = codonDataType.getTripletStates(codonState);
            String triplet = codonDataType.getTriplet(codonState);
            if (stopCodonsList.contains(triplet)) {
                int stopCodonIndex = stopCodonsList.indexOf(triplet);
                stopCodonFreqs[stopCodonIndex]++;
            }
            for (int pos = 0; pos < nPosition; pos++) {
                int nucleotideState = nucleotideStates[pos];
                counts[nucleotideState][pos]++;
                countsPos[pos]++;
            }
        // END: nucleotide positions loop
        }
    // END: sitePatterns loop
    }
    // sites loop
    int total = 0;
    for (int pos = 0; pos < nPosition; pos++) {
        int totalPos = countsPos[pos];
        for (int s = 0; s < nStateCount; s++) {
            counts[s][pos] = counts[s][pos] / totalPos;
        }
        // END: nucleotide states loop
        total += totalPos;
    }
    // END: nucleotide positions loop
    // Utils.printArray(stopCodonFreqs);
    // System.out.println(stopCodonFreqs.length);
    // add stop codon frequencies
    double pi = 0.0;
    for (double stopCodonFreq : stopCodonFreqs) {
        double freq = stopCodonFreq / total;
        pi += freq;
    }
    // System.out.println(pi);
    double[] freqs = new double[cStateCount];
    Arrays.fill(freqs, 1.0);
    for (int codonState = 0; codonState < cStateCount; codonState++) {
        int[] nucleotideStates = codonDataType.getTripletStates(codonState);
        for (int pos = 0; pos < nPosition; pos++) {
            int nucleotide = nucleotideStates[pos];
            freqs[codonState] *= counts[nucleotide][pos];
        }
        // END: nucleotide positions loop
        // TODO: stop codons freqs
        freqs[codonState] = freqs[codonState] / (1 - pi);
    }
    return freqs;
}