[Vegan-commits] r1191 - in pkg/vegan: R inst man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon May 24 19:07:53 CEST 2010


Author: jarioksa
Date: 2010-05-24 19:07:53 +0200 (Mon, 24 May 2010)
New Revision: 1191

Added:
   pkg/vegan/R/ordiR2step.R
Modified:
   pkg/vegan/inst/ChangeLog
   pkg/vegan/man/ordistep.Rd
Log:
ordiR2step: new function to implement R2 forward choice of Blanchet & amates

Added: pkg/vegan/R/ordiR2step.R
===================================================================
--- pkg/vegan/R/ordiR2step.R	                        (rev 0)
+++ pkg/vegan/R/ordiR2step.R	2010-05-24 17:07:53 UTC (rev 1191)
@@ -0,0 +1,59 @@
+### Forward selection to maximize R2.adjusted, but stopping once the
+### R2.adjusted of the maximum model ('scope') is exceeded, after
+### Blanchet, Legendre & Borcard: Ecology 89, 2623--2623; 2008.
+
+`ordiR2step` <-
+    function(object, scope, trace = TRUE)
+{
+    ## Works only for rda(): cca() does not have (yet) R2.adjusted
+    if (!inherits(object, "rda"))
+        stop("can be used only with rda() or capscale()")
+    ## Get R2 of the original object
+    if (is.null(object$CCA))
+        R2.0 <- 0
+    else
+        R2.0 <- RsquareAdj(object)$adj.r.squared
+    ## Get R2 of the scope
+    if (inherits(scope, "rda")) {
+        R2.all <- RsquareAdj(scope)$adj.r.squared
+        scope <- formula(scope)
+    } else {
+        if (!inherits(scope, "formula"))
+            scope <- reformulate(scope)
+        R2.all <- RsquareAdj(update(object, scope))$adj.r.squared
+    }
+    ## Check that the full model can be evaluated
+    if (is.na(R2.all))
+        stop("the upper scope cannot be fitted (too many terms?)")
+    ## Step forward and continue as long as R2.adj improves and R2.adj
+    ## remains below R2.adj < R2.all
+    R2.previous <- R2.0
+    repeat {
+        adds <- add.scope(object, scope)
+        R2.adds <- numeric(length(adds))
+        names(R2.adds) <- adds
+        ## Loop over add scope
+        for (trm in seq_along(adds)) {
+            fla <- paste("~  . + ", adds[trm])
+            R2.adds[trm] <- RsquareAdj(update(object, fla))$adj.r.squared
+        }
+        best <- which.max(R2.adds)
+        if (trace) {
+            names(R2.adds) <- paste("+", names(R2.adds))
+            out <- sort(c("<All variables>" = R2.all, R2.adds), decreasing = TRUE)
+            out <- as.matrix(out)
+            colnames(out) <- "R2.adjusted"
+            print(out)
+            cat("\n")
+        }
+        ## See if the best should be kept
+        if (R2.adds[best] > R2.previous && R2.adds[best] <= R2.all) {
+            fla <- paste("~  . +", adds[best])
+            object <- update(object, fla)
+            R2.previous <- RsquareAdj(object)$adj.r.squared
+        } else {
+            break
+        }
+    }
+    object
+}

Modified: pkg/vegan/inst/ChangeLog
===================================================================
--- pkg/vegan/inst/ChangeLog	2010-05-21 11:58:11 UTC (rev 1190)
+++ pkg/vegan/inst/ChangeLog	2010-05-24 17:07:53 UTC (rev 1191)
@@ -4,6 +4,11 @@
 
 Version 1.18-5 (opened May 21, 2010) 
 
+	* ordiR2step: new function to perform forward model selection
+	following Blanchet, Legendre & Borcard (Ecology 89, 2623-2632;
+	2008) or by adjusted R2 until the adjusted R2 of the full model
+	(scope) is exceeded or adjusted R2 starts to decrease.
+
 	* RsquareAdj.default: handles vector arguments.
 
 	* ordiellipse: works only with 2D data, and now uses only first

Modified: pkg/vegan/man/ordistep.Rd
===================================================================
--- pkg/vegan/man/ordistep.Rd	2010-05-21 11:58:11 UTC (rev 1190)
+++ pkg/vegan/man/ordistep.Rd	2010-05-24 17:07:53 UTC (rev 1191)
@@ -1,6 +1,7 @@
 \name{ordistep}
 \Rdversion{1.1}
 \alias{ordistep}
