Search in sources :

Example 6 with AbstractCase

use of dr.evomodel.epidemiology.casetocase.AbstractCase in project beast-mcmc by beast-dev.

the class TransmissionTreeOperator method doOperation.

public double doOperation() {
    TreeModel tree = c2cLikelihood.getTreeModel();
    BranchMapModel branchMap = c2cLikelihood.getBranchMap();
    AbstractCase[] newBranchMap = branchMap.getArrayCopy();
    int[] oldParents = getParentsArray(tree);
    double[] oldHeights = getHeightsArray(tree);
    double hr = innerOperator.doOperation();
    int[] newParents = getParentsArray(tree);
    ArrayList<Integer> changedNodes = new ArrayList<Integer>();
    for (int i = 0; i < tree.getNodeCount(); i++) {
        if (oldParents[i] != newParents[i]) {
            changedNodes.add(i);
        }
    }
    if (changedNodes.size() != 0) {
        if (innerOperator instanceof ExchangeOperator) {
            //this is a node swap operator
            AbstractCase[] nodePaintings = new AbstractCase[2];
            AbstractCase[] parentPaintings = new AbstractCase[2];
            for (int i = 0; i < 2; i++) {
                nodePaintings[i] = branchMap.get(changedNodes.get(i));
                parentPaintings[i] = branchMap.get(oldParents[changedNodes.get(i)]);
            }
            if (nodePaintings[0] == parentPaintings[0] || nodePaintings[1] == parentPaintings[1]) {
                //If this is not true there is nothing to do - the result is already a valid painting
                for (int i = 0; i < 2; i++) {
                    paintUp(tree, nodePaintings[1 - i], nodePaintings[i], branchMap, newBranchMap, tree.getNode(changedNodes.get(i)).getNumber(), newParents);
                }
            }
        } else if (innerOperator instanceof SubtreeSlideOperator || innerOperator instanceof WilsonBalding) {
            //this is a node transplantation operator
            int movedNode = -1;
            int oldChild = -1;
            int newChild = -1;
            for (int i = 0; i < tree.getNodeCount(); i++) {
                // exception later, though.
                if (tree.getNodeHeight(tree.getNode(i)) != oldHeights[i]) {
                    movedNode = i;
                    break;
                }
            }
            for (int j : changedNodes) {
                if (j != movedNode) {
                    if (tree.getParent(tree.getNode(j)) == tree.getNode(movedNode)) {
                        newChild = j;
                    } else {
                        oldChild = j;
                    }
                }
            }
            if (movedNode == -1 || oldChild == -1 || newChild == -1) {
                // is this a bug or should the move be rejected (i.e., return -Inf HR)?
                throw new RuntimeException("Failed to establish relationship between relocated node and others");
            }
            NodeRef movedNodeObject = tree.getNode(movedNode);
            NodeRef oldChildObject = tree.getNode(oldChild);
            NodeRef newChildObject = tree.getNode(newChild);
            int otherChild = -1;
            NodeRef otherChildObject;
            //Find the other child of the moved node (the root of the transplanted subtree)
            for (int i = 0; i < tree.getChildCount(movedNodeObject); i++) {
                if (tree.getChild(movedNodeObject, i) != newChildObject) {
                    otherChildObject = tree.getChild(movedNodeObject, i);
                    otherChild = otherChildObject.getNumber();
                }
            }
            //If the child of the moved node is the earliest node with its painting:
            if (branchMap.get(otherChild) != branchMap.get(movedNode)) {
                newBranchMap[movedNode] = branchMap.get(newChild);
            } else {
                //Change all paintings up the tree from the old child that used to match the moved node to match
                //the old child
                paintUp(tree, branchMap.get(movedNode), branchMap.get(oldChild), branchMap, newBranchMap, oldChild, oldParents);
                //This may have resulted in the moved node being recoloured wrong.
                newBranchMap[movedNode] = branchMap.get(movedNode);
                branchMap.setAll(newBranchMap, true);
                //Change all paintings up the tree from the moved node that used to match the new child to match
                //the moved node
                paintUp(tree, branchMap.get(newChild), branchMap.get(movedNode), branchMap, newBranchMap, movedNode, newParents);
            }
        } else {
            //I don't know what this is
            throw new UnsupportedOperationException("Operator class " + innerOperator.getOperatorName() + " not yet " + "supported");
        }
    }
    branchMap.setAll(newBranchMap, true);
    c2cLikelihood.makeDirty();
    return hr;
}
Also used : ArrayList(java.util.ArrayList) ExchangeOperator(dr.evomodel.operators.ExchangeOperator) SubtreeSlideOperator(dr.evomodel.operators.SubtreeSlideOperator) AbstractCase(dr.evomodel.epidemiology.casetocase.AbstractCase) TreeModel(dr.evomodel.tree.TreeModel) NodeRef(dr.evolution.tree.NodeRef) BranchMapModel(dr.evomodel.epidemiology.casetocase.BranchMapModel) WilsonBalding(dr.evomodel.operators.WilsonBalding)

