Search in sources :

Example 21 with Locus

use of bacter.Locus in project bacter by tgvaughan.

the class SimulatedACG method generateConversions.

private void generateConversions() {
    // Draw number of conversions:
    int Nconv = (int) Randomizer.nextPoisson(rho * getClonalFrameLength() * (getTotalSequenceLength() + (delta - 1.0) * getLoci().size()));
    // Generate conversions:
    for (int i = 0; i < Nconv; i++) {
        // Choose alignment
        double u = Randomizer.nextDouble() * (getTotalSequenceLength() + (delta - 1.0) * getLoci().size());
        Locus affectedLocus = null;
        for (Locus locus : getLoci()) {
            if (u < locus.getSiteCount() + delta - 1.0) {
                affectedLocus = locus;
                break;
            } else
                u -= locus.getSiteCount() + delta - 1.0;
        }
        if (affectedLocus == null)
            throw new IllegalStateException("Programmer error: " + "locus choice loop fell through.");
        int startSite, endSite;
        if (u < delta) {
            startSite = 0;
        } else {
            startSite = (int) Math.ceil(u - delta);
        }
        endSite = startSite + (int) Randomizer.nextGeometric(1.0 / delta);
        endSite = Math.min(endSite, affectedLocus.getSiteCount() - 1);
        Conversion conv = new Conversion();
        conv.setLocus(affectedLocus);
        conv.setStartSite(startSite);
        conv.setEndSite(endSite);
        associateConversionWithCF(conv);
        addConversion(conv);
    }
}
Also used : Locus(bacter.Locus) Conversion(bacter.Conversion)

Example 22 with Locus

use of bacter.Locus in project bacter by tgvaughan.

the class SimulatedAlignment method initAndValidate.

@Override
public void initAndValidate() {
    acg = acgInput.get();
    siteModel = siteModelInput.get();
    Locus locus;
    if (acg.getLoci().size() == 1)
        locus = acg.getLoci().get(0);
    else {
        locus = locusInput.get();
        if (locus == null)
            throw new IllegalArgumentException("Must specify locus" + " when simulating alignment from ACG associated" + " with multiple loci.");
    }
    // We can't wait for Alignment.initAndValidate() to get the
    // data type for us.
    grabDataType();
    // Simulate alignment
    simulate(locus);
    super.initAndValidate();
    // Write simulated alignment to disk if requested:
    if (outputFileNameInput.get() != null) {
        try (PrintStream pstream = new PrintStream(outputFileNameInput.get())) {
            if (useNexusInput.get()) {
                NexusBuilder nb = new NexusBuilder();
                nb.append(new TaxaBlock(acg.getTaxonset()));
                nb.append(new CharactersBlock(this));
                nb.write(pstream);
            } else {
                for (String taxonName : acg.getTaxaNames()) {
                    pstream.print(">" + taxonName + "\n");
                    pstream.print(getSequenceAsString(taxonName) + "\n");
                }
            }
        } catch (FileNotFoundException e) {
            e.printStackTrace();
        }
    }
}
Also used : PrintStream(java.io.PrintStream) CharactersBlock(feast.nexus.CharactersBlock) FileNotFoundException(java.io.FileNotFoundException) TaxaBlock(feast.nexus.TaxaBlock) Locus(bacter.Locus) NexusBuilder(feast.nexus.NexusBuilder)

Example 23 with Locus

use of bacter.Locus in project bacter by tgvaughan.

the class ACGOperator method disconnectEdge.

/**
 * Disconnect edge <node, node.parent> from its sister and
 * grandparent (if the grandparent exists), forming a new edge
 * <sister, grandparent>. All conversions originally above node.parent
 * are re-attached to sister.
 *
 * Conversions on edge above node are not modified.
 *
 * @param node base of edge to detach.
 */
protected void disconnectEdge(Node node) {
    if (node.isRoot())
        throw new IllegalArgumentException("Programmer error: " + "root argument passed to disconnectEdge().");
    Node parent = node.getParent();
    Node sister = getSibling(node);
    if (parent.isRoot()) {
        parent.removeChild(sister);
        sister.setParent(null);
    } else {
        Node grandParent = parent.getParent();
        grandParent.removeChild(parent);
        parent.setParent(null);
        parent.removeChild(sister);
        grandParent.addChild(sister);
    }
    for (Locus locus : acg.getLoci()) {
        for (Conversion conv : acg.getConversions(locus)) {
            if (conv.getNode1() == parent)
                conv.setNode1(sister);
            if (conv.getNode2() == parent)
                conv.setNode2(sister);
        }
    }
}
Also used : Node(beast.evolution.tree.Node) Locus(bacter.Locus) Conversion(bacter.Conversion)

Example 24 with Locus

use of bacter.Locus in project bacter by tgvaughan.

the class AddRemoveDetour method addDetour.

/**
 * Detour addition variant of move.
 *
 * @return log HGF
 */
