Examples with NodeRef dr.evolution.tree.NodeRef used on opensource projects

Example 31 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class MulTreeNodeSlide method operateOneNode.

public void operateOneNode(final double factor) {
    //            #print "operate: tree", ut.treerep(t)
    //   if( verbose)  System.out.println("  Mau at start: " + tree.getSimpleTree());
    final int count = multree.getExternalNodeCount();
    assert count == species.nSpSeqs();
    NodeRef[] order = new NodeRef[2 * count - 1];
    boolean[] swapped = new boolean[count - 1];
    mauCanonical(multree, order, swapped);
    // internal node to change
    // count-1 - number of internal nodes
    int which = MathUtils.nextInt(count - 1);
    FixedBitSet left = new FixedBitSet(count);
    FixedBitSet right = new FixedBitSet(count);
    for (int k = 0; k < 2 * which + 1; k += 2) {
        left.set(multree.speciesIndex(order[k]));
    }
    for (int k = 2 * (which + 1); k < 2 * count; k += 2) {
        right.set(multree.speciesIndex(order[k]));
    }
    double newHeight;
    if (factor > 0) {
        newHeight = multree.getNodeHeight(order[2 * which + 1]) * factor;
    } else {
        final double limit = species.speciationUpperBound(left, right);
        newHeight = MathUtils.nextDouble() * limit;
    }
    multree.beginTreeEdit();
    multree.setPreorderIndices(preOrderIndexBefore);
    final NodeRef node = order[2 * which + 1];
    multree.setNodeHeight(node, newHeight);
    mauReconstruct(multree, order, swapped);
    // restore pre-order of pops -
    {
        multree.setPreorderIndices(preOrderIndexAfter);
        double[] splitPopValues = null;
        for (int k = 0; k < preOrderIndexBefore.length; ++k) {
            final int b = preOrderIndexBefore[k];
            if (b >= 0) {
                final int a = preOrderIndexAfter[k];
                if (a != b) {
                    //if( verbose)  System.out.println("pops: " + a + " <- " + b);
                    final Parameter p1 = multree.sppSplitPopulations;
                    if (splitPopValues == null) {
                        splitPopValues = p1.getParameterValues();
                    }
                    if (multree.constPopulation()) {
                        p1.setParameterValue(count + a, splitPopValues[count + b]);
                    } else {
                        for (int i = 0; i < 2; ++i) {
                            p1.setParameterValue(count + 2 * a + i, splitPopValues[count + 2 * b + i]);
                        }
                    }
                }
            }
        }
    }
    multree.endTreeEdit();
}

Example 32 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class AnyTipObservationProcess method calculateLogTreeWeight.

public double calculateLogTreeWeight() {
    int L = treeModel.getNodeCount();
    if (u0 == null || p == null) {
        // probability that the trait at node i survives to no leaf
        u0 = new double[L];
        // probability of survival on the branch ancestral to i
        p = new double[L];
    }
    int i, j, childNumber;
    NodeRef node;
    double logWeight = 0.0;
    double averageRate = getAverageRate();
    for (i = 0; i < L; ++i) {
        p[i] = 1.0 - getNodeSurvivalProbability(i, averageRate);
    }
    TreeUtils.postOrderTraversalList(treeModel, postOrderNodeList);
    for (int postOrderIndex = 0; postOrderIndex < nodeCount; postOrderIndex++) {
        i = postOrderNodeList[postOrderIndex];
        if (i < treeModel.getExternalNodeCount()) {
            // Is tip
            u0[i] = 0.0;
            logWeight += 1.0 - p[i];
        } else {
            // Is internal node or root
            u0[i] = 1.0;
            node = treeModel.getNode(i);
            for (j = 0; j < treeModel.getChildCount(node); ++j) {
                childNumber = treeModel.getChild(node, j).getNumber();
                u0[i] *= 1.0 - p[childNumber] * (1.0 - u0[childNumber]);
            }
            logWeight += (1.0 - u0[i]) * (1.0 - p[i]);
        }
    }
    return -logWeight * lam.getParameterValue(0) / (getAverageRate() * mu.getParameterValue(0));
}

Example 33 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class DiscretizedLocationOperator method makeLocationList.

private Set<Point2D> makeLocationList() {
    Set<Point2D> uniquePoints = new HashSet<Point2D>();
    for (int i = 0; i < treeModel.getExternalNodeCount(); i++) {
        NodeRef node = treeModel.getExternalNode(i);
        double[] leafTrait = treeModel.getMultivariateNodeTrait(node, traitName);
        Point2D.Double point = new Point2D.Double(leafTrait[0], leafTrait[1]);
        if (!uniquePoints.contains(point)) {
            uniquePoints.add(point);
            savedPt = point;
        }
    }
    return uniquePoints;
}

Example 34 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class DiscretizedLocationOperator method doOperation.

public double doOperation() {
    NodeRef node;
    if (onlyInternalNodes)
        node = treeModel.getInternalNode(MathUtils.nextInt(treeModel.getInternalNodeCount()));
    else
        node = treeModel.getNode(MathUtils.nextInt(treeModel.getNodeCount()));
    double[] trait = treeModel.getMultivariateNodeTrait(node, traitName);
    Point2D currentPt = new Point2D.Double(trait[0], trait[1]);
    List<WeightedPoint2D> neighbors = nearestNeighborMap.get(currentPt);
    if (neighbors == null)
        throw new RuntimeException("Node location outside allowable values: " + currentPt);
    //        Point2D newPt = neighbors.get(MathUtils.nextInt(disk));
    Point2D newPt = neighbors.get(MathUtils.nextInt(convertFromAutoOptimizeValue(autoOptimize)));
    trait[0] = newPt.getX();
    trait[1] = newPt.getY();
    //        treeModel.setMultivariateTrait(node, traitName, trait);
    recursivelySetTrait(node, trait, null);
    return 0;
}

Example 35 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class ExchangeOperator method narrow.

/**
     * WARNING: Assumes strictly bifurcating tree.
     */
public boolean narrow() {
    final int nNodes = tree.getNodeCount();
    final NodeRef root = tree.getRoot();
    NodeRef i = root;
    while (root == i || tree.getParent(i) == root) {
        i = tree.getNode(MathUtils.nextInt(nNodes));
    }
    final NodeRef iParent = tree.getParent(i);
    final NodeRef iGrandParent = tree.getParent(iParent);
    NodeRef iUncle = tree.getChild(iGrandParent, 0);
    if (iUncle == iParent) {
        iUncle = tree.getChild(iGrandParent, 1);
    }
    assert iUncle == getOtherChild(tree, iGrandParent, iParent);
    assert tree.getNodeHeight(i) <= tree.getNodeHeight(iGrandParent);
    if (tree.getNodeHeight(iUncle) < tree.getNodeHeight(iParent)) {
        exchangeNodes(tree, i, iUncle, iParent, iGrandParent);
        //tree.pushTreeChangedEvent(iGrandParent);
        return true;
    }
    return false;
}