Example 7 with AbstractCase

use of dr.evomodel.epidemiology.casetocase.AbstractCase in project beast-mcmc by beast-dev.

the class TransmissionExchangeOperatorA method exchange.

public double exchange() {
    TreeModel tree = c2cLikelihood.getTreeModel();
    final int nodeCount = tree.getNodeCount();
    final NodeRef root = tree.getRoot();
    NodeRef i = root;
    while (root == i) {
        i = tree.getNode(MathUtils.nextInt(nodeCount));
    }
    ArrayList<NodeRef> candidates = getPossibleExchanges(tree, i);
    int candidateCount = candidates.size();
    if (candidateCount == 0) {
        //            throw new OperatorFailedException("No valid exchanges for this node");
        return Double.NEGATIVE_INFINITY;
    }
    NodeRef j = candidates.get(MathUtils.nextInt(candidates.size()));
    int jFirstCandidateCount = getPossibleExchanges(tree, j).size();
    double HRDenom = (1 / ((double) candidateCount)) + (1 / ((double) jFirstCandidateCount));
    NodeRef iP = tree.getParent(i);
    NodeRef jP = tree.getParent(j);
    if (resampleInfectionTimes) {
        BranchMapModel branchMap = c2cLikelihood.getBranchMap();
        AbstractCase iCase = branchMap.get(i.getNumber());
        AbstractCase jCase = branchMap.get(j.getNumber());
        AbstractCase parentCase = branchMap.get(iP.getNumber());
        if (iCase != parentCase) {
            iCase.setInfectionBranchPosition(MathUtils.nextDouble());
        }
        if (jCase != parentCase) {
            jCase.setInfectionBranchPosition(MathUtils.nextDouble());
        }
    }
    exchangeNodes(tree, i, j, iP, jP);
    ArrayList<NodeRef> reverseCandidatesIfirst = getPossibleExchanges(tree, i);
    ArrayList<NodeRef> reverseCandidatesJfirst = getPossibleExchanges(tree, j);
    double HRNum = (1 / (double) reverseCandidatesIfirst.size()) + (1 / (double) reverseCandidatesJfirst.size());
    return Math.log(HRNum / HRDenom);
}
Also used : TreeModel(dr.evomodel.tree.TreeModel) NodeRef(dr.evolution.tree.NodeRef) BranchMapModel(dr.evomodel.epidemiology.casetocase.BranchMapModel) AbstractCase(dr.evomodel.epidemiology.casetocase.AbstractCase)

Aggregations

NodeRef (dr.evolution.tree.NodeRef)7 AbstractCase (dr.evomodel.epidemiology.casetocase.AbstractCase)7 BranchMapModel (dr.evomodel.epidemiology.casetocase.BranchMapModel)7 TreeModel (dr.evomodel.tree.TreeModel)4 ArrayList (java.util.ArrayList)3 PartitionedTreeModel (dr.evomodel.epidemiology.casetocase.PartitionedTreeModel)1 ExchangeOperator (dr.evomodel.operators.ExchangeOperator)1 SubtreeSlideOperator (dr.evomodel.operators.SubtreeSlideOperator)1 WilsonBalding (dr.evomodel.operators.WilsonBalding)1 HashSet (java.util.HashSet)1