Search in sources :

Example 11 with Conversion

use of bacter.Conversion in project bacter by tgvaughan.

the class CFConversionSwap method proposal.

@Override
public double proposal() {
    double logHGF = 0.0;
    // Determine whether we can apply this operator:
    if (acg.getLeafNodeCount() < 3 || acg.getTotalConvCount() == 0)
        return Double.NEGATIVE_INFINITY;
    // Acquire list of conversions compatible with swap:
    List<Conversion> compatible = getCompatibleConversions();
    if (compatible.isEmpty())
        return Double.NEGATIVE_INFINITY;
    // Select conversion to replace
    Conversion replaceConversion = compatible.get(Randomizer.nextInt(compatible.size()));
    logHGF -= Math.log(1.0 / compatible.size());
    acg.deleteConversion(replaceConversion);
    Node srcNode = replaceConversion.getNode1();
    Node destNode = replaceConversion.getNode2();
    Node srcNodeP = srcNode.getParent();
    Node srcNodeS = getSibling(srcNode);
    double t_srcNodeP = srcNodeP.getHeight();
    if (destNode == srcNode.getParent())
        destNode = srcNodeS;
    double newTime = replaceConversion.getHeight2();
    // Conversion modification:
    replaceConversion.setNode2(srcNodeS);
    replaceConversion.setHeight2(t_srcNodeP);
    logHGF += getAffectedRegionProb(replaceConversion);
    logHGF -= drawAffectedRegion(replaceConversion);
    acg.addConversion(replaceConversion);
    if (!includeRootInput.get() && (srcNodeP.isRoot() || destNode.isRoot()))
        return Double.NEGATIVE_INFINITY;
    // Perform necessary conversion expansions/collapses:
    if (newTime < t_srcNodeP) {
        logHGF += collapseConversions(srcNode, destNode, newTime);
    } else {
        logHGF -= expandConversions(srcNode, destNode, newTime);
    }
    logHGF += Math.log(1.0 / getCompatibleConversions().size());
    assert !acg.isInvalid() : "CFCS proposed invalid state.";
    return logHGF;
}
Also used : Node(beast.evolution.tree.Node) Conversion(bacter.Conversion)

Example 12 with Conversion

use of bacter.Conversion in project bacter by tgvaughan.

the class DifferenceFromTrueACG method countSampledTrueConversions.

/**
 * Count the number of true conversions which have correspondences on the
 * provided sampled ARG.
 *
 * @param trueACG true arg
 * @param trueClades clades in true arg
 * @param acg sampled arg
 * @param clades clades in sampled arg
 * @param boundaryTol minimum relative error in region boundaries to allow.
 * @param convHist
 * @return number of found conversions
 */
public static int countSampledTrueConversions(ConversionGraph trueACG, Clade[] trueClades, ConversionGraph acg, Clade[] clades, double boundaryTol, double ageTol, Map<Conversion, Integer> convHist) {
    int count = 0;
    for (Conversion trueConv : trueACG.getConversions(trueACG.getLoci().get(0))) {
        Clade trueFromClade = trueClades[trueConv.getNode1().getNr()];
        Clade trueToClade = trueClades[trueConv.getNode1().getNr()];
        for (Conversion conv : acg.getConversions(acg.getLoci().get(0))) {
            Clade fromClade = clades[conv.getNode1().getNr()];
            Clade toClade = clades[conv.getNode1().getNr()];
            if (fromClade.equals(trueFromClade) && toClade.equals(trueToClade)) {
                if (Math.abs(trueConv.getStartSite() - conv.getStartSite()) / trueConv.getSiteCount() <= boundaryTol && Math.abs(trueConv.getEndSite() - conv.getEndSite()) / trueConv.getSiteCount() <= boundaryTol) // && Math.abs(trueConv.getHeight1()-conv.getHeight1())/trueConv.getHeight1() <= ageTol
                // && Math.abs(trueConv.getHeight2()-conv.getHeight2())/trueConv.getHeight2() <= ageTol
                {
                    count += 1;
                    convHist.put(trueConv, convHist.get(trueConv) + 1);
                    break;
                }
            }
        }
    }
    return count;
}
Also used : Conversion(bacter.Conversion)

Example 13 with Conversion

use of bacter.Conversion in project bacter by tgvaughan.

the class DifferenceFromTrueACG method countRecoveredTrueConvs.

public static int countRecoveredTrueConvs(Map<Conversion, Integer> convHist, int nSampledACGs, double minimumSupport) {
    int count = 0;
    int sampleThresh = (int) Math.round(nSampledACGs * minimumSupport);
    for (Conversion conv : convHist.keySet()) {
        if (convHist.get(conv) >= sampleThresh)
            count += 1;
    }
    return count;
}
Also used : Conversion(bacter.Conversion)

Example 14 with Conversion

use of bacter.Conversion in project bacter by tgvaughan.

the class ConvertedRegionSwap method proposal.

