Search in sources :

Example 21 with Node

use of dr.evomodel.arg.ARGModel.Node in project beast-mcmc by beast-dev.

the class ARGTree method toGraphString.

public String toGraphString() {
    StringBuffer sb = new StringBuffer();
    for (Node node : nodes) {
        sb.append(node.number);
        if (node.leftParent != null)
            sb.append(" " + node.leftParent.number);
        else
            sb.append(" 0");
        if (node.rightParent != null)
            sb.append(" " + node.rightParent.number);
        else
            sb.append(" 0");
        if (node.leftChild != null)
            sb.append(" " + node.leftChild.number);
        else
            sb.append(" 0");
        if (node.rightChild != null)
            sb.append(" " + node.rightChild.number);
        else
            sb.append(" 0");
        if (node.taxon != null)
            sb.append(" " + node.taxon.toString());
        sb.append("\n");
    }
    sb.append("Root = " + ((Node) getRoot()).number + "\n");
    return new String(sb);
}
Also used : Node(dr.evomodel.arg.ARGModel.Node)

Example 22 with Node

use of dr.evomodel.arg.ARGModel.Node in project beast-mcmc by beast-dev.

the class ObsoleteARGNewEventOperator method AddOperation.

private double AddOperation() throws ARGOperatorFailedException {
    double logq = 0;
    // Draw attachment point for new bifurcation
    //	    System.err.println("Starting add operation.");
    ArrayList<NodeRef> potentialBifurcationChildren = new ArrayList<NodeRef>();
    ArrayList<NodeRef> potentialReassortmentChildren = new ArrayList<NodeRef>();
    double treeHeight = arg.getNodeHeight(arg.getRoot());
    double newBifurcationHeight = treeHeight * MathUtils.nextDouble();
    double newReassortmentHeight = treeHeight * MathUtils.nextDouble();
    if (newReassortmentHeight > newBifurcationHeight) {
        double temp = newReassortmentHeight;
        newReassortmentHeight = newBifurcationHeight;
        newBifurcationHeight = temp;
    }
    //	    logq += 2.0 * Math.log(treeHeight); // Uniform before sorting
    int totalPotentialBifurcationChildren = findPotentialAttachmentPoints(newBifurcationHeight, potentialBifurcationChildren);
    int totalPotentialReassortmentChildren = findPotentialAttachmentPoints(newReassortmentHeight, potentialReassortmentChildren);
    if (totalPotentialBifurcationChildren == 0 || totalPotentialReassortmentChildren == 0) {
        throw new RuntimeException("Unable to find attachment points.");
    }
    Node recNode = (Node) potentialReassortmentChildren.get(MathUtils.nextInt(totalPotentialReassortmentChildren));
    //	    logq += Math.log(totalPotentialReassortmentChildren);
    Node recParentL = recNode.leftParent;
    Node recParentR = recNode.rightParent;
    Node recParent = recParentL;
    if (recParentL != recParentR) {
        if (arg.getNodeHeight(recParentL) < newReassortmentHeight)
            recParent = recParentR;
        else if (arg.getNodeHeight(recParentR) > newReassortmentHeight && MathUtils.nextDouble() > 0.5)
            recParent = recParentR;
    }
    Node sisNode = (Node) potentialBifurcationChildren.get(MathUtils.nextInt(totalPotentialBifurcationChildren));
    //	    logq += Math.log(totalPotentialBifurcationChildren);
    Node sisParentL = sisNode.leftParent;
    Node sisParentR = sisNode.rightParent;
    Node sisParent = sisParentL;
    if (sisParentL != sisParentR) {
        if (arg.getNodeHeight(sisParentL) < newBifurcationHeight)
            sisParent = sisParentR;
        else if (arg.getNodeHeight(sisParentR) > newBifurcationHeight && MathUtils.nextDouble() > 0.5)
            sisParent = sisParentR;
    }
    double newBifurcationRateCategory = 1.0;
    double newReassortmentRateCategory = 1.0;
    Node newBifurcation = arg.new Node();
    newBifurcation.heightParameter = new Parameter.Default(newBifurcationHeight);
    newBifurcation.rateParameter = new Parameter.Default(newBifurcationRateCategory);
    newBifurcation.setupHeightBounds();
    Node newReassortment = arg.new Node();
    newReassortment.bifurcation = false;
    newReassortment.heightParameter = new Parameter.Default(newReassortmentHeight);
    newReassortment.rateParameter = new Parameter.Default(newReassortmentRateCategory);
    newReassortment.setupHeightBounds();
    arg.beginTreeEdit();
    if (sisParent.bifurcation)
        arg.singleRemoveChild(sisParent, sisNode);
    else
        arg.doubleRemoveChild(sisParent, sisNode);
    if (sisNode != recNode) {
        if (recParent.bifurcation)
            arg.singleRemoveChild(recParent, recNode);
        else
            arg.doubleRemoveChild(recParent, recNode);
    }
    if (sisParent.bifurcation)
        arg.singleAddChild(sisParent, newBifurcation);
    else
        arg.doubleAddChild(sisParent, newBifurcation);
    if (sisNode != recNode)
        arg.singleAddChild(newBifurcation, sisNode);
    arg.doubleAddChild(newReassortment, recNode);
    Parameter partitioning = new Parameter.Default(arg.getNumberOfPartitions());
    //        logq +=
    drawRandomPartitioning(partitioning);
    if (sisNode != recNode) {
        arg.addChildAsRecombinant(newBifurcation, recParent, newReassortment, partitioning);
    } else {
        arg.addChildAsRecombinant(newBifurcation, newBifurcation, newReassortment, partitioning);
    }
    arg.expandARGWithRecombinant(newBifurcation, newReassortment, internalNodeParameters, internalAndRootNodeParameters, nodeRates);
    //		System.err.println("point 2");
    arg.pushTreeSizeChangedEvent();
    arg.endTreeEdit();
    //        try {
    //            arg.checkTreeIsValid();
    //        } catch (MutableTree.InvalidTreeException ite) {
    //            throw new RuntimeException(ite.toString() + "\n" + arg.toString()
    //                    + "\n" + Tree.Utils.uniqueNewick(arg, arg.getRoot()));
    //        }
    //	    double d1 = findPotentialNodesToRemove(null);
    //	    System.err.printf("d1 = %5.4f\n",d1);
    //        logq -= Math.log(findPotentialNodesToRemove(null));
    //        if (!(recNode.leftParent.isBifurcation() && recNode.rightParent.isRoot()) &&
    //                !(recNode.rightParent.isBifurcation() && recNode.leftParent.isRoot()))
    //            logq -= Math.log(2.0);
    //	    System.err.println("End add ARG operation.");
    //        int nodes = arg.getInternalNodeCount() - 1;
    //
    /*	    int i = findPotentialNodesToRemove(null);
	    if (i==0) {
		    System.err.println("why can't i remove this one?");
		    System.err.println("graph:"+arg.toGraphString());
		    System.exit(1);
		    return Double.NEGATIVE_INFINITY ;
	    }*/
    logq = 0;
    //	    System.err.println("After add:\n"+arg.toGraphString());
    return logq;
}
Also used : NodeRef(dr.evolution.tree.NodeRef) Node(dr.evomodel.arg.ARGModel.Node) ArrayList(java.util.ArrayList) CompoundParameter(dr.inference.model.CompoundParameter) Parameter(dr.inference.model.Parameter)

