[Analogue-commits] r109 - in pkg: . R inst man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Apr 13 22:40:25 CEST 2009
Author: gsimpson
Date: 2009-04-13 22:40:25 +0200 (Mon, 13 Apr 2009)
New Revision: 109
Added:
pkg/R/print.residLen.R
pkg/R/residLen.R
pkg/man/residLen.Rd
Modified:
pkg/DESCRIPTION
pkg/inst/ChangeLog
Log:
add squared residual length code
Modified: pkg/DESCRIPTION
===================================================================
--- pkg/DESCRIPTION 2009-04-06 17:54:54 UTC (rev 108)
+++ pkg/DESCRIPTION 2009-04-13 20:40:25 UTC (rev 109)
@@ -1,7 +1,7 @@
Package: analogue
Type: Package
Title: Analogue and weighted averaging methods for palaeoecology
-Version: 0.6-7
+Version: 0.6-8
Date: $Date$
Depends: R (>= 2.5.0), stats, graphics, vegan, lattice, MASS
Author: Gavin L. Simpson, Jari Oksanen
Added: pkg/R/print.residLen.R
===================================================================
--- pkg/R/print.residLen.R (rev 0)
+++ pkg/R/print.residLen.R 2009-04-13 20:40:25 UTC (rev 109)
@@ -0,0 +1,15 @@
+## print method
+`print.residLen` <- function(x,
+ digits = min(3, getOption("digits") - 4),
+ probs = c(0.5, 0.75, 0.9, 0.95, 0.99), ...) {
+ cat("\n")
+ writeLines(strwrap("Squared residual lengths",
+ prefix = "\t"))
+ cat("\nCall:\n")
+ cat(paste(deparse(x$call), "\n\n"))
+ cat("Quantiles of training set lengths:\n")
+ print(with(x, quantile(train, probs = probs)), digits = digits)
+ cat("\nQuantiles of passive sample lengths:\n")
+ print(with(x, quantile(passive, probs = probs)), digits = digits)
+ cat("\n")
+}
Added: pkg/R/residLen.R
===================================================================
--- pkg/R/residLen.R (rev 0)
+++ pkg/R/residLen.R 2009-04-13 20:40:25 UTC (rev 109)
@@ -0,0 +1,82 @@
+## compute the squared residual length statistic for a constrained
+## ordination given a single constraint
+
+`residLen` <- function(train, ...) {
+ UseMethod("residLen")
+}
+
+`residLen.default` <- function(train, env, passive,
+ method = c("cca", "rda"), ...) {
+ ## inline functions
+ ## fitted values for passive samples
+ fittedSuppl <- function(ord, newdata, csum) {
+ #########################################################
+ ## Fitted values for passive samples
+ ## Arguments:
+ ## ord = vegan ordination object (CCA/RDA)
+ ## newdata = matrix of passive species data with same
+ ## columns as that used to fit ord
+ ## csum = colum sums for the training data
+ #########################################################
+ ## species scores
+ b <- predict(ord, type = "sp")
+ ## site scores
+ xi <- predict(ord, newdata = newdata, type = "wa")
+ fik <- sweep(1 + (xi %*% t(b)), 2, csum / sum(csum),
+ "*")
+ ## fitted values
+ fik <- sweep(fik, 1, rowSums(newdata), "*")
+ }
+ sqrl.uni <- function(X, colsum, fik) {
+ #########################################################
+ ## Squared residual length for unimodal methods
+ ## Arguments:
+ ## X = species data (training or passive) to which
+ ## we want to find the residual length
+ ## colsum = column sums for the training species data
+ ## fik = fitted species data for either the training
+ ## or passive samples
+ #########################################################
+ yip <- rowSums(X)
+ ypk <- colsum
+ ypp <- sum(ypk)
+ A <- (sweep(X, 2, ypk, "/") - sweep(fik, 2, ypk, "/"))^2
+ B <- sweep(A, 2, ypk, "*")
+ C <- rowSums(B / ypp)
+ res <- (1 / (yip / ypp))^2 * C
+ return(res)
+ }
+ ## merge train and passive
+ dat <- join(train, passive)
+ train <- dat[[1]]
+ passive <- dat[[2]]
+ ## check env is same length as nrow(train)
+ if(!isTRUE(all.equal(length(env), nrow(train))))
+ stop("'train' and 'env' imply different numbers of observations")
+ ## ordinate
+ if(missing(method))
+ method <- "cca"
+ method <- match.arg(method)
+ ##if(method == "rda")
+ ## .NotYetUsed("method = \"rda\" not currently implemented.")
+ FUN <- match.fun(method)
+ ## ordinate
+ ord <- FUN(X = train, Y = env)
+ ## fitted values
+ fit <- fitted(ord, type = "response")
+ ## colSums of train
+ ypk <- colSums(train)
+ ## predict locations for the passive
+ pred <- fittedSuppl(ord, passive, ypk)
+ ## system.call
+ .call <- match.call()
+ .call[[1]] <- as.name("residLen")
+ ## residual lengths for train
+ res <- list(train = sqrl.uni(train, ypk, fit),
+ passive = sqrl.uni(passive, ypk, pred),
+ ordination = ord,
+ call = .call)
+ class(res) <- "residLen"
+ attr(res, "method") <- method
+ return(res)
+}
Modified: pkg/inst/ChangeLog
===================================================================
--- pkg/inst/ChangeLog 2009-04-06 17:54:54 UTC (rev 108)
+++ pkg/inst/ChangeLog 2009-04-13 20:40:25 UTC (rev 109)
@@ -1,7 +1,14 @@
analogue Change Log
-Version 0.6-7 (Opened Fri 3 Apr 2009)
+Version 0.6-8 (Opened Mon 13 Apr 2009)
+ * residLen: new function to compute squared residual length
+ diagnostic for passive samples in a constrained ordination. Used as
+ a test of whether core samples are well fitted in a transfer function
+ model.
+
+Version 0.6-7 (Closed Mon 13 Apr 2009)
+
* optima, tolerance: new methods to coerce objects of these classes
to data frames.
Added: pkg/man/residLen.Rd
===================================================================
--- pkg/man/residLen.Rd (rev 0)
+++ pkg/man/residLen.Rd 2009-04-13 20:40:25 UTC (rev 109)
@@ -0,0 +1,61 @@
+\name{residLen}
+\alias{residLen}
+\alias{residLen.default}
+\alias{print.residLen}
+\title{Squared residual length diagnostics}
+\description{
+ The squared residual length between the fitted values of a constrained
+ ordination and the original species data is one diagnostic for
+ transfer function models.
+}
+\usage{
+residLen(train, ...)
+
+\method{residLen}{default}(train, env, passive,
+ method = c("cca","rda"), \dots)
+}
+%- maybe also 'usage' for other objects documented here.
+\arguments{
+ \item{train}{data frame; the training set species data.}
+ \item{env}{vector; the training set environmental data.}
+ \item{passive}{data frame; the passive samples species data.}
+ \item{method}{the ordination technique to use. Only \code{"cca"} is
+ currently supported.}
+ \item{\dots}{additional arguments passed to methods.}
+}
+\details{
+ The squared residual lengths are computed for the training set samples
+ and the passive samples separately. Passive samples that are poorly
+ fitted in the transfer function model will have large squared residual
+ distances between the observed species data and the fitted values from
+ the constrained ordination.
+}
+\value{
+ \code{residLen} returns an object of class \code{"residLen"} with the
+ attribute \code{"method"} set to \code{"method"}. This object is a
+ list with the following components:
+
+ \item{train, passive}{The squared residual lengths for the training
+ set and the passive samples.}
+ \item{ordination}{The fitted ordination.}
+ \item{call}{The matched call.}
+}
+\references{
+ Ter Braak C.J.F. and Smilauer P. (2002) CANOCO Reference manual and
+ CanoDraw for Windows User's guide: Software for Canonical Ordination
+ (version 4.5). Microcomputer Power (Ithaca, NY, USA), 500pp.
+}
+\author{Gavin L. Simpson}
+\seealso{
+ \code{\link[vegan]{cca}} and \code{\link[vegan]{predict.cca}} for some
+ of the underlying computations.
+}
+\examples{
+data(swapdiat, swappH, rlgh)
+
+## squared residual lengths for RLGH
+rlens <- residLen(swapdiat, swappH, rlgh)
+rlens
+}
+\keyword{methods}
+\keyword{multivariate}
More information about the Analogue-commits
mailing list