private double addDetour() {
    double logHGF = 0.0;
    if (acg.getTotalConvCount() == 0)
        return Double.NEGATIVE_INFINITY;
    // Select conversion at random
    Conversion conv = chooseConversion();
    logHGF -= Math.log(1.0 / acg.getTotalConvCount());
    // Select detour times:
    double t1 = conv.getHeight1() + Randomizer.nextDouble() * (conv.getHeight2() - conv.getHeight1());
    double t2 = conv.getHeight1() + Randomizer.nextDouble() * (conv.getHeight2() - conv.getHeight1());
    double tLower, tUpper;
    if (t1 < t2) {
        tLower = t1;
        tUpper = t2;
    } else {
        tLower = t2;
        tUpper = t1;
    }
    logHGF -= Math.log(2.0) - 2.0 * Math.log(conv.getHeight2() - conv.getHeight1());
    // Select non-root node below detour edge at random
    Node detour = acg.getNode(Randomizer.nextInt(acg.getNodeCount() - 1));
    logHGF -= Math.log(1.0 / (acg.getNodeCount() - 1));
    // Abort if selected detour edge does not contain tLower and tUpper
    if (detour.getHeight() > tLower || detour.getParent().getHeight() < tUpper)
        return Double.NEGATIVE_INFINITY;
    // Abort if conv is already attached to selected detour edge:
    if (detour == conv.getNode1() || detour == conv.getNode2())
        return Double.NEGATIVE_INFINITY;
    Conversion convA = new Conversion();
    convA.setLocus(conv.getLocus());
    convA.setNode1(conv.getNode1());
    convA.setHeight1(conv.getHeight1());
    convA.setNode2(detour);
    convA.setHeight2(tLower);
    convA.setStartSite(conv.getStartSite());
    convA.setEndSite(conv.getEndSite());
    Conversion convB = new Conversion();
    convB.setNode1(detour);
    convB.setHeight1(tUpper);
    convB.setNode2(conv.getNode2());
    convB.setHeight2(conv.getHeight2());
    logHGF -= drawAffectedRegion(convB);
    acg.deleteConversion(conv);
    acg.addConversion(convA);
    acg.addConversion(convB);
    // Count number of node1s and node2s attached to detour edge
    int node1Count = 0;
    int node2Count = 0;
    for (Locus locus : acg.getLoci()) {
        for (Conversion thisConv : acg.getConversions(locus)) {
            if (thisConv.getNode1() == detour && thisConv.getNode2() != detour)
                node1Count += 1;
            if (thisConv.getNode2() == detour && thisConv.getNode1() != detour)
                node2Count += 1;
        }
    }
    // Incorporate probability of reverse move:
    logHGF += Math.log(1.0 / ((acg.getNodeCount() - 1) * node1Count * node2Count));
    return logHGF;
}
Also used : Node(beast.evolution.tree.Node) Locus(bacter.Locus) Conversion(bacter.Conversion)

Example 25 with Locus

use of bacter.Locus in project bacter by tgvaughan.

the class AddRemoveRedundantConversion method drawAffectedRegion.

private double drawAffectedRegion(Conversion conv) {
    double logP = 0.0;
    Locus locus = acg.getLoci().get(Randomizer.nextInt(acg.getLoci().size()));
    conv.setLocus(locus);
    logP += Math.log(1.0 / acg.getLoci().size());
    if (!acg.wholeLocusModeOn()) {
        int site1 = Randomizer.nextInt(locus.getSiteCount());
        int site2 = Randomizer.nextInt(locus.getSiteCount());
        if (site1 < site2) {
            conv.setStartSite(site1);
            conv.setEndSite(site2);
        } else {
            conv.setStartSite(site2);
            conv.setEndSite(site1);
        }
        logP += 2.0 * Math.log(1.0 / locus.getSiteCount());
        if (site1 != site2)
            logP += Math.log(2.0);
    } else {
        conv.setStartSite(0);
        conv.setEndSite(locus.getSiteCount() - 1);
    }
    return logP;
}
Also used : Locus(bacter.Locus)

Aggregations

Locus (bacter.Locus)37 Conversion (bacter.Conversion)27 Node (beast.evolution.tree.Node)17 ConversionGraph (bacter.ConversionGraph)10 ArrayList (java.util.ArrayList)7 Test (org.junit.Test)7 RealParameter (beast.core.parameter.RealParameter)6 TaxonSet (beast.evolution.alignment.TaxonSet)5 SiteModel (beast.evolution.sitemodel.SiteModel)5 JukesCantor (beast.evolution.substitutionmodel.JukesCantor)5 ConstantPopulation (beast.evolution.tree.coalescent.ConstantPopulation)5 PrintStream (java.io.PrintStream)4 ClusterTree (beast.util.ClusterTree)3 SimulatedACG (bacter.model.SimulatedACG)2 NexusBuilder (feast.nexus.NexusBuilder)2 TaxaBlock (feast.nexus.TaxaBlock)2 FileNotFoundException (java.io.FileNotFoundException)2 Iterator (java.util.Iterator)2 CFEventList (bacter.CFEventList)1 ACGCoalescent (bacter.model.ACGCoalescent)1