Search in sources :

Example 1 with BeagleSequenceSimulator

use of dr.app.beagle.tools.BeagleSequenceSimulator in project beast-mcmc by beast-dev.

the class BeagleTreeLikelihood method main.

public static void main(String[] args) {
    try {
        MathUtils.setSeed(666);
        System.out.println("Test case 1: simulateOnePartition");
        int sequenceLength = 1000;
        ArrayList<Partition> partitionsList = new ArrayList<Partition>();
        // create tree
        NewickImporter importer = new NewickImporter("(SimSeq1:73.7468,(SimSeq2:25.256989999999995,SimSeq3:45.256989999999995):18.48981);");
        Tree tree = importer.importTree(null);
        TreeModel treeModel = new TreeModel(tree);
        // create Frequency Model
        Parameter freqs = new Parameter.Default(new double[] { 0.25, 0.25, 0.25, 0.25 });
        FrequencyModel freqModel = new FrequencyModel(Nucleotides.INSTANCE, freqs);
        // create branch model
        Parameter kappa1 = new Parameter.Default(1, 1);
        Parameter kappa2 = new Parameter.Default(1, 1);
        HKY hky1 = new HKY(kappa1, freqModel);
        HKY hky2 = new HKY(kappa2, freqModel);
        HomogeneousBranchModel homogenousBranchSubstitutionModel = new HomogeneousBranchModel(hky1);
        List<SubstitutionModel> substitutionModels = new ArrayList<SubstitutionModel>();
        substitutionModels.add(hky1);
        substitutionModels.add(hky2);
        List<FrequencyModel> freqModels = new ArrayList<FrequencyModel>();
        freqModels.add(freqModel);
        Parameter epochTimes = new Parameter.Default(1, 20);
        // create branch rate model
        Parameter rate = new Parameter.Default(1, 0.001);
        BranchRateModel branchRateModel = new StrictClockBranchRates(rate);
        // create site model
        GammaSiteRateModel siteRateModel = new GammaSiteRateModel("siteModel");
        BranchModel homogeneousBranchModel = new HomogeneousBranchModel(hky1);
        BranchModel epochBranchModel = new EpochBranchModel(treeModel, substitutionModels, epochTimes);
        // create partition
        Partition partition1 = new //
        Partition(//
        treeModel, //
        homogenousBranchSubstitutionModel, //
        siteRateModel, //
        branchRateModel, //
        freqModel, // from
        0, // to
        sequenceLength - 1, // every
        1);
        partitionsList.add(partition1);
        // feed to sequence simulator and generate data
        BeagleSequenceSimulator simulator = new BeagleSequenceSimulator(partitionsList);
        Alignment alignment = simulator.simulate(false, false);
        BeagleTreeLikelihood nbtl = new BeagleTreeLikelihood(alignment, treeModel, homogeneousBranchModel, siteRateModel, branchRateModel, null, false, PartialsRescalingScheme.DEFAULT, false);
        System.out.println("nBTL(homogeneous) = " + nbtl.getLogLikelihood());
        nbtl = new BeagleTreeLikelihood(alignment, treeModel, epochBranchModel, siteRateModel, branchRateModel, null, false, PartialsRescalingScheme.DEFAULT, false);
        System.out.println("nBTL(epoch) = " + nbtl.getLogLikelihood());
    } catch (Exception e) {
        e.printStackTrace();
        System.exit(-1);
    }
// END: try-catch block
}
Also used : FrequencyModel(dr.evomodel.substmodel.FrequencyModel) Partition(dr.app.beagle.tools.Partition) EpochBranchModel(dr.evomodel.branchmodel.EpochBranchModel) HomogeneousBranchModel(dr.evomodel.branchmodel.HomogeneousBranchModel) GammaSiteRateModel(dr.evomodel.siteratemodel.GammaSiteRateModel) HomogeneousBranchModel(dr.evomodel.branchmodel.HomogeneousBranchModel) BranchModel(dr.evomodel.branchmodel.BranchModel) EpochBranchModel(dr.evomodel.branchmodel.EpochBranchModel) StrictClockBranchRates(dr.evomodel.branchratemodel.StrictClockBranchRates) SubstitutionModel(dr.evomodel.substmodel.SubstitutionModel) MarkovModulatedSubstitutionModel(dr.evomodel.substmodel.MarkovModulatedSubstitutionModel) BeagleSequenceSimulator(dr.app.beagle.tools.BeagleSequenceSimulator) TreeModel(dr.evomodel.tree.TreeModel) Alignment(dr.evolution.alignment.Alignment) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) NewickImporter(dr.evolution.io.NewickImporter) HKY(dr.evomodel.substmodel.nucleotide.HKY) Tree(dr.evolution.tree.Tree) Parameter(dr.inference.model.Parameter)

