Examples with ConversionGraph bacter.ConversionGraph used on opensource projects

Example 1 with ConversionGraph

use of bacter.ConversionGraph in project bacter by tgvaughan.

the class ACGLikelihoodTest method testLikelihoodCaching.

@Test
public void testLikelihoodCaching() throws Exception {
    ConstantPopulation popFunc = new ConstantPopulation();
    popFunc.initByName("popSize", new RealParameter("1.0"));
    Locus locus = new Locus("locus", 10000);
    TaxonSet taxonSet = getTaxonSet(10);
    ConversionGraph acg = new SimulatedACG();
    acg.initByName("rho", 5.0 / locus.getSiteCount(), "delta", 1000.0, "populationModel", popFunc, "locus", locus, "taxonset", taxonSet);
    State state = new State();
    state.initByName("stateNode", acg);
    state.initialise();
    System.out.println(acg);
    // Site model:
    JukesCantor jc = new JukesCantor();
    jc.initByName();
    SiteModel siteModel = new SiteModel();
    siteModel.initByName("mutationRate", new RealParameter("1"), "substModel", jc);
    // Simulate alignment:
    SimulatedAlignment alignment = new SimulatedAlignment();
    alignment.initByName("acg", acg, "siteModel", siteModel, "outputFileName", "simulated_alignment.nexus", "useNexus", true);
    // Calculate likelihood 1:
    ACGLikelihood argLikelihood = new ACGLikelihood();
    argLikelihood.initByName("locus", locus, "data", alignment, "tree", acg, "siteModel", siteModel);
    ACGLikelihoodSlow argLikelihoodSlow = new ACGLikelihoodSlow();
    argLikelihoodSlow.initByName("locus", locus, "data", alignment, "tree", acg, "siteModel", siteModel);
    double logP1 = argLikelihood.calculateLogP();
    double logP1prime = argLikelihoodSlow.calculateLogP();
    double relError = 2.0 * Math.abs(logP1 - logP1prime) / Math.abs(logP1 + logP1prime);
    System.out.format("logP=%g\nlogPprime=%g\nrelError=%g\n", logP1, logP1prime, relError);
    assertTrue(relError < 1e-13);
    // Add a single recombination event
    Node node1 = acg.getExternalNodes().get(0);
    Node node2 = node1.getParent();
    double height1 = 0.5 * (node1.getHeight() + node1.getParent().getHeight());
    double height2 = 0.5 * (node2.getHeight() + node2.getParent().getHeight());
    int startLocus = 500;
    int endLocus = 600;
    Conversion recomb1 = new Conversion(node1, height1, node2, height2, startLocus, endLocus, acg, locus);
    acg.addConversion(recomb1);
    double logP2 = argLikelihood.calculateLogP();
    double logP2prime = argLikelihoodSlow.calculateLogP();
    relError = 2.0 * Math.abs(logP2 - logP2prime) / Math.abs(logP2 + logP2prime);
    System.out.format("logP=%g\nlogPprime=%g\nrelError=%g\n", logP2, logP2prime, relError);
    assertTrue(relError < 1e-13);
    state.restore();
    double logP3 = argLikelihood.calculateLogP();
    double logP3prime = argLikelihoodSlow.calculateLogP();
    relError = 2.0 * Math.abs(logP3 - logP3prime) / Math.abs(logP3 + logP3prime);
    System.out.format("logP=%g\nlogPprime=%g\nrelError=%g\n", logP3, logP3prime, relError);
    assertTrue(relError < 1e-13);
}

Example 2 with ConversionGraph

use of bacter.ConversionGraph in project bacter by tgvaughan.

the class ACGLikelihoodTest method testLikelihoodAmbiguities.

