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

Example 71 with Node

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

the class ACGLikelihood method setPartials.

/**
 * Set leaf partials in likelihood core.
 *
 * @param lhc likelihood core object
 * @param patterns leaf state patterns
 */
protected void setPartials(LikelihoodCore lhc, Multiset<int[]> patterns) {
    for (Node node : acg.getExternalNodes()) {
        int nStates = alignment.getDataType().getStateCount();
        double[] partials = new double[patterns.elementSet().size() * nStates];
        int k = 0;
        int iTaxon = alignment.getTaxonIndex(node.getID());
        for (int[] pattern : patterns.elementSet()) {
            int code = pattern[iTaxon];
            boolean[] stateSet = alignment.getDataType().getStateSet(code);
            for (int iState = 0; iState < nStates; iState++) {
                partials[k++] = (stateSet[iState] ? 1.0 : 0.0);
            }
        }
        lhc.setNodePartials(node.getNr(), partials);
    }
}

Example 72 with Node

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

the class ACGLikelihood method preComputeCFTransitionProbs.

/**
 * Pre-compute transition probabilities for CF edges.
 */
void preComputeCFTransitionProbs() {
    if (cfTransitionProbs == null)
        cfTransitionProbs = new double[acg.getNodeCount() - 1][siteModel.getCategoryCount()][(nStates + 1) * (nStates + 1)];
    for (int ni = 0; ni < acg.getNodeCount() - 1; ni++) {
        Node node = acg.getNode(ni);
        for (int ci = 0; ci < siteModel.getCategoryCount(); ci++) {
            double jointBranchRate = siteModel.getRateForCategory(ci, node) * branchRateModel.getRateForBranch(node);
            double parentHeight = node.getParent().getHeight();
            double nodeHeight = node.getHeight();
            substitutionModel.getTransitionProbabilities(node, parentHeight, nodeHeight, jointBranchRate, cfTransitionProbs[ni][ci]);
        }
    }
}

Example 73 with Node

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

the class ACGLikelihoodApprox method getCoalescenceHeights.

/**
 * @return map from heights of coalescences to objects describing
 * the sites and samples they involve.
 */
Map<Double, Coalescence> getCoalescenceHeights() {
    Map<Double, Coalescence> heightMap = new HashMap<>();
    Map<Node, SiteAncestry> activeCFNodes = new HashMap<>();
    Map<Conversion, SiteAncestry> activeConversions = new HashMap<>();
    ACGEventList acgEventList = new ACGEventList(acg, locus);
    for (ACGEventList.Event event : acgEventList.getACGEvents()) {
        switch(event.type) {
            case CF_LEAF:
                activeCFNodes.put(event.node, new SiteAncestry(event.node, locus));
                break;
            case CF_COALESCENCE:
                Node node1 = event.node.getLeft();
                Node node2 = event.node.getRight();
                SiteAncestry ancestryCF = new SiteAncestry();
                Coalescence coalescenceCF = new Coalescence();
                activeCFNodes.get(node1).merge(activeCFNodes.get(node2), coalescenceCF, ancestryCF);
                activeCFNodes.remove(node1);
                activeCFNodes.remove(node2);
                activeCFNodes.put(event.node, ancestryCF);
                if (coalescenceCF.getIntervalCount() > 0)
                    heightMap.put(event.t, coalescenceCF);
                break;
            case CONV_DEPART:
                SiteAncestry inside = new SiteAncestry();
                SiteAncestry outside = new SiteAncestry();
                activeCFNodes.get(event.node).split(event.conversion.getStartSite(), event.conversion.getEndSite() + 1, inside, outside);
                if (inside.getIntervalCount() > 0) {
                    activeCFNodes.put(event.node, outside);
                    activeConversions.put(event.conversion, inside);
                }
                break;
            case CONV_ARRIVE:
                if (!activeConversions.containsKey(event.conversion))
                    continue;
                SiteAncestry ancestry = new SiteAncestry();
                Coalescence coalescence = new Coalescence();
                activeCFNodes.get(event.node).merge(activeConversions.get(event.conversion), coalescence, ancestry);
                activeCFNodes.put(event.node, ancestry);
                activeConversions.remove(event.conversion);
                if (coalescence.getIntervalCount() > 0)
                    heightMap.put(event.t, coalescence);
                break;
        }
    }
    return heightMap;
}

