Examples with Alignment jebl.evolution.alignments.Alignment used on opensource projects

Example 1 with Alignment

use of jebl.evolution.alignments.Alignment 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;
        }
    }
}