@Test
public void testLikelihoodAmbiguities() throws Exception {
    Locus locus = new Locus("locus", getAlignment());
    // ConversionGraph
    ConversionGraph acg = new ConversionGraph();
    ClusterTree tree = new ClusterTree();
    tree.initByName("clusterType", "upgma", "taxa", locus.getAlignment());
    acg.assignFrom(tree);
    acg.initByName("locus", locus);
    // Site model:
    JukesCantor jc = new JukesCantor();
    jc.initByName();
    SiteModel siteModel = new SiteModel();
    siteModel.initByName("substModel", jc);
    // Likelihood
    ACGLikelihood argLikelihood = new ACGLikelihood();
    argLikelihood.initByName("locus", locus, "tree", acg, "siteModel", siteModel, "useAmbiguities", true);
    ACGLikelihoodSlow argLikelihoodSlow = new ACGLikelihoodSlow();
    argLikelihoodSlow.initByName("locus", locus, "tree", acg, "siteModel", siteModel, "useAmbiguities", true);
    acg.setEverythingDirty(true);
    double logP = argLikelihood.calculateLogP();
    double logPtrue = argLikelihoodSlow.calculateLogP();
    double relativeDiff = Math.abs(2.0 * (logPtrue - logP) / (logPtrue + logP));
    System.out.println("Test 1.  Truth: " + logPtrue + " Value: " + logP);
    assertTrue(relativeDiff < 1e-14);
    // Add a single recombination event
    Node node1 = acg.getExternalNodes().get(0);
    Node node2 = node1.getParent();
    double height1 = 0.5 * (node1.getHeight() + node1.getParent().getHeight());
    double height2 = 0.5 * (node2.getHeight() + node2.getParent().getHeight());
    int startLocus = 100;
    int endLocus = 200;
    Conversion recomb1 = new Conversion(node1, height1, node2, height2, startLocus, endLocus, acg, locus);
    acg.addConversion(recomb1);
    logP = argLikelihood.calculateLogP();
    logPtrue = argLikelihoodSlow.calculateLogP();
    relativeDiff = Math.abs(2.0 * (logPtrue - logP) / (logPtrue + logP));
    assertTrue(relativeDiff < 1e-14);
    // Add another recombination event
    node1 = acg.getExternalNodes().get(0);
    node2 = acg.getNode(20);
    height1 = 0.75 * (node1.getHeight() + node1.getParent().getHeight());
    height2 = 0.5 * (node2.getHeight() + node2.getParent().getHeight());
    startLocus = 250;
    endLocus = 300;
    Conversion recomb2 = new Conversion(node1, height1, node2, height2, startLocus, endLocus, acg, locus);
    acg.addConversion(recomb2);
    logP = argLikelihood.calculateLogP();
    logPtrue = argLikelihoodSlow.calculateLogP();
    relativeDiff = Math.abs(2.0 * (logPtrue - logP) / (logPtrue + logP));
    assertTrue(relativeDiff < 1e-14);
}

Example 3 with ConversionGraph

use of bacter.ConversionGraph in project bacter by tgvaughan.

the class AddRemoveConversion method main.

public static void main(String[] args) throws Exception {
    ConversionGraph acg = new ConversionGraph();
    ConstantPopulation popFunc = new ConstantPopulation();
    AddRemoveConversion operator = new AddRemoveConversion();
    operator.initByName("weight", 1.0, "acg", acg, "populationModel", popFunc, "rho", new RealParameter(Double.toString(1.0 / 10000.0)), "delta", new RealParameter("50.0"));
    popFunc.initByName("popSize", new RealParameter("1.0"));
    TaxonSet taxonSet = new TaxonSet();
    taxonSet.taxonsetInput.setValue(new Taxon("t1"), taxonSet);
    taxonSet.taxonsetInput.setValue(new Taxon("t2"), taxonSet);
    Locus locus = new Locus("locus", 10000);
    try (PrintStream ps = new PrintStream("out.txt")) {
        for (int i = 0; i < 100000; i++) {
            acg.initByName("locus", locus, "taxonset", taxonSet, "fromString", "(0:1.0,1:1.0)2:0.0;");
            operator.drawNewConversion();
            ps.println(acg.getConversions(locus).get(0).getStartSite() + " " + acg.getConversions(locus).get(0).getEndSite());
        }
    }
}

Example 4 with ConversionGraph

use of bacter.ConversionGraph in project bacter by tgvaughan.

the class COACGLogFileReader method iterator.

