Search in sources :

Example 31 with FrequencyModel

use of dr.oldevomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class EmpiricalCodonModelParser method createNewFreqModel.

// creates new FrequencyModel from XML frequencies
private FrequencyModel createNewFreqModel(DataType codons, EmpiricalCodonRateMatrix type) throws XMLParseException {
    double[] freqs = type.getFrequencies();
    double sum = 0;
    for (int j = 0; j < freqs.length; j++) {
        sum += freqs[j];
    }
    if (Math.abs(sum - 1.0) > 1e-8) {
        throw new XMLParseException("Frequencies do not sum to 1 (they sum to " + sum + ")");
    }
    FrequencyModel fm = new FrequencyModel(codons, freqs);
    Logger.getLogger("dr.evomodel").info("Using frequencies from data file");
    return fm;
}
Also used : FrequencyModel(dr.oldevomodel.substmodel.FrequencyModel)

Example 32 with FrequencyModel

use of dr.oldevomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class SeqGen method main.

public static void main(String[] argv) {
    String treeFileName = argv[0];
    String outputFileStem = argv[1];
    int length = 500;
    double[] frequencies = new double[] { 0.25, 0.25, 0.25, 0.25 };
    double kappa = 10.0;
    double alpha = 0.5;
    double substitutionRate = argv.length < 3 ? 1.0E-3 : Double.parseDouble(argv[2]);
    int categoryCount = argv.length < 4 ? 8 : Integer.parseInt(argv[3]);
    //1.56E-6;
    double damageRate = argv.length < 5 ? 0 : Double.parseDouble(argv[4]);
    System.out.println("substitutionRate = " + substitutionRate + "; categoryCount = " + categoryCount + "; damageRate = " + damageRate);
    FrequencyModel freqModel = new FrequencyModel(dr.evolution.datatype.Nucleotides.INSTANCE, frequencies);
    HKY hkyModel = new HKY(kappa, freqModel);
    SiteModel siteModel = null;
    if (categoryCount > 1) {
        siteModel = new GammaSiteModel(hkyModel, alpha, categoryCount);
    } else {
        // no rate heterogeneity
        siteModel = new GammaSiteModel(hkyModel);
    }
    List<Tree> trees = new ArrayList<Tree>();
    FileReader reader = null;
    try {
        reader = new FileReader(treeFileName);
        //            TreeImporter importer = new NexusImporter(reader);
        TreeImporter importer = new NewickImporter(reader);
        while (importer.hasTree()) {
            Tree tree = importer.importNextTree();
            trees.add(tree);
            System.out.println("tree height = " + tree.getNodeHeight(tree.getRoot()) + "; leave nodes = " + tree.getExternalNodeCount());
        }
    } catch (FileNotFoundException e) {
        e.printStackTrace();
        return;
    } catch (Importer.ImportException e) {
        e.printStackTrace();
        return;
    } catch (IOException e) {
        e.printStackTrace();
        return;
    }
    SeqGen seqGen = new SeqGen(length, substitutionRate, freqModel, hkyModel, siteModel, damageRate);
    int i = 1;
    for (Tree tree : trees) {
        Alignment alignment = seqGen.simulate(tree);
        FileWriter writer = null;
        try {
            //                writer = new FileWriter(outputFileStem + (i < 10 ? "00" : (i < 100 ? "0" : "")) + i + ".nex");
            //                NexusExporter exporter = new NexusExporter(writer);
            //
            //                exporter.exportAlignment(alignment);
            //
            //                writer.close();
            String outputFileName = outputFileStem + "-" + substitutionRate + ".fasta";
            writer = new FileWriter(outputFileName);
            BufferedWriter bf = new BufferedWriter(writer);
            FastaExporter exporter = new FastaExporter(bf);
            exporter.exportSequences(alignment.getSequenceList());
            bf.close();
            System.out.println("Write " + i + "th sequence file : " + outputFileName);
            i++;
        } catch (IOException e) {
            e.printStackTrace();
            return;
        }
    }
}
Also used : FrequencyModel(dr.oldevomodel.substmodel.FrequencyModel) FastaExporter(jebl.evolution.io.FastaExporter) ArrayList(java.util.ArrayList) GammaSiteModel(dr.oldevomodel.sitemodel.GammaSiteModel) SiteModel(dr.oldevomodel.sitemodel.SiteModel) Alignment(jebl.evolution.alignments.Alignment) BasicAlignment(jebl.evolution.alignments.BasicAlignment) GammaSiteModel(dr.oldevomodel.sitemodel.GammaSiteModel) HKY(dr.oldevomodel.substmodel.HKY) NewickImporter(dr.evolution.io.NewickImporter) TreeImporter(dr.evolution.io.TreeImporter) Tree(dr.evolution.tree.Tree) NewickImporter(dr.evolution.io.NewickImporter) Importer(dr.evolution.io.Importer) TreeImporter(dr.evolution.io.TreeImporter)

