Examples with SiteModel dr.oldevomodel.sitemodel.SiteModel used on opensource projects

Example 6 with SiteModel

use of dr.oldevomodel.sitemodel.SiteModel in project beast-mcmc by beast-dev.

the class TreeLikelihoodParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    boolean useAmbiguities = xo.getAttribute(USE_AMBIGUITIES, false);
    boolean allowMissingTaxa = xo.getAttribute(ALLOW_MISSING_TAXA, false);
    boolean storePartials = xo.getAttribute(STORE_PARTIALS, true);
    boolean forceJavaCore = xo.getAttribute(FORCE_JAVA_CORE, false);
    if (Boolean.valueOf(System.getProperty("java.only"))) {
        forceJavaCore = true;
    }
    PatternList patternList = (PatternList) xo.getChild(PatternList.class);
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    SiteModel siteModel = (SiteModel) xo.getChild(SiteModel.class);
    BranchRateModel branchRateModel = (BranchRateModel) xo.getChild(BranchRateModel.class);
    TipStatesModel tipStatesModel = (TipStatesModel) xo.getChild(TipStatesModel.class);
    if (tipStatesModel != null && tipStatesModel.getPatternList() != null) {
        throw new XMLParseException("The same sequence error model cannot be used for multiple partitions");
    }
    if (tipStatesModel != null && tipStatesModel.getModelType() == TipStatesModel.Type.STATES) {
        throw new XMLParseException("The state emitting TipStateModel requires BEAGLE");
    }
    boolean forceRescaling = xo.getAttribute(FORCE_RESCALING, false);
    return new TreeLikelihood(patternList, treeModel, siteModel, branchRateModel, tipStatesModel, useAmbiguities, allowMissingTaxa, storePartials, forceJavaCore, forceRescaling);
}

Example 7 with SiteModel

use of dr.oldevomodel.sitemodel.SiteModel in project beast-mcmc by beast-dev.

the class ALSTreeLikelihoodParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    boolean useAmbiguities = false;
    boolean storePartials = true;
    if (xo.hasAttribute(TreeLikelihoodParser.USE_AMBIGUITIES)) {
        useAmbiguities = xo.getBooleanAttribute(TreeLikelihoodParser.USE_AMBIGUITIES);
    }
    if (xo.hasAttribute(TreeLikelihoodParser.STORE_PARTIALS)) {
        storePartials = xo.getBooleanAttribute(TreeLikelihoodParser.STORE_PARTIALS);
    }
    boolean integrateGainRate = xo.getBooleanAttribute(INTEGRATE_GAIN_RATE);
    //AbstractObservationProcess observationProcess = (AbstractObservationProcess) xo.getChild(AbstractObservationProcess.class);
    PatternList patternList = (PatternList) xo.getChild(PatternList.class);
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    SiteModel siteModel = (SiteModel) xo.getChild(SiteModel.class);
    BranchRateModel branchRateModel = (BranchRateModel) xo.getChild(BranchRateModel.class);
    Parameter mu = ((MutationDeathModel) siteModel.getSubstitutionModel()).getDeathParameter();
    Parameter lam;
    if (!integrateGainRate) {
        lam = (Parameter) xo.getElementFirstChild(IMMIGRATION_RATE);
    } else {
        lam = new Parameter.Default("gainRate", 1.0, 0.001, 1.999);
    }
    AbstractObservationProcess observationProcess = null;
    Logger.getLogger("dr.evolution").info("\n ---------------------------------\nCreating ALSTreeLikelihood model.");
    for (int i = 0; i < xo.getChildCount(); ++i) {
        Object cxo = xo.getChild(i);
        if (cxo instanceof XMLObject && ((XMLObject) cxo).getName().equals(OBSERVATION_PROCESS)) {
            if (((XMLObject) cxo).getStringAttribute(OBSERVATION_TYPE).equals("singleTip")) {
                String taxonName = ((XMLObject) cxo).getStringAttribute(OBSERVATION_TAXON);
                Taxon taxon = treeModel.getTaxon(treeModel.getTaxonIndex(taxonName));
                observationProcess = new SingleTipObservationProcess(treeModel, patternList, siteModel, branchRateModel, mu, lam, taxon);
                Logger.getLogger("dr.evolution").info("All traits are assumed extant in " + taxonName);
            } else {
                // "anyTip" observation process
                observationProcess = new AnyTipObservationProcess(ANY_TIP, treeModel, patternList, siteModel, branchRateModel, mu, lam);
                Logger.getLogger("dr.evolution").info("Observed traits are assumed to be extant in at least one tip node.");
            }
            observationProcess.setIntegrateGainRate(integrateGainRate);
        }
    }
    Logger.getLogger("dr.evolution").info("\tIf you publish results using Acquisition-Loss-Mutation (ALS) Model likelihood, please reference Alekseyenko, Lee and Suchard (2008) Syst. Biol 57: 772-784.\n---------------------------------\n");
    boolean forceRescaling = xo.getAttribute(FORCE_RESCALING, false);
    return new ALSTreeLikelihood(observationProcess, patternList, treeModel, siteModel, branchRateModel, useAmbiguities, storePartials, forceRescaling);
}

