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

Example 1 with SiteModel

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

the class RecomboGen method main.

public static void main(String[] argv) {
    // Simulate sequences on this tree to generate sequences at the top of the
    // recombination process.
    Parameter kappa = new Parameter.Default(1, 2);
    Parameter freqs = new Parameter.Default(new double[] { 0.25, 0.25, 0.25, 0.25 });
    FrequencyModel freqModel = new FrequencyModel(Nucleotides.INSTANCE, freqs);
    HKY hky = new HKY(kappa, freqModel);
    // create site model
    SiteModel siteModel = new GammaSiteModel(hky);
    // create branch rate model
    Parameter rate = new Parameter.Default(1, 1.0E-4);
    BranchRateModel branchRateModel = new StrictClockBranchRates(rate);
//        RecomboGen recomboGen = new RecomboGen(8.0E-2, 100, 1000, 1.0E-4,
//                new String[] { "taxon1_0", "taxon2_0", "taxon3_10", "taxon4_10", "taxon5_20", "taxon6_20"},
//                new double[] { 0, 0, 10, 10, 20, 20} );
//
//        recomboGen.generate();
}

Example 2 with SiteModel

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

the class Coevolve method runModel2.

private void runModel2(PatternList patternList, PrintWriter pw, Tree tree, SubstitutionModel substModel, final Parameter betaParameter) {
    final Parameter muParameter = new Parameter.Default(1.0);
    muParameter.setId("mu");
    SiteModel siteModel = new GammaSiteModel(substModel, muParameter, null, 1, null);
    TreeModel treeModel = new TreeModel(tree);
    final TreeLikelihood treeLikelihood = new TreeLikelihood(patternList, treeModel, siteModel, null, null, false, false, true, false, false);
    treeLikelihood.setId("likelihood");
    MultivariateFunction function = new MultivariateFunction() {

        public double evaluate(double[] argument) {
            betaParameter.setParameterValue(0, argument[0]);
            betaParameter.setParameterValue(1, argument[1]);
            muParameter.setParameterValue(0, argument[2]);
            double lnL = -treeLikelihood.getLogLikelihood();
            // System.err.println("" + argument[0] + "\t" + argument[1] + "\t" + argument[2] + "\t" + lnL);
            return lnL;
        }

        public int getNumArguments() {
            return 3;
        }

        public double getLowerBound(int n) {
            return 0.0;
        }

        public double getUpperBound(int n) {
            return 100.0;
        }
    };
    MultivariateMinimum optimizer = new ConjugateGradientSearch();
    double lnL = optimizer.findMinimum(function, new double[] { 1.0, 1.0, 1.0 }, 6, 6);
    pw.write(betaParameter.getParameterValue(0) + "\t");
    pw.write(betaParameter.getParameterValue(1) + "\t");
    pw.write(muParameter.getParameterValue(0) + "\t");
    pw.write(lnL + "\n");
    pw.flush();
    System.out.println("" + betaParameter.getParameterValue(0) + "\t" + betaParameter.getParameterValue(1) + "\t" + muParameter.getParameterValue(0) + "\t" + lnL);
}

Example 3 with SiteModel

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

the class SingleTipObservationProcessParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    Parameter mu = (Parameter) xo.getElementFirstChild(AnyTipObservationProcessParser.DEATH_RATE);
    Parameter lam = (Parameter) xo.getElementFirstChild(AnyTipObservationProcessParser.IMMIGRATION_RATE);
    TreeModel treeModel = (TreeModel) xo.getChild(TreeModel.class);
    PatternList patterns = (PatternList) xo.getChild(PatternList.class);
    Taxon sourceTaxon = (Taxon) xo.getChild(Taxon.class);
    SiteModel siteModel = (SiteModel) xo.getChild(SiteModel.class);
    BranchRateModel branchRateModel = (BranchRateModel) xo.getChild(BranchRateModel.class);
    Logger.getLogger("dr.evomodel.MSSD").info("Creating SingleTipObservationProcess model. All traits are assumed extant in " + sourceTaxon.getId() + "Initial mu = " + mu.getParameterValue(0) + " initial lam = " + lam.getParameterValue(0));
    return new SingleTipObservationProcess(treeModel, patterns, siteModel, branchRateModel, mu, lam, sourceTaxon);
}

Example 4 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;
        }
    }
}

Example 5 with SiteModel

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

the class SequenceSimulator method main.

// getDefaultSiteModel
public static void main(String[] args) {
    try {
        int nReplications = 10;
        // create tree
        NewickImporter importer = new NewickImporter("((A:1.0,B:1.0)AB:1.0,(C:1.0,D:1.0)CD:1.0)ABCD;");
        Tree tree = importer.importTree(null);
        // create site model
        SiteModel siteModel = getDefaultSiteModel();
        // create branch rate model
        BranchRateModel branchRateModel = new DefaultBranchRateModel();
        // feed to sequence simulator and generate leaves
        SequenceSimulator treeSimulator = new SequenceSimulator(tree, siteModel, branchRateModel, nReplications);
        Sequence ancestralSequence = new Sequence();
        ancestralSequence.appendSequenceString("TCAGGTCAAG");
        treeSimulator.setAncestralSequence(ancestralSequence);
        System.out.println(treeSimulator.simulate().toString());
    } catch (Exception e) {
        e.printStackTrace();
    }
//END: try-catch block
}