[Vegan-commits] r1356 - in pkg/vegan: . R inst

[noreply at r-forge.r-project.org]

Author: jarioksa
Date: 2010-11-09 16:44:31 +0100 (Tue, 09 Nov 2010)
New Revision: 1356

Modified:
   pkg/vegan/DESCRIPTION
   pkg/vegan/R/prc.R
   pkg/vegan/R/print.summary.prc.R
   pkg/vegan/R/summary.prc.R
   pkg/vegan/inst/ChangeLog
Log:
Cajo's new direct implementation of PRC

Modified: pkg/vegan/DESCRIPTION
===================================================================
--- pkg/vegan/DESCRIPTION	2010-11-09 09:23:48 UTC (rev 1355)
+++ pkg/vegan/DESCRIPTION	2010-11-09 15:44:31 UTC (rev 1356)
@@ -1,7 +1,7 @@
 Package: vegan
 Title: Community Ecology Package
-Version: 1.18-15
-Date: November 1, 2010
+Version: 1.18-16
+Date: November 9, 2010
 Author: Jari Oksanen, F. Guillaume Blanchet, Roeland Kindt, Pierre Legendre, 
    R. B. O'Hara, Gavin L. Simpson, Peter Solymos, M. Henry H. Stevens, 
    Helene Wagner  

Modified: pkg/vegan/R/prc.R
===================================================================
--- pkg/vegan/R/prc.R	2010-11-09 09:23:48 UTC (rev 1355)
+++ pkg/vegan/R/prc.R	2010-11-09 15:44:31 UTC (rev 1356)
@@ -1,7 +1,7 @@
-`prc` <-
-    function (response, treatment, time,...) 
+`prc`  <-
+    function (response, treatment, time, ...)
 {
-    extras <- match.call(expand.dots=FALSE)$...
+    extras <- match.call(expand.dots = FALSE)$...
     if (is.null(extras$data))
         data <- parent.frame()
     else
@@ -15,11 +15,16 @@
     mf <- model.frame(fla, data)
     if (!all(sapply(mf, is.factor)))
         stop(x, " and ", z, " must be factors")
-    if (any(sapply(mf, is.ordered))) 
+    if (any(sapply(mf, is.ordered)))
         stop(x, " or ", z, " cannot be ordered factors")
-    fla <- as.formula(paste(y, "~", z, "*", x, "+ Condition(", 
-                            z, ")"))
-    mod <- rda(fla, ...)
+    fla.zx <- as.formula(paste("~", z, ":", x))
+    fla.z <- as.formula(paste("~", z))
+    # delete first (control) level from the design matrix
+    X = model.matrix(fla.zx, data)[,-c(seq_len(nlevels(time)+1))]
+    Z = model.matrix(fla.z, data)[,-1]
+    mod <- rda(response, X, Z, ...)
+    mod$terminfo$xlev = list(levels(time), levels(treatment))
+    names(mod$terminfo$xlev) = c(paste(z), paste(x))
     mod$call <- match.call()
     class(mod) <- c("prc", class(mod))
     mod

Modified: pkg/vegan/R/print.summary.prc.R
===================================================================
--- pkg/vegan/R/print.summary.prc.R	2010-11-09 09:23:48 UTC (rev 1355)
+++ pkg/vegan/R/print.summary.prc.R	2010-11-09 15:44:31 UTC (rev 1356)
@@ -5,7 +5,9 @@
     cat(deparse(x$call), "\n")
     cat("Species scores:\n")
     print(x$sp, digits=x$digits, ...)
-    cat("\nCoefficients for", paste(x$names, collapse=":"), "interaction\n")
+    cat("\nCoefficients for",
+        paste(x$names[1], "+", paste(x$names, collapse=":")),
+        "interaction\n")
     cat(paste("which are contrasts to", x$names[2], x$corner, "\n"))
     cat(paste(c("rows are",", columns are"), x$names[2:1], collapse=""))
     cat("\n")

Modified: pkg/vegan/R/summary.prc.R
===================================================================
--- pkg/vegan/R/summary.prc.R	2010-11-09 09:23:48 UTC (rev 1355)
+++ pkg/vegan/R/summary.prc.R	2010-11-09 15:44:31 UTC (rev 1356)
@@ -1,24 +1,19 @@
 `summary.prc` <-
-    function (object, axis = 1, scaling = 3, digits = 4, ...) 
+    function (object, axis = 1, scaling = 3, digits = 4, ...)
 {
-    species <- drop(scores(object, scaling = scaling, display="sp",
-                           choices=axis, ...))
-    sites <- drop(scores(object, scaling = scaling, display="lc",
-                         choices=axis, ... )) 
+    sc = scores(object, scaling = scaling, display = c("sp", "lc"),
+                           choices=axis, ...)
     ## coef for scaled sites (coef(object) gives for orthonormal)
-    b <- qr.coef(object$CCA$QR, sites)
+    b <- qr.coef(object$CCA$QR, sc$constraints)
     prnk <- object$pCCA$rank
     lentreat <- length(object$terminfo$xlev[[2]])
-    lenb <- length(b)
-    b <- b[-(1:(2 * prnk))]
-    bx <- b[1:(lentreat - 1)]
-    by <- b[lentreat:length(b)] + rep(bx, each = length(object$terminfo$xlev[[1]])-1)
-    b <- cbind(bx, matrix(by, nrow = lentreat - 1, byrow = TRUE))
+    b = matrix(b[-(1:prnk)], nrow = lentreat-1, byrow = TRUE)
     rownames(b) <- (object$terminfo$xlev[[2]])[-1]
     colnames(b) <- object$terminfo$xlev[[1]]
-    out <- list(sp = species, coefficients = b, names = names(object$terminfo$xlev), 
-                corner = (object$terminfo$xlev[[2]])[1], call = object$call, 
-                digits = digits)
+    out <- list(sp = drop(sc$species), coefficients = b,
+                names = names(object$terminfo$xlev),
+                corner = (object$terminfo$xlev[[2]])[1],
+                call = object$call, digits = digits)
     class(out) <- "summary.prc"
     out
 }

Modified: pkg/vegan/inst/ChangeLog
===================================================================
--- pkg/vegan/inst/ChangeLog	2010-11-09 09:23:48 UTC (rev 1355)
+++ pkg/vegan/inst/ChangeLog	2010-11-09 15:44:31 UTC (rev 1356)
@@ -2,8 +2,14 @@
 
 VEGAN DEVEL VERSIONS at http://r-forge.r-project.org/
 
-Version 1.18-15 (opened November 1, 2010)
+Version 1.18-16 (opened November 9, 20109
 
+	* prc: Cajo ter Braak wrote new code for more direct
+	implementation of PRC. This drops formula interface and directly
+	finds model matrices avoiding aliased terms.
+	
+Version 1.18-15 (closed November 9, 2010)
+
 	* procrustes: 'translation' needs to take into account 'scale',
 	although Mardia et al. omit it. Reported, analysed and fix
 	suggested by Christian Dudel (Bochum).