Example 2 with BeagleSequenceSimulator

use of dr.app.beagle.tools.BeagleSequenceSimulator in project beast-mcmc by beast-dev.

the class MainFrame method generateNumberOfSimulations.

// END: doExport
// threading, UI, exceptions handling
private void generateNumberOfSimulations(final File outFile) {
    setBusy();
    SwingWorker<Void, Void> worker = new SwingWorker<Void, Void>() {

        ArrayList<TreeModel> simulatedTreeModelList = new ArrayList<TreeModel>();

        // Executed in background thread
        public Void doInBackground() {
            try {
                if (BeagleSequenceSimulatorApp.VERBOSE) {
                    Utils.printPartitionDataList(dataList);
                    System.out.println();
                }
                long startingSeed = dataList.startingSeed;
                for (int i = 0; i < dataList.simulationsCount; i++) {
                    String fullPath = Utils.getMultipleWritePath(outFile, dataList.outputFormat.toString().toLowerCase(), i);
                    PrintWriter writer = new PrintWriter(new FileWriter(fullPath));
                    ArrayList<Partition> partitionsList = new ArrayList<Partition>();
                    for (PartitionData data : dataList) {
                        if (data.record == null) {
                            writer.close();
                            throw new RuntimeException("Set data in Partitions tab for " + (partitionsList.size() + 1) + " partition.");
                        } else {
                            TreeModel treeModel = data.createTreeModel();
                            simulatedTreeModelList.add(treeModel);
                            // create partition
                            Partition partition = new Partition(// 
                            treeModel, // 
                            data.createBranchModel(), // 
                            data.createSiteRateModel(), // 
                            data.createClockRateModel(), // 
                            data.createFrequencyModel(), // from
                            data.from - 1, // to
                            data.to - 1, // every
                            data.every);
                            if (data.ancestralSequenceString != null) {
                                partition.setRootSequence(data.createAncestralSequence());
                            }
                            partitionsList.add(partition);
                        }
                    }
                    if (dataList.setSeed) {
                        MathUtils.setSeed(startingSeed);
                        startingSeed += 1;
                    }
                    beagleSequenceSimulator = new BeagleSequenceSimulator(partitionsList);
                    SimpleAlignment alignment = beagleSequenceSimulator.simulate(dataList.useParallel, dataList.outputAncestralSequences);
                    alignment.setOutputType(dataList.outputFormat);
                    // if (dataList.outputFormat == SimpleAlignment.OutputType.NEXUS) {
                    // alignment.setOutputType(dataList.outputFormat);
                    // } else if(dataList.outputFormat == SimpleAlignment.OutputType.XML) {
                    // alignment.setOutputType(dataList.outputFormat);
                    // }else {
                    // //
                    // }
                    writer.println(alignment.toString());
                    writer.close();
                }
            // END: simulationsCount loop
            } catch (Exception e) {
                Utils.handleException(e);
                setStatus("Exception occured.");
                setIdle();
            }
            return null;
        }

        // END: doInBackground
        // Executed in event dispatch thread
        public void done() {
            // LinkedHashMap<Integer, LinkedHashMap<NodeRef, int[]>> partitionSequencesMap = beagleSequenceSimulator.getPartitionSequencesMap();
            terminalPanel.setText(Utils.partitionDataListToString(dataList, simulatedTreeModelList));
            setStatus("Generated " + Utils.getSiteCount(dataList) + " sites.");
            setIdle();
        }
    };
    worker.execute();
}
Also used : Partition(dr.app.beagle.tools.Partition) FileWriter(java.io.FileWriter) ArrayList(java.util.ArrayList) BeagleSequenceSimulator(dr.app.beagle.tools.BeagleSequenceSimulator) IOException(java.io.IOException) FileNotFoundException(java.io.FileNotFoundException) TreeModel(dr.evomodel.tree.TreeModel) SimpleAlignment(dr.evolution.alignment.SimpleAlignment) SwingWorker(javax.swing.SwingWorker) PrintWriter(java.io.PrintWriter)