Example 8 with SiteModel

use of dr.oldevomodel.sitemodel.SiteModel in project beast-mcmc by beast-dev.

the class ALSTreeLikelihoodParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    boolean useAmbiguities = false;
    boolean storePartials = true;
    if (xo.hasAttribute(TreeLikelihoodParser.USE_AMBIGUITIES)) {
        useAmbiguities = xo.getBooleanAttribute(TreeLikelihoodParser.USE_AMBIGUITIES);
    }
    if (xo.hasAttribute(TreeLikelihoodParser.STORE_PARTIALS)) {
        storePartials = xo.getBooleanAttribute(TreeLikelihoodParser.STORE_PARTIALS);
    }
    boolean integrateGainRate = xo.getBooleanAttribute(INTEGRATE_GAIN_RATE);
    //AbstractObservationProcess observationProcess = (AbstractObservationProcess) xo.getChild(AbstractObservationProcess.class);
    PatternList patternList = (PatternList) xo.getChild(PatternList.class);
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    SiteModel siteModel = (SiteModel) xo.getChild(SiteModel.class);
    BranchRateModel branchRateModel = (BranchRateModel) xo.getChild(BranchRateModel.class);
    Parameter mu = ((MutationDeathModel) siteModel.getSubstitutionModel()).getDeathParameter();
    Parameter lam;
    if (!integrateGainRate) {
        lam = (Parameter) xo.getElementFirstChild(IMMIGRATION_RATE);
    } else {
        lam = new Parameter.Default("gainRate", 1.0, 0.001, 1.999);
    }
    AbstractObservationProcess observationProcess = null;
    Logger.getLogger("dr.evolution").info("\n ---------------------------------\nCreating ALSTreeLikelihood model.");
    for (int i = 0; i < xo.getChildCount(); ++i) {
        Object cxo = xo.getChild(i);
        if (cxo instanceof XMLObject && ((XMLObject) cxo).getName().equals(OBSERVATION_PROCESS)) {
            if (((XMLObject) cxo).getStringAttribute(OBSERVATION_TYPE).equals("singleTip")) {
                String taxonName = ((XMLObject) cxo).getStringAttribute(OBSERVATION_TAXON);
                Taxon taxon = treeModel.getTaxon(treeModel.getTaxonIndex(taxonName));
                observationProcess = new SingleTipObservationProcess(treeModel, patternList, siteModel, branchRateModel, mu, lam, taxon);
                Logger.getLogger("dr.evolution").info("All traits are assumed extant in " + taxonName);
            } else {
                // "anyTip" observation process
                observationProcess = new AnyTipObservationProcess(ANY_TIP, treeModel, patternList, siteModel, branchRateModel, mu, lam);
                Logger.getLogger("dr.evolution").info("Observed traits are assumed to be extant in at least one tip node.");
            }
            observationProcess.setIntegrateGainRate(integrateGainRate);
        }
    }
    Logger.getLogger("dr.evolution").info("\tIf you publish results using Acquisition-Loss-Mutation (ALS) Model likelihood, please reference Alekseyenko, Lee and Suchard (2008) Syst. Biol 57: 772-784.\n---------------------------------\n");
    boolean forceRescaling = xo.getAttribute(FORCE_RESCALING, false);
    return new ALSTreeLikelihood(observationProcess, patternList, treeModel, siteModel, branchRateModel, useAmbiguities, storePartials, forceRescaling);
}

Example 9 with SiteModel

use of dr.oldevomodel.sitemodel.SiteModel in project beast-mcmc by beast-dev.

the class AnyTipObservationProcessParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    Parameter mu = (Parameter) xo.getElementFirstChild(DEATH_RATE);
    Parameter lam = (Parameter) xo.getElementFirstChild(IMMIGRATION_RATE);
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    PatternList patterns = (PatternList) xo.getChild(PatternList.class);
    SiteModel siteModel = (SiteModel) xo.getChild(SiteModel.class);
    BranchRateModel branchRateModel = (BranchRateModel) xo.getChild(BranchRateModel.class);
    Logger.getLogger("dr.evomodel.MSSD").info("Creating AnyTipObservationProcess model. Observed traits are assumed to be extant in at least one tip node. Initial mu = " + mu.getParameterValue(0) + " initial lam = " + lam.getParameterValue(0));
    return new AnyTipObservationProcess(MODEL_NAME, treeModel, patterns, siteModel, branchRateModel, mu, lam);
}

Example 10 with SiteModel

use of dr.oldevomodel.sitemodel.SiteModel 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;
        }
    }
}