@Override
public Iterator<ConversionGraph> iterator() {
    try {
        reset();
        skipBurnin();
    } catch (XMLStreamException | IOException e) {
        throw new IllegalStateException(e.getMessage());
    }
    return new Iterator<ConversionGraph>() {

        int current = 0;

        @Override
        public boolean hasNext() {
            return current < nACGs;
        }

        @Override
        public ConversionGraph next() {
            current += 1;
            ConversionGraph acg = new ConversionGraph();
            for (Locus locus : loci) acg.lociInput.setValue(locus, acg);
            try {
                acg.initAndValidate();
            } catch (Exception e) {
                throw new IllegalStateException(e.getMessage());
            }
            String newick = null;
            List<Integer> rStarts = new ArrayList<>();
            List<Integer> rEnds = new ArrayList<>();
            List<Integer> rEFroms = new ArrayList<>();
            List<Integer> rETos = new ArrayList<>();
            List<Double> rAFroms = new ArrayList<>();
            List<Double> rATos = new ArrayList<>();
            try {
                skipUntil("iteration");
                while (xmlStreamReader.hasNext()) {
                    int type = xmlStreamReader.next();
                    if (type == XMLStreamConstants.END_ELEMENT && xmlStreamReader.getLocalName().toLowerCase().equals("iteration"))
                        break;
                    if (type != XMLStreamConstants.START_ELEMENT)
                        continue;
                    switch(xmlStreamReader.getLocalName().toLowerCase()) {
                        case "tree":
                            xmlStreamReader.next();
                            newick = xmlStreamReader.getText().trim();
                            break;
                        case "recedge":
                            while ((type = xmlStreamReader.next()) != XMLStreamConstants.END_ELEMENT || !xmlStreamReader.getLocalName().toLowerCase().equals("recedge")) {
                                if (type != XMLStreamConstants.START_ELEMENT)
                                    continue;
                                switch(xmlStreamReader.getLocalName().toLowerCase()) {
                                    case "start":
                                        xmlStreamReader.next();
                                        rStarts.add(Integer.valueOf(xmlStreamReader.getText().trim()));
                                        break;
                                    case "end":
                                        xmlStreamReader.next();
                                        rEnds.add(Integer.valueOf(xmlStreamReader.getText().trim()) - 1);
                                        break;
                                    case "efrom":
                                        xmlStreamReader.next();
                                        rEFroms.add(Integer.valueOf(xmlStreamReader.getText().trim()));
                                        break;
                                    case "eto":
                                        xmlStreamReader.next();
                                        rETos.add(Integer.valueOf(xmlStreamReader.getText().trim()));
                                        break;
                                    case "afrom":
                                        xmlStreamReader.next();
                                        rAFroms.add(Double.valueOf(xmlStreamReader.getText().trim()));
                                        break;
                                    case "ato":
                                        xmlStreamReader.next();
                                        rATos.add(Double.valueOf(xmlStreamReader.getText().trim()));
                                        break;
                                    default:
                                }
                            }
                            break;
                        default:
                    }
                }
            } catch (XMLStreamException e) {
                throw new IllegalStateException(e.getMessage());
            }
            acg.fromExtendedNewick(newick, true, 0);
            for (int i = 0; i < rStarts.size(); i++) {
                Node fromNode = acg.getNode(rEFroms.get(i));
                Node toNode = acg.getNode(rETos.get(i));
                Conversion conv = new Conversion(toNode, toNode.getHeight() + rATos.get(i), fromNode, fromNode.getHeight() + rAFroms.get(i), rStarts.get(i), rEnds.get(i), acg, loci.get(0));
                acg.addConversion(conv);
            }
            current += 1;
            return acg;
        }
    };
}

Example 5 with ConversionGraph

use of bacter.ConversionGraph in project bacter by tgvaughan.

the class COACGLogFileReader method main.

public static void main(String[] args) throws IOException, XMLStreamException {
    COACGLogFileReader reader = new COACGLogFileReader(new File("/home/tvaughan/articles/bacter-paper/simulation_studies/robustness/blah.xml"), 0);
    System.out.println("Iterations: " + reader.nACGs);
    for (ConversionGraph acg : reader) System.out.println(acg);
}