Example 3 with BeagleSequenceSimulator

use of dr.app.beagle.tools.BeagleSequenceSimulator in project beast-mcmc by beast-dev.

the class BeagleSeqSimTest method simulateAminoAcid.

// END: simulateThreePartitions
static void simulateAminoAcid() {
    try {
        System.out.println("Test case 4: simulateAminoAcid");
        MathUtils.setSeed(666);
        int sequenceLength = 10;
        ArrayList<Partition> partitionsList = new ArrayList<Partition>();
        // create tree
        NewickImporter importer = new NewickImporter("(SimSeq1:73.7468,(SimSeq2:25.256989999999995,SimSeq3:45.256989999999995):18.48981);");
        Tree tree = importer.importTree(null);
        DefaultTreeModel treeModel = new DefaultTreeModel(tree);
        // create site model
        GammaSiteRateModel siteRateModel = new GammaSiteRateModel("siteModel");
        // create branch rate model
        BranchRateModel branchRateModel = new DefaultBranchRateModel();
        // create Frequency Model
        Parameter freqs = new Parameter.Default(new double[] { 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05 });
        FrequencyModel freqModel = new FrequencyModel(AminoAcids.INSTANCE, freqs);
        // create substitution model
        EmpiricalRateMatrix rateMatrix = Blosum62.INSTANCE;
        EmpiricalAminoAcidModel empiricalAminoAcidModel = new EmpiricalAminoAcidModel(rateMatrix, freqModel);
        HomogeneousBranchModel substitutionModel = new HomogeneousBranchModel(empiricalAminoAcidModel);
        // create partition
        Partition partition1 = new // 
        Partition(// 
        treeModel, // 
        substitutionModel, // 
        siteRateModel, // 
        branchRateModel, // 
        freqModel, // from
        0, // to
        sequenceLength - 1, // every
        1);
        partitionsList.add(partition1);
        // feed to sequence simulator and generate data
        BeagleSequenceSimulator simulator = new BeagleSequenceSimulator(partitionsList);
        System.out.println(simulator.simulate(simulateInPar, false).toString());
    } catch (Exception e) {
        e.printStackTrace();
        System.exit(-1);
    }
// END: try-catch
}
Also used : FrequencyModel(dr.evomodel.substmodel.FrequencyModel) Partition(dr.app.beagle.tools.Partition) EmpiricalRateMatrix(dr.evomodel.substmodel.EmpiricalRateMatrix) EmpiricalAminoAcidModel(dr.evomodel.substmodel.aminoacid.EmpiricalAminoAcidModel) ArrayList(java.util.ArrayList) HomogeneousBranchModel(dr.evomodel.branchmodel.HomogeneousBranchModel) GammaSiteRateModel(dr.evomodel.siteratemodel.GammaSiteRateModel) DefaultTreeModel(dr.evomodel.tree.DefaultTreeModel) BeagleSequenceSimulator(dr.app.beagle.tools.BeagleSequenceSimulator) ImportException(dr.evolution.io.Importer.ImportException) IOException(java.io.IOException) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) NewickImporter(dr.evolution.io.NewickImporter) Tree(dr.evolution.tree.Tree) Parameter(dr.inference.model.Parameter)

Example 4 with BeagleSequenceSimulator

use of dr.app.beagle.tools.BeagleSequenceSimulator in project beast-mcmc by beast-dev.

the class BeagleSeqSimTest method simulateTwoPartitions.