Example 33 with FrequencyModel

use of dr.oldevomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class SequenceSimulator method simulate.

//END: sequence2intArray
/**
	 * perform the actual sequence generation
	 * @return alignment containing randomly generated sequences for the nodes in the
	 * leaves of the tree
	 */
public Alignment simulate() {
    NodeRef root = m_tree.getRoot();
    double[] categoryProbs = m_siteModel.getCategoryProportions();
    int[] category = new int[m_sequenceLength];
    for (int i = 0; i < m_sequenceLength; i++) {
        category[i] = MathUtils.randomChoicePDF(categoryProbs);
    }
    int[] seq = new int[m_sequenceLength];
    if (has_ancestralSequence) {
        seq = sequence2intArray(ancestralSequence);
    } else {
        FrequencyModel frequencyModel = m_siteModel.getFrequencyModel();
        for (int i = 0; i < m_sequenceLength; i++) {
            seq[i] = MathUtils.randomChoicePDF(frequencyModel.getFrequencies());
        }
    }
    if (DEBUG) {
        synchronized (this) {
            System.out.println();
            System.out.println("root Sequence:");
            Utils.printArray(seq);
        }
    }
    SimpleAlignment alignment = new SimpleAlignment();
    alignment.setReportCountStatistics(false);
    alignment.setDataType(m_siteModel.getFrequencyModel().getDataType());
    traverse(root, seq, category, alignment);
    return alignment;
}
Also used : FrequencyModel(dr.oldevomodel.substmodel.FrequencyModel) NodeRef(dr.evolution.tree.NodeRef) SimpleAlignment(dr.evolution.alignment.SimpleAlignment)

Example 34 with FrequencyModel

use of dr.oldevomodel.substmodel.FrequencyModel 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();
    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);
}
Also used : FrequencyModel(dr.oldevomodel.substmodel.FrequencyModel) GeneralizedLinearModel(dr.inference.distribution.GeneralizedLinearModel) DataType(dr.evolution.datatype.DataType) LogLinearModel(dr.inference.distribution.LogLinearModel) GLMSubstitutionModel(dr.oldevomodel.substmodel.GLMSubstitutionModel)

Example 35 with FrequencyModel

use of dr.oldevomodel.substmodel.FrequencyModel in project beast-mcmc by beast-dev.

the class GeneralF81ModelParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    FrequencyModel freqModel = null;
    if (xo.hasChildNamed(FREQUENCIES)) {
        XMLObject cxo = xo.getChild(FREQUENCIES);
        freqModel = (FrequencyModel) cxo.getChild(FrequencyModel.class);
    }
    Logger.getLogger("dr.evomodel").info("  General F81 Model from frequencyModel '" + freqModel.getId() + "' (stateCount=" + freqModel.getFrequencyCount() + ")");
    return new GeneralF81Model(freqModel);
}
Also used : FrequencyModel(dr.oldevomodel.substmodel.FrequencyModel) GeneralF81Model(dr.oldevomodel.substmodel.GeneralF81Model)

Aggregations

FrequencyModel (dr.oldevomodel.substmodel.FrequencyModel)59 Parameter (dr.inference.model.Parameter)44 GammaSiteModel (dr.oldevomodel.sitemodel.GammaSiteModel)20 HKY (dr.oldevomodel.substmodel.HKY)18 SitePatterns (dr.evolution.alignment.SitePatterns)16 TreeLikelihood (dr.oldevomodel.treelikelihood.TreeLikelihood)16 DataType (dr.evolution.datatype.DataType)11 GeneralSubstitutionModel (dr.oldevomodel.substmodel.GeneralSubstitutionModel)6 ExchangeOperator (dr.evomodel.operators.ExchangeOperator)5 SubtreeSlideOperator (dr.evomodel.operators.SubtreeSlideOperator)5 WilsonBalding (dr.evomodel.operators.WilsonBalding)5 ArrayLogFormatter (dr.inference.loggers.ArrayLogFormatter)5 MCLogger (dr.inference.loggers.MCLogger)5 TabDelimitedFormatter (dr.inference.loggers.TabDelimitedFormatter)5 MCMC (dr.inference.mcmc.MCMC)5 MCMCOptions (dr.inference.mcmc.MCMCOptions)5 ArrayTraceList (dr.inference.trace.ArrayTraceList)5 Trace (dr.inference.trace.Trace)5 TraceCorrelation (dr.inference.trace.TraceCorrelation)5 GTR (dr.oldevomodel.substmodel.GTR)5