Example 74 with Node

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

the class ACGLikelihoodBeagle method setPartials.

/**
 * Set leaf partials in a Beagle instance
 *
 * @param beagle beagle instance object
 * @param patterns leaf state patterns
 */
protected void setPartials(Beagle beagle, Multiset<int[]> patterns) {
    for (Node node : acg.getExternalNodes()) {
        Alignment data = dataInput.get();
        int nStates = data.getDataType().getStateCount();
        double[] partials = new double[patterns.elementSet().size() * nStates * siteModel.getCategoryCount()];
        int k = 0;
        int iTaxon = alignment.getTaxonIndex(node.getID());
        for (int[] pattern : patterns.elementSet()) {
            int code = pattern[iTaxon];
            boolean[] stateSet = alignment.getDataType().getStateSet(code);
            for (int iState = 0; iState < nStates; iState++) {
                partials[k++] = (stateSet[iState] ? 1.0 : 0.0);
            }
        }
        int n = patterns.elementSet().size() * siteModel.getCategoryCount();
        for (int cIdx = 1; cIdx < siteModel.getCategoryCount(); cIdx++) {
            System.arraycopy(partials, 0, partials, n * cIdx, n);
        }
        beagle.setTipPartials(node.getNr(), partials);
    }
}

Example 75 with Node

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

the class SimulatedACG method associateConversionWithCF.

/**
 * Associates recombination with the clonal frame, selecting points of
 * departure and coalescence.
 *
 * @param conv recombination to associate
 */
private void associateConversionWithCF(Conversion conv) {
    List<CFEventList.Event> eventList = getCFEvents();
    // Select departure point
    double u = Randomizer.nextDouble() * getClonalFrameLength();
    boolean started = false;
    for (int eidx = 0; eidx < eventList.size(); eidx++) {
        CFEventList.Event event = eventList.get(eidx);
        if (!started) {
            double interval = eventList.get(eidx + 1).getHeight() - event.getHeight();
            if (u < interval * event.getLineageCount()) {
                for (Node node : getNodesAsArray()) {
                    if (!node.isRoot() && node.getHeight() <= event.getHeight() && node.getParent().getHeight() > event.getHeight()) {
                        if (u < interval) {
                            conv.setNode1(node);
                            conv.setHeight1(event.getHeight() + u);
                            break;
                        } else
                            u -= interval;
                    }
                }
                started = true;
                u = Randomizer.nextExponential(1.0);
            } else
                u -= interval * event.getLineageCount();
        }
        if (started) {
            double t = Math.max(event.getHeight(), conv.getHeight1());
            double intervalArea;
            if (eidx < eventList.size() - 1) {
                intervalArea = popFunc.getIntegral(t, eventList.get(eidx + 1).getHeight());
            } else
                intervalArea = Double.POSITIVE_INFINITY;
            if (u < intervalArea * event.getLineageCount()) {
                // Fix height of attachment point
                double tauEnd = popFunc.getIntensity(t) + u / event.getLineageCount();
                double tEnd = popFunc.getInverseIntensity(tauEnd);
                conv.setHeight2(tEnd);
                // Choose particular lineage to attach to
                int nodeNumber = Randomizer.nextInt(event.getLineageCount());
                for (Node node : getNodesAsArray()) {
                    if (node.getHeight() <= event.getHeight() && (node.isRoot() || node.getParent().getHeight() > event.getHeight())) {
                        if (nodeNumber == 0) {
                            conv.setNode2(node);
                            break;
                        } else
                            nodeNumber -= 1;
                    }
                }
                break;
            } else
                u -= intervalArea * event.getLineageCount();
        }
    }
}