]> git.donarmstrong.com Git - ape.git/blobdiff - R/compar.ou.R
a collection of bug fixes
[ape.git] / R / compar.ou.R
index 83e0d1cf5ab5eb0d939caac0a370a9884bc153d6..ada36d72cc0ffed14934ae9677666ae12ddd2c75 100644 (file)
@@ -1,15 +1,15 @@
-## compar.ou.R (2006-10-05)
+## compar.ou.R (2009-05-10)
 
 ##   Ornstein--Uhlenbeck Model for Continuous Characters
 
-## Copyright 2005-2006 Emmanuel Paradis
+## Copyright 2005-2009 Emmanuel Paradis
 
 ## This file is part of the R-package `ape'.
 ## See the file ../COPYING for licensing issues.
 
 compar.ou <- function(x, phy, node = NULL, alpha = NULL)
 {
-    if (class(phy) != "phylo")
+    if (!inherits(phy, "phylo"))
       stop('object "phy" is not of class "phylo".')
     if (!is.numeric(x)) stop("'x' must be numeric.")
     if (!is.null(names(x))) {
@@ -40,7 +40,9 @@ compar.ou <- function(x, phy, node = NULL, alpha = NULL)
     }
     W <- cophenetic.phylo(phy)
     dev <- function(p) {
-        M <- rowSums(exp(-p[1] * Wstart) - exp(-p[1] * Wend) * p[-(1:2)])
+        ##M <- rowSums(exp(-p[1] * Wstart) - exp(-p[1] * Wend) * p[-(1:2)])
+        ##M <- -rowSums(exp(-p[1] * Wstart) - exp(-p[1] * Wend) * p[-(1:2)])
+        M <- rowSums((exp(-p[1] * Wend) - exp(-p[1] * Wstart)) * p[-(1:2)])
         V <- exp(-p[1]*W) * (1 - exp(-2*p[1]*(Tmax - W/2)))
         nb.tip*log(2*pi*p[2]) + log(det(V)) +
           (t(x - M) %*% chol2inv(V) %*% (x - M)) / p[2]