// END: simulateOnePartition
static void simulateTwoPartitions() {
    try {
        System.out.println("Test case 3: simulateTwoPartitions");
        MathUtils.setSeed(666);
        int sequenceLength = 11;
        ArrayList<Partition> partitionsList = new ArrayList<Partition>();
        // create tree
        NewickImporter importer = new NewickImporter("(SimSeq1:73.7468,(SimSeq2:25.256989999999995,SimSeq3:45.256989999999995):18.48981);");
        Tree tree = importer.importTree(null);
        DefaultTreeModel treeModel = new DefaultTreeModel(tree);
        // create Frequency Model
        Parameter freqs = new Parameter.Default(new double[] { 0.25, 0.25, 0.25, 0.25 });
        FrequencyModel freqModel = new FrequencyModel(Nucleotides.INSTANCE, freqs);
        // create substitution model
        Parameter kappa = new Parameter.Default(1, 10);
        HKY hky = new HKY(kappa, freqModel);
        HomogeneousBranchModel substitutionModel = new HomogeneousBranchModel(hky);
        // create site model
        GammaSiteRateModel siteRateModel = new GammaSiteRateModel("siteModel");
        // create branch rate model
        BranchRateModel branchRateModel = new DefaultBranchRateModel();
        // create partition
        Partition partition1 = new // 
        Partition(// 
        treeModel, // 
        substitutionModel, // 
        siteRateModel, // 
        branchRateModel, // 
        freqModel, // from
        0, // to
        3, // every
        1);
        // create partition
        Partition Partition = new // 
        Partition(// 
        treeModel, // 
        substitutionModel, // 
        siteRateModel, // 
        branchRateModel, // 
        freqModel, // from
        4, // to
        sequenceLength - 1, // every
        1);
        Sequence ancestralSequence = new Sequence();
        ancestralSequence.appendSequenceString("TCAAGTG");
        Partition.setRootSequence(ancestralSequence);
        partitionsList.add(partition1);
        partitionsList.add(Partition);
        // feed to sequence simulator and generate data
        BeagleSequenceSimulator simulator = new BeagleSequenceSimulator(partitionsList);
        System.out.println(simulator.simulate(simulateInPar, false).toString());
    } catch (Exception e) {
        e.printStackTrace();
        System.exit(-1);
    }
// END: try-catch block
}
Also used : FrequencyModel(dr.evomodel.substmodel.FrequencyModel) Partition(dr.app.beagle.tools.Partition) ArrayList(java.util.ArrayList) HomogeneousBranchModel(dr.evomodel.branchmodel.HomogeneousBranchModel) GammaSiteRateModel(dr.evomodel.siteratemodel.GammaSiteRateModel) DefaultTreeModel(dr.evomodel.tree.DefaultTreeModel) Sequence(dr.evolution.sequence.Sequence) BeagleSequenceSimulator(dr.app.beagle.tools.BeagleSequenceSimulator) ImportException(dr.evolution.io.Importer.ImportException) IOException(java.io.IOException) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) NewickImporter(dr.evolution.io.NewickImporter) HKY(dr.evomodel.substmodel.nucleotide.HKY) Tree(dr.evolution.tree.Tree) Parameter(dr.inference.model.Parameter)

Example 5 with BeagleSequenceSimulator

use of dr.app.beagle.tools.BeagleSequenceSimulator in project beast-mcmc by beast-dev.

the class BeagleSeqSimTest method simulateRandomBranchAssignment.