+\alias{ordiR2step}
 
 \title{
   Choose a Model by Permutation Tests in Constrained Ordination 
@@ -8,22 +9,27 @@
 \description{
   Automatic stepwise model building for constrained ordination methods 
   (\code{\link{cca}}, \code{\link{rda}}, \code{\link{capscale}}). 
-  The function is modelled after \code{\link{step}} and can do forward,
-  backward and stepwise model selection.
+  The function \code{ordistep} is modelled after \code{\link{step}} and 
+  can do forward,backward and stepwise model selection using permutation tests.
+  Function \code{ordiR2step} peforms forward model choice solely on adjusted 
+  \eqn{R^2}{R2} in \code{\link{rda}} or \code{\link{capscale}}.
 }
 \usage{
 ordistep(object, scope, direction = c("both", "backward", "forward"), Pin = 0.05, Pout = 0.1, pstep = 100, perm.max = 1000, steps = 50, trace = TRUE, ...)
+ordiR2step(object, scope, trace = TRUE)
 }
 %- maybe also 'usage' for other objects documented here.
 \arguments{
   \item{object}{
-  An ordination object inheriting from \code{\link{cca}}.
+  An ordination object inheriting from \code{\link{cca}} (in \code{ordiR2step} 
+  this must be from \code{\link{rda}} or \code{\link{capscale}}).
 }
   \item{scope}{
   Defines the range of models examined in the stepwise search. 
   This should be either a single formula, or a list containing
   components \code{upper} and \code{lower}, both formulae. 
-  See \code{\link{step}} for details. 
+  See \code{\link{step}} for details. In \code{ordiR2step} this defines the
+  upper scope and can also be an ordination object from with the model is extracted. 
 }
   \item{direction}{
   The mode of stepwise search, can be one of \code{"both"},
@@ -75,6 +81,21 @@
   is selected adaptively with respect to the defined decision limit. It
   is often sensible to have \code{Pout} \eqn{>} \code{Pin} in stepwise
   models to avoid cyclic adds and drops of single terms. 
+
+  Function \code{ordiR2step} implements the forward selection to
+  maximize adjusted \eqn{R^2}{R2} (function \code{\link{RsquareAdj}})
+  at every step, and stopping when the adjusted \eqn{R^2}{R2} starts
+  to decrease, or the adjusted \eqn{R^2}{R2} of the \code{scope} is
+  exceeded (Blanchet et al. 2008). There are no permutation tests, but
+  the adusted \eqn{R^2}{R2} is the only criterion. If tests are
+  wished, the final model can be analysed using
+  \code{\link{anova.cca}} with argument \code{by = "term"} which
+  performs sequential test in the same order as the model was built in
+  \code{ordiR2step}.
+
+  Functions \code{ordistep} (based on \eqn{P} values) and \code{ordiR2step}
+  (based on adjusted \eqn{R^2}{R2} and hence on eigenvalues) can select
+  variables in different order.
 }
 
 \value{ 
@@ -85,6 +106,11 @@
   in the \code{anova} item.  
 }
 
+\references{
+  Blanchet, F. G., Legendre, P. & Borcard, D. (2008) Forward selection
+  of explanatory variables. \emph{Ecology} 89, 2623--2632. 
+}
+
 \author{
   Jari Oksanen
 }
@@ -95,7 +121,8 @@
   are \code{\link{add1.cca}} and \code{\link{drop1.cca}}, and the 
   function is modelled after standard \code{\link{step}} (which also can 
   be used directly but uses AIC for model choice, see 
-  \code{\link{extractAIC.cca}}). 
+  \code{\link{extractAIC.cca}}). Function \code{ordiR2step} builds upon
+  \code{\link{RsquareAdj}}.  
 }
 \examples{
 ## See add1.cca for another example
@@ -104,6 +131,7 @@
 mod1 <- rda(dune ~ ., dune.env)
 ordistep(mod1)
 ordistep(rda(dune ~ 1, dune.env), scope = formula(mod1))
+ordiR2step(rda(dune ~ 1, dune.env), mod1)
 }
 
 \keyword{ multivariate }



More information about the Vegan-commits mailing list