[Vegan-commits] r315 - in pkg: . R inst

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

Author: jarioksa
Date: 2008-04-14 02:02:15 +0200 (Mon, 14 Apr 2008)
New Revision: 315

Modified:
   pkg/DESCRIPTION
   pkg/R/adonis.R
   pkg/inst/ChangeLog
Log:
adonis handles properly aliased variables and rank deficit models

Modified: pkg/DESCRIPTION
===================================================================
--- pkg/DESCRIPTION	2008-04-10 16:32:30 UTC (rev 314)
+++ pkg/DESCRIPTION	2008-04-14 00:02:15 UTC (rev 315)
@@ -1,7 +1,7 @@
 Package: vegan
 Title: Community Ecology Package
-Version: 1.12-9
-Date: Apr 6, 2008
+Version: 1.12-10
+Date: Apr 13, 2008
 Author: Jari Oksanen, Roeland Kindt, Pierre Legendre, Bob O'Hara, Gavin L. Simpson, 
   M. Henry H. Stevens  
 Maintainer: Jari Oksanen <jari.oksanen at oulu.fi>

Modified: pkg/R/adonis.R
===================================================================
--- pkg/R/adonis.R	2008-04-10 16:32:30 UTC (rev 314)
+++ pkg/R/adonis.R	2008-04-14 00:02:15 UTC (rev 315)
@@ -16,9 +16,13 @@
     rhs <- model.matrix(formula, rhs.frame) # and finally the model.matrix
     options(contrasts=op.c)
     grps <- attr(rhs, "assign")
+    qrhs <- qr(rhs)
+    ## Take care of aliased variables and pivoting in rhs
+    rhs <- rhs[, qrhs$pivot, drop=FALSE]
+    rhs <- rhs[, 1:qrhs$rank, drop=FALSE]
+    grps <- grps[qrhs$pivot][1:qrhs$rank]
     u.grps <- unique(grps)
-    nterms <- max(u.grps)
-    
+    nterms <- length(u.grps) - 1
     H.s <- lapply(2:length(u.grps),
                   function(j) {Xj <- rhs[, grps %in% u.grps[1:j] ]
                                qrX <- qr(Xj, tol=1e-7)
@@ -41,8 +45,8 @@
     SS.Exp.each <- c(SS.Exp.comb - c(0,SS.Exp.comb[-nterms]) )
     SS.Res <- sum(diag( ( G %*% (I-H.s[[nterms]]) )))
     df.Exp <- sapply(u.grps[-1], function(i) sum(grps==i) )
-    df.Res <- n - dim(rhs)[2]
-    beta <-  qr.coef( qr(rhs), as.matrix(lhs) )
+    df.Res <- n - qrhs$rank
+    beta <-  qr.coef(qrhs, as.matrix(lhs) )
     colnames(beta) <- colnames(lhs)
     F.Mod <- (SS.Exp.each/df.Exp) / (SS.Res/df.Res)
 
@@ -75,7 +79,7 @@
                       F.Model = c(F.Mod, NA,NA),
                       R2 = SumsOfSqs/SumsOfSqs[length(SumsOfSqs)],
                       P = c(rowSums(t(f.perms) > F.Mod)/permutations, NA, NA))
-    rownames(tab) <- c(attr(attr(rhs.frame, "terms"), "term.labels"),
+    rownames(tab) <- c(attr(attr(rhs.frame, "terms"), "term.labels")[u.grps],
                        "Residuals", "Total")
     colnames(tab)[ncol(tab)] <- "Pr(>F)"
     out <- list(aov.tab = tab, call = match.call(), 

Modified: pkg/inst/ChangeLog
===================================================================
--- pkg/inst/ChangeLog	2008-04-10 16:32:30 UTC (rev 314)
+++ pkg/inst/ChangeLog	2008-04-14 00:02:15 UTC (rev 315)
@@ -2,8 +2,14 @@
 
 VEGAN DEVEL VERSIONS at http://r-forge.r-project.org/
 
-Version 1.12-9 (opened April 6, 2008)
+Version 1.12-10 (opened April 13, 2008)
 
+	* adonis: the fix for the unused factor levels (r312) revealed
+	that adonis did not handle properly aliased variables and
+	pivoting. Results need checking.
+	
+Version 1.12-9 (closed April 13, 2008)
+
 	* ordiArrowMul: added function 'ordiArrowMul' to improve automatic
 	scaling of biplot arrows (in plot.cca, text.cca, points.cca) and
 	fitted vectors (in plot.envfit). Reversed axes 'xlim=c(1,-1)' work
@@ -22,8 +28,8 @@
 	can be deduced. This was the intended behaviour but allPerms failed
 	in such cases.
 
-	* adonis: Was not dropping empty levels of factors, which was 
-	leading to incorrect results (d.f. was inflated). Reported by
+	* adonis: Was not dropping empty levels of factors, which was
+	leading to incorrect results (d.f. was inflated).  Reported by
 	Tyler Smith.
 	
 Version 1.12-8 (closed April 6, 2008)