]> git.donarmstrong.com Git - ape.git/blobdiff - R/rTrait.R
final commit for ape 3.0-8
[ape.git] / R / rTrait.R
index 2c18d0940af08ada1fcaf1a67c8cba9b5402b1ed..1b5852d102024a01a9a1a312643e82f57f9c0d76 100644 (file)
@@ -1,8 +1,8 @@
-## rTrait.R (2011-04-02)
+## rTrait.R (2012-02-09)
 
 ##   Trait Evolution
 
-## Copyright 2010-2011 Emmanuel Paradis
+## Copyright 2010-2012 Emmanuel Paradis
 
 ## This file is part of the R-package `ape'.
 ## See the file ../COPYING for licensing issues.
@@ -64,7 +64,7 @@ rTraitDisc <-
         diag(Q) <- -rowSums(Q)
         for (i in N:1) {
             p <- matexpo(Q * el[i])[x[anc[i]], ]
-            x[des[i]] <- .Internal(sample(k, size = 1, FALSE, prob = p))
+            x[des[i]] <- sample.int(k, size = 1, FALSE, prob = p)
         }
     }
 
@@ -82,7 +82,7 @@ rTraitDisc <-
 
 rTraitCont <-
     function(phy, model = "BM", sigma = 0.1, alpha = 1, theta = 0,
-             ancestor = FALSE, root.value = 0, linear = TRUE, ...)
+             ancestor = FALSE, root.value = 0, ...)
 {
     if (is.null(phy$edge.length))
         stop("tree has no branch length")
@@ -115,7 +115,6 @@ rTraitCont <-
             if (length(theta) == 1) theta <- rep(theta, N)
             else if (length(theta) != N)
                 stop("'theta' must have one or Nedge(phy) elements")
-            if (!linear) model <- model + 1L
         }
         .C("rTraitCont", as.integer(model), as.integer(N),
            as.integer(anc - 1L), as.integer(des - 1L), el,
@@ -135,7 +134,7 @@ rTraitCont <-
 
 rTraitMult <-
     function(phy, model, p = 1, root.value = rep(0, p), ancestor = FALSE,
-             as.factor = NULL, ...)
+             asFactor = NULL, trait.labels = paste("x", 1:p, sep = ""), ...)
 {
     phy <- reorder(phy, "pruningwise")
     n <- length(phy$tip.label)
@@ -152,6 +151,8 @@ rTraitMult <-
     el <- phy$edge.length
     if (is.null(el)) el <- numeric(N)
 
+    environment(model) <- environment() # to find 'p'
+
     for (i in N:1) x[des[i], ] <- model(x[anc[i], ], el[i], ...)
 
     if (ancestor) {
@@ -162,8 +163,9 @@ rTraitMult <-
         rownames(x) <- phy$tip.label
     }
     x <- as.data.frame(x)
-    if (!is.null(as.factor)) {
-        for (i in as.factor) {
+    names(x) <- trait.labels
+    if (!is.null(asFactor)) {
+        for (i in asFactor) {
             y <- x[, i]
             x[, i] <- factor(y, labels = LETTERS[1:length(unique(y))])
         }