[Vegan-commits] r1617 - in pkg/vegan: . R inst
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon May 30 14:43:32 CEST 2011
Author: jarioksa
Date: 2011-05-30 14:43:32 +0200 (Mon, 30 May 2011)
New Revision: 1617
Modified:
pkg/vegan/DESCRIPTION
pkg/vegan/R/capscale.R
pkg/vegan/R/rda.default.R
pkg/vegan/R/rda.formula.R
pkg/vegan/inst/ChangeLog
Log:
rda & capscale handle zero-rank pCCA, CCA & CA components as cca after r1517
Modified: pkg/vegan/DESCRIPTION
===================================================================
--- pkg/vegan/DESCRIPTION 2011-05-29 14:46:57 UTC (rev 1616)
+++ pkg/vegan/DESCRIPTION 2011-05-30 12:43:32 UTC (rev 1617)
@@ -1,7 +1,7 @@
Package: vegan
Title: Community Ecology Package
-Version: 1.18-31
-Date: May 26, 2011
+Version: 1.18-32
+Date: May 30, 2011
Author: Jari Oksanen, F. Guillaume Blanchet, Roeland Kindt, Pierre Legendre,
Peter R. Minchin, R. B. O'Hara, Gavin L. Simpson, Peter Solymos,
M. Henry H. Stevens, Helene Wagner
Modified: pkg/vegan/R/capscale.R
===================================================================
--- pkg/vegan/R/capscale.R 2011-05-29 14:46:57 UTC (rev 1616)
+++ pkg/vegan/R/capscale.R 2011-05-30 12:43:32 UTC (rev 1617)
@@ -95,9 +95,10 @@
colnames(sol$CCA$wa) <- colnames(sol$CCA$v) <-
paste("CAP", 1:ncol(sol$CCA$u), sep = "")
}
- if (!is.null(sol$CA)) {
+ if (!is.null(sol$CA) && sol$CA$rank > 0) {
colnames(sol$CA$u) <- names(sol$CA$eig) <- colnames(sol$CA$v) <-
paste("MDS", 1:ncol(sol$CA$u), sep = "")
+ }
## update for negative eigenvalues
poseig <- length(sol$CA$eig)
if (any(X$eig < 0)) {
@@ -106,7 +107,6 @@
sol$tot.chi <- sol$tot.chi + sol$CA$imaginary.chi
sol$CA$imaginary.rank <- length(negax)
sol$CA$imaginary.u.eig <- X$negaxes
- }
}
if (!is.null(comm)) {
comm <- scale(comm, center = TRUE, scale = FALSE)
@@ -125,7 +125,7 @@
"/")
comm <- qr.resid(sol$CCA$QR, comm)
}
- if (!is.null(sol$CA)) {
+ if (!is.null(sol$CA) && sol$CA$rank > 0) {
sol$CA$v.eig <- t(comm) %*% sol$CA$u/sqrt(k)
sol$CA$v <- sweep(sol$CA$v.eig, 2, sqrt(sol$CA$eig[1:poseig]),
"/")
Modified: pkg/vegan/R/rda.default.R
===================================================================
--- pkg/vegan/R/rda.default.R 2011-05-29 14:46:57 UTC (rev 1616)
+++ pkg/vegan/R/rda.default.R 2011-05-30 12:43:32 UTC (rev 1617)
@@ -23,11 +23,12 @@
Fit = Z, envcentre = attr(Z.r, "scaled:center"))
Xbar <- qr.resid(Q, Xbar)
}
+ if (tmp < ZERO)
+ pCCA$tot.chi <- 0
}
else Z.r <- NULL
if (!missing(Y) && !is.null(Y)) {
Y <- as.matrix(Y)
- rawmat <- Y
Y.r <- scale(Y, center = TRUE, scale = FALSE)
Q <- qr(cbind(Z.r, Y.r), tol = ZERO)
if (is.null(pCCA))
@@ -72,31 +73,43 @@
CCA$envcentre <- attr(Y.r, "scaled:center")
CCA$Xbar <- Xbar
Xbar <- qr.resid(Q, Xbar)
+ } else {
+ CCA <- list(eig = 0, rank = rank, qrank = qrank, tot.chi = 0,
+ QR = Q, Xbar = Xbar)
+ u <- matrix(0, nrow=nrow(sol$u), ncol=0)
+ v <- matrix(0, nrow=nrow(sol$v), ncol=0)
+ CCA$u <- CCA$u.eig <- CCA$wa <- CCA$wa.eig <- u
+ CCA$v <- CCA$v.eig <- v
+ CCA$biplot <- matrix(0, 0, 0)
+ CCA$alias <- colnames(Y.r)
}
}
Q <- qr(Xbar)
- if (Q$rank > 0) {
- sol <- svd(Xbar)
- sol$d <- sol$d/sqrt(NR)
- ax.names <- paste("PC", 1:length(sol$d), sep = "")
- colnames(sol$u) <- ax.names
- colnames(sol$v) <- ax.names
- names(sol$d) <- ax.names
- rownames(sol$u) <- rownames(X)
- rownames(sol$v) <- colnames(X)
- rank <- min(Q$rank, sum(sol$d > ZERO))
- if (rank) {
- CA <- list(eig = (sol$d[1:rank]^2))
- CA$u <- as.matrix(sol$u)[, 1:rank, drop = FALSE]
- CA$v <- as.matrix(sol$v)[, 1:rank, drop = FALSE]
- CA$u.eig <- sweep(as.matrix(CA$u), 2, sol$d[1:rank],
- "*")
- CA$v.eig <- sweep(as.matrix(CA$v), 2, sol$d[1:rank],
- "*")
- CA$rank <- rank
- CA$tot.chi <- sum(CA$eig)
- CA$Xbar <- Xbar
- }
+ sol <- svd(Xbar)
+ sol$d <- sol$d/sqrt(NR)
+ ax.names <- paste("PC", 1:length(sol$d), sep = "")
+ colnames(sol$u) <- ax.names
+ colnames(sol$v) <- ax.names
+ names(sol$d) <- ax.names
+ rownames(sol$u) <- rownames(X)
+ rownames(sol$v) <- colnames(X)
+ rank <- min(Q$rank, sum(sol$d > ZERO))
+ if (rank) {
+ CA <- list(eig = (sol$d[1:rank]^2))
+ CA$u <- as.matrix(sol$u)[, 1:rank, drop = FALSE]
+ CA$v <- as.matrix(sol$v)[, 1:rank, drop = FALSE]
+ CA$u.eig <- sweep(as.matrix(CA$u), 2, sol$d[1:rank],
+ "*")
+ CA$v.eig <- sweep(as.matrix(CA$v), 2, sol$d[1:rank],
+ "*")
+ CA$rank <- rank
+ CA$tot.chi <- sum(CA$eig)
+ CA$Xbar <- Xbar
+ } else { # zero rank: no residual component
+ CA <- list(eig = 0, rank = rank, tot.chi = 0,
+ Xbar = Xbar)
+ CA$u <- CA$u.eig <- matrix(0, nrow(sol$u), 0)
+ CA$v <- CA$v.eig <- matrix(0, nrow(sol$v), 0)
}
call <- match.call()
call[[1]] <- as.name("rda")
Modified: pkg/vegan/R/rda.formula.R
===================================================================
--- pkg/vegan/R/rda.formula.R 2011-05-29 14:46:57 UTC (rev 1616)
+++ pkg/vegan/R/rda.formula.R 2011-05-30 12:43:32 UTC (rev 1617)
@@ -1,6 +1,6 @@
-"rda.formula" <-
-function (formula, data, scale = FALSE, na.action = na.fail,
- subset = NULL, ...)
+`rda.formula` <-
+ function (for§mula, data, scale = FALSE, na.action = na.fail,
+ subset = NULL, ...)
{
if (missing(data)) {
data <- parent.frame()
@@ -10,7 +10,7 @@
d <- ordiParseFormula(formula, data, na.action = na.action,
subset = substitute(subset))
sol <- rda.default(d$X, d$Y, d$Z, scale)
- if (!is.null(sol$CCA))
+ if (!is.null(sol$CCA) && sol$CCA$rank > 0)
sol$CCA$centroids <- centroids.cca(sol$CCA$wa, d$modelframe)
if (!is.null(sol$CCA$alias))
sol$CCA$centroids <- unique(sol$CCA$centroids)
Modified: pkg/vegan/inst/ChangeLog
===================================================================
--- pkg/vegan/inst/ChangeLog 2011-05-29 14:46:57 UTC (rev 1616)
+++ pkg/vegan/inst/ChangeLog 2011-05-30 12:43:32 UTC (rev 1617)
@@ -2,8 +2,17 @@
VEGAN DEVEL VERSIONS at http://r-forge.r-project.org/
-Version 1.18-31 (opened May 26, 2011)
+Version 1.18-32 (opened May 30, 2011)
+ * rda & capscale: similar changes as in cca in r1517 (version
+ 1.18-24). CA component is always returned, even with zero rank,
+ and CCA and pCCA components are NULL only if they were not
+ specified originally, and if they become zero rank because of
+ aliasing etc, they are returned as zero components (instead of
+ NULL).
+
+Version 1.18-31 (closed May 30, 2011)
+
* ordiR2step: more informative about rejecting 'scope': partial
models are not (currently) accepted due to the design of
RsquareAdj.rda().
More information about the Vegan-commits
mailing list