Example 23 with Node

use of dr.evomodel.arg.ARGModel.Node in project beast-mcmc by beast-dev.

the class ObsoleteARGNewEventOperator method findPotentialNodesToRemove.

private int findPotentialNodesToRemove(ArrayList<NodeRef> list) {
    int count = 0;
    int n = arg.getNodeCount();
    //	    int max = arg.getReassortmentNodeCount();
    Node root = (Node) arg.getRoot();
    for (int i = 0; i < n; i++) {
        Node node = (Node) arg.getNode(i);
        if (node.isReassortment() && ((node.leftParent != root && node.leftParent.isBifurcation()) || (node.rightParent != root && node.rightParent.isBifurcation()))) //	                && !(
        //		                    (node.leftParent == root && node.rightParent.isReassortment()) ||
        //				            (node.rightParent== root && node.leftParent.isReassortment()))
        //		            && !(node.leftParent == root && node.rightParent == root)
        //		            && (node.leftParent.isBifurcation() || node.rightParent.isBifurcation())
        {
            if (list != null)
                list.add(node);
            count++;
        }
    }
    //	    System.err.printf("max = %d, available = %d\n",max,count);
    return count;
}
Also used : Node(dr.evomodel.arg.ARGModel.Node)

Example 24 with Node

use of dr.evomodel.arg.ARGModel.Node in project beast-mcmc by beast-dev.

