[Vegan-commits] r406 - in branches/1.13: R inst man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Jun 9 08:24:24 CEST 2008


Author: jarioksa
Date: 2008-06-09 08:24:24 +0200 (Mon, 09 Jun 2008)
New Revision: 406

Modified:
   branches/1.13/R/alias.cca.R
   branches/1.13/R/calibrate.cca.R
   branches/1.13/R/cca.default.R
   branches/1.13/R/extractAIC.cca.R
   branches/1.13/R/permutest.cca.R
   branches/1.13/R/predict.cca.R
   branches/1.13/R/predict.rda.R
   branches/1.13/R/print.cca.R
   branches/1.13/R/rda.default.R
   branches/1.13/inst/ChangeLog
   branches/1.13/man/cca.object.Rd
Log:
merged r400:405 to branches/1.13 (fixes to rda/cca rank issue)

Modified: branches/1.13/R/alias.cca.R
===================================================================
--- branches/1.13/R/alias.cca.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/alias.cca.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -14,7 +14,7 @@
     R <- R[1:min(dim(R)), , drop = FALSE]
     R[lower.tri(R)] <- 0
     d <- dim(R)
-    rank <- object$CCA$rank
+    rank <- object$CCA$QR$rank
     p <- d[2]
     value$Complete <- if (is.null(p) || rank == p) 
         NULL

Modified: branches/1.13/R/calibrate.cca.R
===================================================================
--- branches/1.13/R/calibrate.cca.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/calibrate.cca.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -5,6 +5,8 @@
         stop("does not work with conditioned (partial) models")
     if (is.null(object$CCA))
         stop("needs constrained model")
+    if (object$CCA$rank < object$CCA$qrank)
+        stop("rank of constraints is higher than rank of dependent data")
     if (rank != "full")
         rank <- min(rank, object$CCA$rank)
     else
@@ -13,11 +15,12 @@
         wa <- object$CCA$wa        
     else
         wa <- predict(object, type="wa", newdata=newdata)
-    b <- (coef(object))[object$CCA$QR$pivot[1:object$CCA$rank], , drop=FALSE]
+    qrank <- object$CCA$qrank
+    b <- (coef(object))[object$CCA$QR$pivot[1:qrank], , drop=FALSE]
     b <- solve(b)
-    pred <- wa[ , 1:rank, drop=FALSE]  %*% b[1:rank, , drop =FALSE]
+    pred <- wa[ , 1:rank, drop=FALSE]  %*% b[1:qrank, , drop =FALSE]
     envcen <- object$CCA$envcentre[object$CCA$QR$pivot]
-    envcen <- envcen[1:object$CCA$rank]
+    envcen <- envcen[1:object$CCA$qrank]
     pred <- sweep(pred, 2, envcen, "+")
     pred
 }

Modified: branches/1.13/R/cca.default.R
===================================================================
--- branches/1.13/R/cca.default.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/cca.default.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -79,6 +79,7 @@
             CCA$biplot <- cor(Y.r[, oo, drop = FALSE], sol$u[, 
                                         1:rank, drop = FALSE])
             CCA$rank <- rank
+            CCA$qrank <- qrank
             CCA$tot.chi <- sum(CCA$eig)
             CCA$QR <- Q
             CCA$envcentre <- attr(Y.r, "centre")

Modified: branches/1.13/R/extractAIC.cca.R
===================================================================
--- branches/1.13/R/extractAIC.cca.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/extractAIC.cca.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -3,7 +3,7 @@
 {
    n <- nrow(fit$CA$Xbar)
    edf <- 1
-   if (!is.null(fit$CCA$rank)) edf <- edf + fit$CCA$rank
+   if (!is.null(fit$CCA$rank)) edf <- edf + fit$CCA$qrank
    if (!is.null(fit$pCCA$rank)) edf <- edf + fit$pCCA$rank
    #edf <- n - fit$CA$rank
    RSS <- deviance(fit)

Modified: branches/1.13/R/permutest.cca.R
===================================================================
--- branches/1.13/R/permutest.cca.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/permutest.cca.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -18,7 +18,7 @@
     else {
         Chi.z <- x$CCA$tot.chi
         names(Chi.z) <- "Model"
-        q <- x$CCA$rank
+        q <- x$CCA$qrank
     }
     Chi.xz <- x$CA$tot.chi
     names(Chi.xz) <- "Residual"

Modified: branches/1.13/R/predict.cca.R
===================================================================
--- branches/1.13/R/predict.cca.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/predict.cca.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -39,7 +39,8 @@
             }
             E <- sweep(E, 2, c(object$pCCA$envcentre, object$CCA$envcentre), 
                        "-")