// END: main
static void simulateRandomBranchAssignment() {
    try {
        System.out.println("Test case I dunno which: simulate random branch assignments");
        MathUtils.setSeed(666);
        int sequenceLength = 10;
        ArrayList<Partition> partitionsList = new ArrayList<Partition>();
        File treeFile = new File("/home/filip/Dropbox/BeagleSequenceSimulator/SimTree/SimTree.figtree");
        Tree tree = Utils.importTreeFromFile(treeFile);
        DefaultTreeModel treeModel = new DefaultTreeModel(tree);
        // create Frequency Model
        Parameter freqs = new Parameter.Default(Utils.UNIFORM_CODON_FREQUENCIES);
        FrequencyModel freqModel = new FrequencyModel(Codons.UNIVERSAL, freqs);
        // create base subst model
        Parameter omegaParameter = new Parameter.Default("omega", 1, 1.0);
        Parameter kappaParameter = new Parameter.Default("kappa", 1, 1.0);
        GY94CodonModel baseSubModel = new GY94CodonModel(Codons.UNIVERSAL, omegaParameter, kappaParameter, freqModel);
        RandomBranchModel substitutionModel = new RandomBranchModel(treeModel, baseSubModel, 0.25, false, -1);
        // create site model
        GammaSiteRateModel siteRateModel = new GammaSiteRateModel("siteModel");
        // create branch rate model
        BranchRateModel branchRateModel = new DefaultBranchRateModel();
        // create partition
        Partition partition1 = new // 
        Partition(// 
        treeModel, // 
        substitutionModel, // 
        siteRateModel, // 
        branchRateModel, // 
        freqModel, // from
        0, // to
        sequenceLength - 1, // every
        1);
        // Sequence ancestralSequence = new Sequence();
        // ancestralSequence.appendSequenceString("TCAAGTGAGG");
        // partition1.setRootSequence(ancestralSequence);
        partitionsList.add(partition1);
        // feed to sequence simulator and generate data
        BeagleSequenceSimulator simulator = new BeagleSequenceSimulator(partitionsList);
        SimpleAlignment alignment = simulator.simulate(simulateInPar, true);
        // alignment.setOutputType(SimpleAlignment.OutputType.NEXUS);
        alignment.setOutputType(SimpleAlignment.OutputType.XML);
        System.out.println(alignment.toString());
    } catch (Exception e) {
        e.printStackTrace();
    }
// END: try-catch
}
Also used : FrequencyModel(dr.evomodel.substmodel.FrequencyModel) Partition(dr.app.beagle.tools.Partition) ArrayList(java.util.ArrayList) GammaSiteRateModel(dr.evomodel.siteratemodel.GammaSiteRateModel) DefaultTreeModel(dr.evomodel.tree.DefaultTreeModel) BeagleSequenceSimulator(dr.app.beagle.tools.BeagleSequenceSimulator) ImportException(dr.evolution.io.Importer.ImportException) IOException(java.io.IOException) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) RandomBranchModel(dr.evomodel.branchmodel.RandomBranchModel) SimpleAlignment(dr.evolution.alignment.SimpleAlignment) DefaultBranchRateModel(dr.evomodel.branchratemodel.DefaultBranchRateModel) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel) Tree(dr.evolution.tree.Tree) Parameter(dr.inference.model.Parameter) GY94CodonModel(dr.evomodel.substmodel.codon.GY94CodonModel) File(java.io.File)

Aggregations

BeagleSequenceSimulator (dr.app.beagle.tools.BeagleSequenceSimulator)13 Partition (dr.app.beagle.tools.Partition)13 Tree (dr.evolution.tree.Tree)10 BranchRateModel (dr.evomodel.branchratemodel.BranchRateModel)10 GammaSiteRateModel (dr.evomodel.siteratemodel.GammaSiteRateModel)10 FrequencyModel (dr.evomodel.substmodel.FrequencyModel)10 Parameter (dr.inference.model.Parameter)10 ArrayList (java.util.ArrayList)10 NewickImporter (dr.evolution.io.NewickImporter)9 HomogeneousBranchModel (dr.evomodel.branchmodel.HomogeneousBranchModel)9 DefaultBranchRateModel (dr.evomodel.branchratemodel.DefaultBranchRateModel)9 IOException (java.io.IOException)9 ImportException (dr.evolution.io.Importer.ImportException)8 DefaultTreeModel (dr.evomodel.tree.DefaultTreeModel)8 SimpleAlignment (dr.evolution.alignment.SimpleAlignment)6 HKY (dr.evomodel.substmodel.nucleotide.HKY)6 Alignment (dr.evolution.alignment.Alignment)4 TreeModel (dr.evomodel.tree.TreeModel)4 SubstitutionModel (dr.evomodel.substmodel.SubstitutionModel)3 BeagleTreeLikelihood (dr.evomodel.treelikelihood.BeagleTreeLikelihood)3