private final double calculateLogProbabilityDensity(TreeColouring colouring,Tree tree,NodeRef node,ColourChangeMatrix mm,double[] N){
  double p=1.0;
  if (tree.isRoot(node)) {
    p=mm.getEquilibrium(colouring.getNodeColour(node));
  }
 else {
    NodeRef parent=tree.getParent(node);
    BranchColouring history=colouring.getBranchColouring(node);
    int fromColour=colouring.getNodeColour(parent);
    double fromHeight=tree.getNodeHeight(parent);
    for (int i=1; i <= history.getNumEvents(); i++) {
      int toColour=history.getForwardColourBelow(i);
      double toHeight=history.getForwardTime(i);
      p*=Math.exp(-(fromHeight - toHeight) * (-mm.getForwardRate(fromColour,fromColour)));
      p*=mm.getForwardRate(fromColour,toColour);
      fromHeight=toHeight;
      fromColour=toColour;
    }
    double toHeight=tree.getNodeHeight(node);
    p*=Math.exp(-(fromHeight - toHeight) * (-mm.getForwardRate(fromColour,fromColour)));
  }
  double logP=Math.log(p);
  for (int i=0; i < tree.getChildCount(node); i++) {
    logP+=calculateLogProbabilityDensity(colouring,tree,tree.getChild(node,i),mm,N);
  }
  return logP;
}

private final double calculateLogProbabilityDensity ( TreeColouring colouring , Tree tree , NodeRef node , ColourChangeMatrix mm , double [ ] N ) { double p = 1.0 ; if ( tree . isRoot ( node ) ) { p = mm . getEquilibrium ( colouring . getNodeColour ( node ) ) ; } else { NodeRef parent = tree . getParent ( node ) ; BranchColouring history = colouring . getBranchColouring ( node ) ; int fromColour = colouring . getNodeColour ( parent ) ; double fromHeight = tree . getNodeHeight ( parent ) ; for ( int i = 1 ; i <= history . getNumEvents ( ) ; i ++ ) { int toColour = history . getForwardColourBelow ( i ) ; double toHeight = history . getForwardTime ( i ) ; p *= Math . exp ( - ( fromHeight - toHeight ) * ( - mm . getForwardRate ( fromColour , fromColour ) ) ) ; p *= mm . getForwardRate ( fromColour , toColour ) ; fromHeight = toHeight ; fromColour = toColour ; } double toHeight = tree . getNodeHeight ( node ) ; p *= Math . exp ( - ( fromHeight - toHeight ) * ( - mm . getForwardRate ( fromColour , fromColour ) ) ) ; } double logP = Math . log ( p ) ; for ( int i = 0 ; i < tree . getChildCount ( node ) ; i ++ ) { logP += calculateLogProbabilityDensity ( colouring , tree , tree . getChild ( node , i ) , mm , N ) ; } return logP ; }
calculates log probability density of the proposal colouring of the tree on the branch leading to this node , and everything descending from it .
