[Vegan-commits] r1517 - in pkg/vegan: . R inst
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Feb 28 17:00:36 CET 2011
Author: jarioksa
Date: 2011-02-28 17:00:35 +0100 (Mon, 28 Feb 2011)
New Revision: 1517
Modified:
pkg/vegan/DESCRIPTION
pkg/vegan/R/cca.default.R
pkg/vegan/R/cca.formula.R
pkg/vegan/R/print.cca.R
pkg/vegan/inst/ChangeLog
Log:
cca handling of zero-rank pCCA, CCA and CA items (experimental code)
Modified: pkg/vegan/DESCRIPTION
===================================================================
--- pkg/vegan/DESCRIPTION 2011-02-28 15:53:47 UTC (rev 1516)
+++ pkg/vegan/DESCRIPTION 2011-02-28 16:00:35 UTC (rev 1517)
@@ -1,7 +1,7 @@
Package: vegan
Title: Community Ecology Package
-Version: 1.18-23
-Date: February 16, 2011
+Version: 1.18-24
+Date: February 28, 2011
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/cca.default.R
===================================================================
--- pkg/vegan/R/cca.default.R 2011-02-28 15:53:47 UTC (rev 1516)
+++ pkg/vegan/R/cca.default.R 2011-02-28 16:00:35 UTC (rev 1517)
@@ -39,11 +39,12 @@
Fit = Z, envcentre = attr(Z.r, "centre"))
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 <- weight.centre(Y, rowsum)
Q <- qr(cbind(Z.r, Y.r), tol = ZERO)
if (is.null(pCCA))
@@ -86,42 +87,59 @@
CCA$QR <- Q
CCA$envcentre <- attr(Y.r, "centre")
CCA$Xbar <- Xbar
- Xbar <- qr.resid(Q, Xbar)
- if (exists("exclude.spec")) {
- attr(CCA$v, "na.action") <- exclude.spec
- attr(CCA$v.eig, "na.action") <- exclude.spec
- }
+ } else { # zero rank
+ 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)
}
+ Xbar <- qr.resid(Q, Xbar)
+ if (exists("exclude.spec")) {
+ attr(CCA$v, "na.action") <- exclude.spec
+ attr(CCA$v.eig, "na.action") <- exclude.spec
+ }
+
}
Q <- qr(Xbar)
- if (Q$rank > 0) {
- sol <- svd(Xbar)
- ax.names <- paste("CA", 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 <- sweep(as.matrix(sol$u[, 1:rank, drop = FALSE]),
- 1, 1/sqrt(rowsum), "*")
- CA$v <- sweep(as.matrix(sol$v[, 1:rank, drop = FALSE]),
- 1, 1/sqrt(colsum), "*")
- CA$u.eig <- sweep(CA$u, 2, sol$d[1:rank], "*")
- CA$v.eig <- sweep(CA$v, 2, sol$d[1:rank], "*")
- CA$rank <- rank
- CA$tot.chi <- sum(CA$eig)
- CA$Xbar <- Xbar
- if (exists("exclude.spec")) {
- attr(CA$v, "na.action") <- exclude.spec
- attr(CA$v.eig, "na.action") <- exclude.spec
- }
- }
+ sol <- svd(Xbar)
+ ax.names <- paste("CA", 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 <- sweep(as.matrix(sol$u[, 1:rank, drop = FALSE]),
+ 1, 1/sqrt(rowsum), "*")
+ CA$v <- sweep(as.matrix(sol$v[, 1:rank, drop = FALSE]),
+ 1, 1/sqrt(colsum), "*")
+ CA$u.eig <- sweep(CA$u, 2, sol$d[1:rank], "*")
+ CA$v.eig <- sweep(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)
}
+ if (exists("exclude.spec")) {
+ attr(CA$v, "na.action") <- exclude.spec
+ attr(CA$v.eig, "na.action") <- exclude.spec
+ }
call <- match.call()
call[[1]] <- as.name("cca")
+ ## computed pCCA$rank was needed before, but zero it here
+ if (!is.null(pCCA) && pCCA$tot.chi == 0)
+ pCCA$rank <- 0
sol <- list(call = call, grand.total = gran.tot, rowsum = rowsum,
colsum = colsum, tot.chi = tot.chi, pCCA = pCCA, CCA = CCA,
CA = CA)
Modified: pkg/vegan/R/cca.formula.R
===================================================================
--- pkg/vegan/R/cca.formula.R 2011-02-28 15:53:47 UTC (rev 1516)
+++ pkg/vegan/R/cca.formula.R 2011-02-28 16:00:35 UTC (rev 1517)
@@ -9,7 +9,7 @@
d <- ordiParseFormula(formula, data, na.action = na.action,
subset = substitute(subset))
sol <- cca.default(d$X, d$Y, d$Z)
- if (!is.null(sol$CCA))
+ if (!is.null(sol$CCA) && sol$CCA$rank > 0)
sol$CCA$centroids <- centroids.cca(sol$CCA$wa, d$modelframe,
sol$rowsum)
if (!is.null(sol$CCA$alias))
Modified: pkg/vegan/R/print.cca.R
===================================================================
--- pkg/vegan/R/print.cca.R 2011-02-28 15:53:47 UTC (rev 1516)
+++ pkg/vegan/R/print.cca.R 2011-02-28 16:00:35 UTC (rev 1517)
@@ -37,11 +37,11 @@
cat(length(sp.na), "species",
ifelse(length(sp.na)==1, "(variable)", "(variables)"),
"deleted due to missingness\n")
- if (!is.null(x$CCA)) {
+ if (!is.null(x$CCA) && x$CCA$rank > 0) {
cat("\nEigenvalues for constrained axes:\n")
print(x$CCA$eig, digits = digits, ...)
}
- if (!is.null(x$CA)) {
+ if (!is.null(x$CA) && x$CA$rank > 0) {
ax.lim <- 8
ax.trig <- 16
cat("\nEigenvalues for unconstrained axes:\n")
Modified: pkg/vegan/inst/ChangeLog
===================================================================
--- pkg/vegan/inst/ChangeLog 2011-02-28 15:53:47 UTC (rev 1516)
+++ pkg/vegan/inst/ChangeLog 2011-02-28 16:00:35 UTC (rev 1517)
@@ -3,8 +3,28 @@
VEGAN DEVEL VERSIONS at http://r-forge.r-project.org/
-Version 1.18-23 (opened February 16, 2011)
+Version 1.18-24 (opened February 28, 2011)
+ * cca: will return NULL item for CCA only if constraints were not
+ given. If the CCA component has zero rank (= constraints were
+ aliased or were orthogonal to the data), will still return a CCA
+ item with zero eigenvalue, rank, scores etc plus info on
+ "alias". The residual component will never be NULL, but similar
+ zero-containing component is returned even when there is no
+ residual variation. Usually these changes only manifest in
+ arbitrary data sets, but we have had email queries in R-sig-eco on
+ completely aliased constrainsts when users assumed that model
+ cca(y ~ A + Condition(A)) would be similar to a random effects
+ model and A could be analysed both as a Condition and as a
+ Constraint instead of being aliased. Naturally, overfitted models
+ with empty residual components do appear often. Currently the
+ printed output includes these zero items with Inertia and Rank 0,
+ whereas NULL components are not displayed. The changes so far only
+ concern cca(); rda() and capscale() are only changed after (and
+ if) this change appears sensible and correct in cca().
+
+Version 1.18-23 (closed February 28, 2011)
+
* opened with the release of 1.17-7 based on version 1.18-22 on
February 16, 2011.
More information about the Vegan-commits
mailing list