Examples with Node beast.evolution.tree.Node used on opensource projects

Example 76 with Node

use of beast.evolution.tree.Node in project bacter by tgvaughan.

the class SimulatedACG method simulateClonalFrame.

/**
 * Use coalescent model to simulate clonal frame.
 */
private void simulateClonalFrame() {
    // Initialize leaf nodes
    List<Node> leafNodes = new ArrayList<>();
    for (int i = 0; i < m_taxonset.get().getTaxonCount(); i++) {
        Node leaf = new Node();
        leaf.setNr(i);
        leaf.setID(m_taxonset.get().getTaxonId(i));
        if (hasDateTrait())
            leaf.setHeight(getDateTrait().getValue(leaf.getID()));
        else
            leaf.setHeight(0.0);
        leafNodes.add(leaf);
    }
    // Create and sort list of inactive nodes
    List<Node> inactiveNodes = new ArrayList<>(leafNodes);
    Collections.sort(inactiveNodes, (Node n1, Node n2) -> {
        if (n1.getHeight() < n2.getHeight())
            return -1;
        if (n1.getHeight() > n2.getHeight())
            return 1;
        return 0;
    });
    List<Node> activeNodes = new ArrayList<>();
    double tau = 0.0;
    int nextNr = leafNodes.size();
    while (true) {
        // Calculate coalescence propensity
        int k = activeNodes.size();
        double chi = 0.5 * k * (k - 1);
        // Draw scaled coalescent time
        if (chi > 0.0)
            tau += Randomizer.nextExponential(chi);
        else
            tau = Double.POSITIVE_INFINITY;
        // Convert to real time
        double t = popFunc.getInverseIntensity(tau);
        // If new time takes us past next sample time, insert that sample
        if (!inactiveNodes.isEmpty() && t > inactiveNodes.get(0).getHeight()) {
            Node nextActive = inactiveNodes.remove(0);
            activeNodes.add(nextActive);
            tau = popFunc.getIntensity(nextActive.getHeight());
            continue;
        }
        // Coalesce random pair of active nodes.
        Node node1 = activeNodes.remove(Randomizer.nextInt(k));
        Node node2 = activeNodes.remove(Randomizer.nextInt(k - 1));
        Node parent = new Node();
        parent.addChild(node1);
        parent.addChild(node2);
        parent.setHeight(t);
        parent.setNr(nextNr++);
        activeNodes.add(parent);
        if (inactiveNodes.isEmpty() && activeNodes.size() < 2)
            break;
    }
    // Remaining active node is root
    setRoot(activeNodes.get(0));
}

Example 77 with Node

use of beast.evolution.tree.Node in project bacter by tgvaughan.

the class SimulatedAlignment method simulate.

/**
 * Perform actual sequence simulation.
 */
private void simulate(Locus locus) {
    Node cfRoot = acg.getRoot();
    int nTaxa = acg.getLeafNodeCount();
    double[] categoryProbs = siteModel.getCategoryProportions(cfRoot);
    int nCategories = siteModel.getCategoryCount();
    int nStates = dataType.getStateCount();
    double[][] transitionProbs = new double[nCategories][nStates * nStates];
    int[][] alignment = new int[nTaxa][locus.getSiteCount()];
    for (Region region : acg.getRegions(locus)) {
        int thisLength = region.getRegionLength();
        MarginalTree marginalTree = new MarginalTree(acg, region);
        int[] categories = new int[thisLength];
        for (int i = 0; i < thisLength; i++) categories[i] = Randomizer.randomChoicePDF(categoryProbs);
        int[][] regionAlignment = new int[nTaxa][thisLength];
        Node thisRoot = marginalTree.getRoot();
        int[] parentSequence = new int[region.getRegionLength()];
        double[] frequencies = siteModel.getSubstitutionModel().getFrequencies();
        for (int i = 0; i < parentSequence.length; i++) parentSequence[i] = Randomizer.randomChoicePDF(frequencies);
        traverse(thisRoot, parentSequence, categories, transitionProbs, regionAlignment);
        // DEBUG: Count differences
        int segsites = 0;
        for (int s = 0; s < regionAlignment[0].length; s++) {
            int state = regionAlignment[0][s];
            for (int l = 1; l < nTaxa; l++) {
                if (state != regionAlignment[l][s]) {
                    segsites += 1;
                    break;
                }
            }
        }
        System.out.println(segsites + " segregating sites in region " + region);
        copyToAlignment(alignment, regionAlignment, region);
    }
    for (int leafIdx = 0; leafIdx < nTaxa; leafIdx++) {
        String sSeq = dataType.encodingToString(alignment[leafIdx]);
        String sTaxon = acg.getNode(leafIdx).getID();
        sequenceInput.setValue(new Sequence(sTaxon, sSeq), this);
    }
}

Example 78 with Node

use of beast.evolution.tree.Node 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);
        }
    }
}

Example 79 with Node

use of beast.evolution.tree.Node 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;
}

Example 80 with Node

use of beast.evolution.tree.Node in project bacter by tgvaughan.

the class ACGCladeSystem method applyToClades.

/**
 * Apply a function to each sub-clade.
 *
 * @param node MRCA of clade
 * @param function function to apply. Given sub-clade parent node
 *                 and bitset as arguments.
 * @return BitSet representing clade.
 */
public BitSet applyToClades(Node node, BiFunction<Node, BitSet, Void> function) {
    BitSet bits = new BitSet();
    if (node.isLeaf()) {
        bits.set(2 * getTaxonIndex(node));
    } else {
        for (Node child : node.getChildren()) bits.or(applyToClades(child, function));
    }
    function.apply(node, bits);
    return bits;
}