[beast-mcmc] r6754 committed - Hastings ratio implemented

[beast...@googlecode.com]

Revision: 6754
Author: rambaut
Date: Fri Apr 24 15:19:15 2015 UTC
Log: Hastings ratio implemented
https://code.google.com/p/beast-mcmc/source/detail?r=6754

Modified:
/trunk/src/dr/evomodel/operators/SubtreeJumpOperator.java
/trunk/src/dr/evomodelxml/operators/SubtreeJumpOperatorParser.java

=======================================
--- /trunk/src/dr/evomodel/operators/SubtreeJumpOperator.java Fri Apr 24
14:44:00 2015 UTC
+++ /trunk/src/dr/evomodel/operators/SubtreeJumpOperator.java Fri Apr 24
15:19:15 2015 UTC
@@ -127,7 +127,9 @@

tree.endTreeEdit();

- double reverseProbability = getReverseProbability(tree, j, i,
height, destinations);
+ final List<NodeRef> reverseDestinations =
getIntersectingEdges(tree, height);
+
+ double reverseProbability = getReverseProbability(tree, CiP, j,
height, reverseDestinations);

logq = Math.log(forwardProbability) - Math.log(reverseProbability);

@@ -147,6 +149,8 @@
for (int i = 0; i < tree.getNodeCount(); i++) {
final NodeRef node = tree.getNode(i);
final NodeRef parent = tree.getParent(node);
+
+ // The original node and its sibling will not be included
because their height is exactly equal to the target height
if (parent != null && tree.getNodeHeight(node) < height &&
tree.getNodeHeight(parent) > height) {
intersectingEdges.add(node);
}
@@ -159,13 +163,11 @@
double sum = 0.0;
int i = 0;
for (NodeRef node1 : intersectingEdges) {
- if (node1 != node0) {
- double age =
tree.getNodeHeight(Tree.Utils.getCommonAncestor(tree, node0, node1)) -
height;
- weights[i] = 1.0 / Math.pow(age, size);
- sum += weights[i];
- } else {
- weights[i] = 0.0;
- }
+ assert(node1 != node0);
+
+ double age =
tree.getNodeHeight(Tree.Utils.getCommonAncestor(tree, node0, node1)) -
height;
+ weights[i] = 1.0 / Math.pow(age, 1.0 / size);
+ sum += weights[i];
i++;
}
for (int j = 0; j < weights.length; j++) {
@@ -183,7 +185,7 @@
for (NodeRef node1 : intersectingEdges) {
if (node1 != targetNode) {
double age =
tree.getNodeHeight(Tree.Utils.getCommonAncestor(tree, targetNode, node1)) -
height;
- weights[i] = 1.0 / Math.pow(age, size);
+ weights[i] = 1.0 / Math.pow(age, 1.0 / size);
sum += weights[i];

if (node1 == originalNode) {
=======================================
--- /trunk/src/dr/evomodelxml/operators/SubtreeJumpOperatorParser.java Thu
Apr 23 18:51:55 2015 UTC
+++ /trunk/src/dr/evomodelxml/operators/SubtreeJumpOperatorParser.java Fri
Apr 24 15:19:15 2015 UTC
@@ -41,7 +41,7 @@
" for tree " + treeModel.getId() );
}

- SubtreeJumpOperator operator = new SubtreeJumpOperator(treeModel,
size, weight, mode);
+ SubtreeJumpOperator operator = new SubtreeJumpOperator(treeModel,
weight, size, mode);
// operator.setTargetAcceptanceProbability(targetAcceptance);

return operator;