@Override
public double proposal() {
    if (acg.getTotalConvCount() < 2)
        return Double.NEGATIVE_INFINITY;
    // Select a random pair of recombinations
    Conversion conv1 = chooseConversion();
    Conversion conv2;
    do {
        conv2 = chooseConversion();
    } while (conv2 == conv1);
    // Switch edges corresponding to recombinations.  (Can't switch
    // loci, as this would break ordering constraint.)
    double tmpHeight1 = conv1.getHeight1();
    double tmpHeight2 = conv1.getHeight2();
    Node tmpNode1 = conv1.getNode1();
    Node tmpNode2 = conv1.getNode2();
    conv1.setHeight1(conv2.getHeight1());
    conv1.setHeight2(conv2.getHeight2());
    conv1.setNode1(conv2.getNode1());
    conv1.setNode2(conv2.getNode2());
    conv2.setHeight1(tmpHeight1);
    conv2.setHeight2(tmpHeight2);
    conv2.setNode1(tmpNode1);
    conv2.setNode2(tmpNode2);
    return 0.0;
}
Also used : Node(beast.evolution.tree.Node) Conversion(bacter.Conversion)

Example 15 with Conversion

use of bacter.Conversion in project bacter by tgvaughan.

the class MergeSplitConversion method splitProposal.

/**
 * Perform split component of merge/split move.
 *
 * @param locus locus on which to apply move
 * @return log of move HR
 */
private double splitProposal(Locus locus) {
    double logHGF = 0.0;
    if (acg.getConvCount(locus) == 0)
        return Double.NEGATIVE_INFINITY;
    Conversion conv1 = acg.getConversions(locus).get(Randomizer.nextInt(acg.getConvCount(locus)));
    logHGF -= Math.log(1.0 / acg.getConvCount(locus));
    Conversion conv2 = new Conversion();
    conv2.setLocus(locus);
    conv2.setNode1(conv1.getNode1());
    conv2.setNode2(conv1.getNode2());
    int s1, s2, e1, e2;
    int m1 = conv1.getStartSite() + Randomizer.nextInt(conv1.getSiteCount());
    int m2 = conv1.getStartSite() + Randomizer.nextInt(conv1.getSiteCount());
    if (Randomizer.nextBoolean()) {
        s1 = conv1.getStartSite();
        s2 = m1;
    } else {
        s1 = m1;
        s2 = conv1.getStartSite();
    }
    if (Randomizer.nextBoolean()) {
        e1 = conv1.getEndSite();
        e2 = m2;
    } else {
        e1 = m2;
        e2 = conv1.getEndSite();
    }
    if (e1 < s1 || e2 < s2)
        return Double.NEGATIVE_INFINITY;
    logHGF -= 2.0 * Math.log(0.5 / (conv1.getSiteCount()));
    conv1.setStartSite(s1);
    conv1.setEndSite(e1);
    conv2.setStartSite(s2);
    conv2.setEndSite(e2);
    conv2.setHeight1(conv1.getNode1().getHeight() + Randomizer.nextDouble() * conv1.getNode1().getLength());
    logHGF -= Math.log(1.0 / conv1.getNode1().getLength());
    if (conv1.getNode2().isRoot()) {
        double lambda = 1.0 / (conv1.getHeight2() - conv1.getNode2().getHeight());
        conv2.setHeight1(conv1.getNode2().getHeight() + Randomizer.nextExponential(lambda));
        logHGF -= -lambda * (conv2.getHeight2() - conv2.getNode2().getHeight()) + Math.log(lambda);
    } else {
        conv2.setHeight2(conv1.getNode2().getHeight() + Randomizer.nextDouble() * conv1.getNode2().getLength());
        logHGF -= Math.log(1.0 / conv1.getNode2().getLength());
    }
    if (conv2.getHeight2() < conv2.getHeight1())
        return Double.NEGATIVE_INFINITY;
    acg.addConversion(conv2);
    logHGF += Math.log(1.0 / (acg.getConvCount(locus) * (acg.getConvCount(locus) - 1)));
    return logHGF;
}
Also used : Conversion(bacter.Conversion)

Aggregations

Conversion (bacter.Conversion)49 Locus (bacter.Locus)27 Node (beast.evolution.tree.Node)24 ConversionGraph (bacter.ConversionGraph)7 ArrayList (java.util.ArrayList)7 Test (org.junit.Test)6 RealParameter (beast.core.parameter.RealParameter)4 SiteModel (beast.evolution.sitemodel.SiteModel)4 JukesCantor (beast.evolution.substitutionmodel.JukesCantor)4 TaxonSet (beast.evolution.alignment.TaxonSet)3 ConstantPopulation (beast.evolution.tree.coalescent.ConstantPopulation)3 ClusterTree (beast.util.ClusterTree)3 PrintStream (java.io.PrintStream)3 SimulatedACG (bacter.model.SimulatedACG)2 CFEventList (bacter.CFEventList)1 ACGCoalescent (bacter.model.ACGCoalescent)1 ACGLogReader (bacter.util.ACGLogReader)1 BacterACGLogReader (bacter.util.BacterACGLogReader)1 COACGLogFileReader (bacter.util.COACGLogFileReader)1 State (beast.core.State)1