the class ObsoleteARGNewEventOperator method sanityCheck.

//	if( !recParent1.bifurcation || recParent1.isRoot() ) { // One orientation is valid
public void sanityCheck() {
    int len = arg.getNodeCount();
    for (int i = 0; i < len; i++) {
        Node node = (Node) arg.getNode(i);
        if (node.bifurcation) {
            boolean equalChild = (node.leftChild == node.rightChild);
            if ((equalChild && node.leftChild != null)) {
                if (!node.leftChild.bifurcation && ((node.leftChild).leftParent == node))
                    ;
                else {
                    System.err.println("Node " + (i + 1) + " is insane.");
                    System.err.println(arg.toGraphString());
                    System.exit(-1);
                }
            }
        } else {
            if ((node.leftChild != node.rightChild)) {
                System.err.println("Node " + (i + 1) + " is insane.");
                System.err.println(arg.toGraphString());
                System.exit(-1);
            }
        }
        if (!node.isRoot()) {
            double d;
            d = node.getHeight();
        }
    }
}
Also used : Node(dr.evomodel.arg.ARGModel.Node)

Example 25 with Node

use of dr.evomodel.arg.ARGModel.Node in project beast-mcmc by beast-dev.

the class ARGSwapOperator method doNarrowSwap.

private void doNarrowSwap(NarrowSwap swap) {
    arg.beginTreeEdit();
    String before = arg.toARGSummary();
    if (swap.c == swap.pb) {
        Node c = (Node) swap.c;
        Node p = (Node) swap.p;
        Node gp = (Node) swap.gp;
        if (c.leftParent == p) {
            c.leftParent = gp;
            c.rightParent = p;
        } else {
            c.leftParent = p;
            c.rightParent = gp;
        }
    } else if (arg.getChild(swap.p, 0) == arg.getChild(swap.p, 1)) {
        Node p = (Node) swap.p;
        Node c = (Node) swap.c;
        if (MathUtils.nextBoolean())
            p.leftChild = c.leftParent = null;
        else
            p.rightChild = c.rightParent = null;
        arg.removeChild(swap.gp, swap.pb);
        arg.singleAddChild(swap.gp, swap.c);
        arg.singleAddChild(swap.p, swap.pb);
    } else {
        arg.removeChild(swap.gp, swap.pb);
        arg.removeChild(swap.p, swap.c);
        arg.singleAddChild(swap.gp, swap.c);
        arg.singleAddChild(swap.p, swap.pb);
    }
    assert nodeCheck() : swap + " " + before + " " + arg.toARGSummary();
    arg.pushTreeChangedEvent(swap.gp);
    arg.pushTreeChangedEvent(swap.p);
    arg.endTreeEdit();
    try {
        arg.checkTreeIsValid();
    } catch (MutableTree.InvalidTreeException ite) {
        System.out.println(swap);
        System.out.println(before);
        System.err.println(ite.getMessage());
        System.exit(-1);
    } catch (NullPointerException e) {
        System.out.println(swap);
        System.out.println(before);
        System.err.println(e.getMessage());
        System.exit(-1);
    }
}
Also used : MutableTree(dr.evolution.tree.MutableTree) Node(dr.evomodel.arg.ARGModel.Node)

Aggregations

Node (dr.evomodel.arg.ARGModel.Node)25 NodeRef (dr.evolution.tree.NodeRef)9 ArrayList (java.util.ArrayList)8 MutableTree (dr.evolution.tree.MutableTree)5 CompoundParameter (dr.inference.model.CompoundParameter)3 Parameter (dr.inference.model.Parameter)3