-            p1 <- object[[model]]$QR$pivot[1:object[[model]]$rank]
+            Q <- object[[model]]$QR
+            p1 <- Q$pivot[1:Q$rank]
             u <- E[, p1, drop = FALSE] %*% coef(object)[p1, , 
                          drop = FALSE]
             u <- u[, 1:take, drop = FALSE]

Modified: branches/1.13/R/predict.rda.R
===================================================================
--- branches/1.13/R/predict.rda.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/predict.rda.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -44,7 +44,8 @@
                                       object$terminfo$xlev)
                 E <- cbind(d$Z, d$Y)
             }
-            p1 <- object[[model]]$QR$pivot[1:object[[model]]$rank]
+            Q <- object[[model]]$QR
+            p1 <- Q$pivot[1:Q$rank]
             E <- sweep(E, 2, c(object$pCCA$envcentre, object$CCA$envcentre), 
                        "-")
             u <- E[, p1, drop = FALSE] %*% coef(object)[p1, , 

Modified: branches/1.13/R/print.cca.R
===================================================================
--- branches/1.13/R/print.cca.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/print.cca.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -4,13 +4,6 @@
     if (inherits(x, "pcaiv")) {
         stop("this is an ade4 object which vegan cannot handle")
     }
-    if (!is.null(x$CCA)) {
-        qrank <- x$CCA$QR$rank
-        if (!is.null(x$pCCA))
-            qrank <- qrank - x$pCCA$rank
-        if (x$CCA$rank < qrank)
-            warning("rank of constraints higher than the rank of dependent data\nvegan may not handle this -- wait for a bug fix release", call. = FALSE)
-    }
     cat("\nCall:\n")
     cat(deparse(x$call), "\n\n")
     chi <- rbind(x$tot.chi, x$pCCA$tot.chi, x$CCA$tot.chi, x$CA$tot.chi)

Modified: branches/1.13/R/rda.default.R
===================================================================
--- branches/1.13/R/rda.default.R	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/R/rda.default.R	2008-06-09 06:24:24 UTC (rev 406)
@@ -63,6 +63,7 @@
             CCA$biplot <- cor(Y.r[, oo, drop = FALSE], sol$u[, 
                                         1:rank, drop = FALSE])
             CCA$rank <- rank
+            CCA$qrank <- qrank
             CCA$tot.chi <- sum(CCA$eig)
             CCA$QR <- Q
             CCA$envcentre <- attr(Y.r, "scaled:center")

Modified: branches/1.13/inst/ChangeLog
===================================================================
--- branches/1.13/inst/ChangeLog	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/inst/ChangeLog	2008-06-09 06:24:24 UTC (rev 406)
@@ -5,6 +5,8 @@
 
 Version 1.13-1 (opened May 22, 2008)
 
+	* merged r400:405 -- should fix rank problems in cca/rda.
+
 	* merged r396: cca/rda/capscale print issues a warning when rank
 	of constraints+conditions < rank of dependent data, because vegan
 	cannot handle these cases properly.

Modified: branches/1.13/man/cca.object.Rd
===================================================================
--- branches/1.13/man/cca.object.Rd	2008-06-09 04:39:11 UTC (rev 405)
+++ branches/1.13/man/cca.object.Rd	2008-06-09 06:24:24 UTC (rev 406)
@@ -69,7 +69,10 @@
       \code{\link{predict.cca}}, \code{\link{predict.rda}},
       \code{\link{calibrate.cca}}.  For possible uses of this component,
       see \code{\link{qr}}. In \code{pCCA} and \code{CCA}.} 
-    \item{\code{rank}}{The rank of the component.}
+    \item{\code{rank}}{The rank of the ordination component.}
+    \item{\code{qrank}}{The rank of the constraints which is the
+     difference of the ranks of QR decompositions in \code{pCCA} and
+     \code{CCA} components. Only in \code{CCA}.}
     \item{\code{tot.chi}}{Total inertia or the sum of all eigenvalues of the
       component.}
     \item{\code{u}}{(Weighted) orthonormal site scores.  Please note that



More information about the Vegan-commits mailing list