[Vegan-commits] r852 - in branches/1.15: R inst man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Jun 2 18:35:56 CEST 2009
Author: jarioksa
Date: 2009-06-02 18:35:56 +0200 (Tue, 02 Jun 2009)
New Revision: 852
Added:
branches/1.15/R/meandist.R
branches/1.15/R/plot.meandist.R
branches/1.15/R/print.summary.meandist.R
branches/1.15/R/summary.meandist.R
Modified:
branches/1.15/R/ordisurf.R
branches/1.15/R/tsallis.R
branches/1.15/inst/ChangeLog
branches/1.15/man/mrpp.Rd
branches/1.15/man/ordisurf.Rd
branches/1.15/man/tsallis.Rd
Log:
merged 818,820 (tsallis gained arg hill), 775, 789 (meandist), r796,801 (ordisurf can do linear and quadratic surfaces) to branches/1.15
Copied: branches/1.15/R/meandist.R (from rev 775, pkg/vegan/R/meandist.R)
===================================================================
--- branches/1.15/R/meandist.R (rev 0)
+++ branches/1.15/R/meandist.R 2009-06-02 16:35:56 UTC (rev 852)
@@ -0,0 +1,27 @@
+`meandist` <-
+ function(dist, grouping, ...)
+{
+ ## merge levels so that lower is always first (filling lower triangle)
+ mergenames <- function(X, Y, ...) {
+ xy <- cbind(X, Y)
+ xy <- apply(xy, 1, sort)
+ apply(xy, 2, paste, collapse = " ")
+ }
+ grouping <- factor(grouping, exclude = NULL)
+ cl <- outer(grouping, grouping, mergenames)
+ cl <- cl[lower.tri(cl)]
+ ## Cannot have within-group dissimilarity for group size 1
+ n <- table(grouping)
+ take <- matrix(TRUE, nlevels(grouping), nlevels(grouping))
+ diag(take) <- n > 1
+ take[upper.tri(take)] <- FALSE
+ ## Get output matrix
+ out <- matrix(NA, nlevels(grouping), nlevels(grouping))
+ out[take] <- tapply(dist, cl, mean)
+ out[upper.tri(out)] <- t(out)[upper.tri(out)]
+ rownames(out) <- colnames(out) <- levels(grouping)
+ class(out) <- c("meandist", "matrix")
+ attr(out, "n") <- table(grouping)
+ out
+}
+
Modified: branches/1.15/R/ordisurf.R
===================================================================
--- branches/1.15/R/ordisurf.R 2009-06-02 13:28:33 UTC (rev 851)
+++ branches/1.15/R/ordisurf.R 2009-06-02 16:35:56 UTC (rev 852)
@@ -14,7 +14,15 @@
X <- scores(x, choices = choices, display = display, ...)
x1 <- X[, 1]
x2 <- X[, 2]
- if (thinplate)
+ if (knots <= 0)
+ mod <- gam(y ~ x1 + x2, family = family, weights = w)
+ else if (knots == 1)
+ mod <- gam(y ~ poly(x1, 1) + poly(x2, 1),
+ family = family, weights = w)
+ else if (knots == 2)
+ mod <- gam(y ~ poly(x1, 2) + poly(x2, 2) + poly(x1, 1):poly(x2, 1),
+ family = family, weights = w)
+ else if (thinplate)
mod <- gam(y ~ s(x1, x2, k = knots), family = family,
weights = w)
else mod <- gam(y ~ s(x1, k = knots) + s(x2, k = knots),
Copied: branches/1.15/R/plot.meandist.R (from rev 775, pkg/vegan/R/plot.meandist.R)
===================================================================
--- branches/1.15/R/plot.meandist.R (rev 0)
+++ branches/1.15/R/plot.meandist.R 2009-06-02 16:35:56 UTC (rev 852)
@@ -0,0 +1,17 @@
+`plot.meandist` <-
+ function(x, cluster = "average", ...)
+{
+ n <- attr(x, "n")
+ cl <- hclust(as.dist(x), method = cluster, members = n)
+ cl <- as.dendrogram(cl, hang = 0)
+ w <- diag(x)[labels(cl)]
+ tr <- unlist(dendrapply(cl, function(n) attr(n, "height")))
+ root <- attr(cl, "height")
+ plot(cl, ylim = range(c(w, tr, root), na.rm = TRUE), leaflab = "none", ...)
+ for (i in 1:length(w)) segments(i, tr[i], i, w[i])
+ pos <- ifelse(w < tr, 1, 3)
+ pos[is.na(pos)] <- 1
+ w[is.na(w)] <- tr[is.na(w)]
+ text(1:length(w), w, labels = labels(cl), pos = pos, srt = 0)
+}
+
Copied: branches/1.15/R/print.summary.meandist.R (from rev 775, pkg/vegan/R/print.summary.meandist.R)
===================================================================
--- branches/1.15/R/print.summary.meandist.R (rev 0)
+++ branches/1.15/R/print.summary.meandist.R 2009-06-02 16:35:56 UTC (rev 852)
@@ -0,0 +1,19 @@
+`print.summary.meandist` <-
+ function(x, ...)
+{
+ cat("\nMean distances:\n")
+ tab <- rbind("within groups" = x$W,
+ "between groups" = x$B,
+ "overall" = x$D)
+ colnames(tab) <- "Average"
+ print(tab, ...)
+ cat("\nSummary statistics:\n")
+ tab <- rbind("MRPP A weights n" = x$A1,
+ "MRPP A weights n-1" = x$A2,
+ "MRPP A weights n(n-1)"= x$A3,
+ "Classification strength"=x$CS)
+ colnames(tab) <- "Statistic"
+ print(tab, ...)
+ invisible(x)
+}
+
Copied: branches/1.15/R/summary.meandist.R (from rev 775, pkg/vegan/R/summary.meandist.R)
===================================================================
--- branches/1.15/R/summary.meandist.R (rev 0)
+++ branches/1.15/R/summary.meandist.R 2009-06-02 16:35:56 UTC (rev 852)
@@ -0,0 +1,22 @@
+`summary.meandist` <-
+ function(object, ...)
+{
+ n <- attr(object, "n")
+ wmat <- n %o% n
+ diag(wmat) <- diag(wmat) - n
+ ## mean distances within, between groups and in total
+ W <- weighted.mean(diag(object), w = diag(wmat), na.rm = TRUE)
+ B <- weighted.mean(object[lower.tri(object)],
+ w = wmat[lower.tri(wmat)], na.rm = TRUE)
+ D <- weighted.mean(object, w = wmat, na.rm = TRUE)
+ ## Variants of MRPP statistics
+ A1 <- weighted.mean(diag(object), w = n, na.rm = TRUE)
+ A2 <- weighted.mean(diag(object), w = n - 1, na.rm = TRUE)
+ A3 <- weighted.mean(diag(object), w = n * (n - 1), na.rm = TRUE)
+ ##
+ out <- list(W = W, B = B, D = D, CS = B-W,
+ A1 = 1 - A1/D, A2 = 1 - A2/D, A3 = 1 - A3/D)
+ class(out) <- "summary.meandist"
+ out
+}
+
Modified: branches/1.15/R/tsallis.R
===================================================================
--- branches/1.15/R/tsallis.R 2009-06-02 13:28:33 UTC (rev 851)
+++ branches/1.15/R/tsallis.R 2009-06-02 16:35:56 UTC (rev 852)
@@ -1,6 +1,8 @@
tsallis <-
-function (x, scales = seq(0, 2, 0.2), norm=FALSE)
+function (x, scales = seq(0, 2, 0.2), norm=FALSE, hill=FALSE)
{
+ if (norm && hill)
+ stop("'norm = TRUE' and 'hill = TRUE' should not be used at the same time")
x <- as.matrix(x)
n <- nrow(x)
p <- ncol(x)
@@ -26,6 +28,13 @@
if (scales[a] == 1) result[, a] <- result[, a] / log(ST)
else result[, a] <- result[, a] / ((ST^(1-scales[a]) - 1) / (1 - scales[a]))
}
+ if (hill) {
+ result[, a] <- if (scales[a] == 1) {
+ exp(result[, a])
+ } else {
+ (1 - (scales[a] - 1) * result[, a])^(1/(1-scales[a]))
+ }
+ }
}
result <- as.data.frame(result)
if (any(dim(result) == 1))
Modified: branches/1.15/inst/ChangeLog
===================================================================
--- branches/1.15/inst/ChangeLog 2009-06-02 13:28:33 UTC (rev 851)
+++ branches/1.15/inst/ChangeLog 2009-06-02 16:35:56 UTC (rev 852)
@@ -5,6 +5,14 @@
Version 1.15-3 (opened 28 May, 2009)
+ * merged r818, 820: tsallis got argument hill.
+
+ * merged r796, 801: ordisurf can do linear or quadratic surfaces
+ with knots 2, 1 or 0.
+
+ * merged r775, r789: meandist with support function (sister
+ function to mrpp).
+
* merged r834, r837:r842 and r846 that put the observed test
statistic among permutations like should be done with order
statistics. Concerns functions anosim, mantel, mantel.partial,
Modified: branches/1.15/man/mrpp.Rd
===================================================================
--- branches/1.15/man/mrpp.Rd 2009-06-02 13:28:33 UTC (rev 851)
+++ branches/1.15/man/mrpp.Rd 2009-06-02 16:35:56 UTC (rev 852)
@@ -1,6 +1,10 @@
\name{mrpp}
\alias{mrpp}
\alias{print.mrpp}
+\alias{meandist}
+\alias{summary.meandist}
+\alias{print.summary.meandist}
+\alias{plot.meandist}
\title{ Multi Response Permutation Procedure of Within- versus
Among-Group Dissimilarities}
@@ -12,6 +16,9 @@
\usage{
mrpp(dat, grouping, permutations = 999, distance = "euclidean",
weight.type = 1, strata)
+meandist(dist, grouping, ...)
+\method{summary}{meandist}(object, ...)
+\method{plot}{meandist}(x, cluster = "average", ...)
}
\arguments{
@@ -29,6 +36,14 @@
\item{strata}{An integer vector or factor specifying the strata for
permutation. If supplied, observations are permuted only within the
specified strata.}
+ \item{dist}{A \code{\link{dist}} object of dissimilarities, such as
+ produced by functions \code{\link{dist}}, \code{\link{vegdist}} or
+ \code{\link{designdist}}.}.
+ \item{object, x}{A \code{meandist} result object.}
+ \item{cluster}{A clustering method for the \code{\link{hclust}}
+ function. Any \code{hclust} method can be used, but perhaps only
+ \code{"average"} and \code{"single"} make sense.}
+ \item{\dots}{Further arguments passed to functions.}
}
\details{ Multiple Response Permutation Procedure (MRPP) provides a test
@@ -75,8 +90,23 @@
\code{dat} as observations, and uses \code{\link{vegdist}} to find
the dissimilarities. The default \code{distance} is Euclidean as in the
traditional use of the method, but other dissimilarities in
-\code{\link{vegdist}} also are available. }
+\code{\link{vegdist}} also are available.
+Function \code{meandist} calculates a matrix of mean within-cluster
+dissimilarities (diagonal) and between-cluster dissimilarites
+(off-diagonal elements), and an attribute \code{n} of \code{grouping
+counts}. Function \code{summary} finds the within-class, between-class
+and overall means of these dissimilarities, and the MRPP statistics
+with all \code{weight.type} options and the classification
+strength. The function does not allow significance tests for these
+statistics, but you must use \code{mrpp} with appropriate
+\code{weight.type}. Function \code{plot} draws a dendrogram of the
+result matrix with given \code{cluster} method (see
+\code{\link{hclust}}). The terminal segments hang to within-cluster
+dissimilarity. If some of the clusters is more heterogeneous than the
+combined class, the leaf segment is reversed.
+}
+
\value{
The function returns a list of class mrpp with following items:
\item{call }{ Function call.}
@@ -154,6 +184,11 @@
expression(bold(delta)), cex=1.5 ) }
)
par(def.par)
+## meandist
+dune.md <- meandist(vegdist(dune), dune.env$Management)
+dune.md
+summary(dune.md)
+plot(dune.md)
}
\keyword{ multivariate }
\keyword{ nonparametric }
Modified: branches/1.15/man/ordisurf.Rd
===================================================================
--- branches/1.15/man/ordisurf.Rd 2009-06-02 13:28:33 UTC (rev 851)
+++ branches/1.15/man/ordisurf.Rd 2009-06-02 16:35:56 UTC (rev 852)
@@ -19,7 +19,10 @@
\item{y}{ Variable to be plotted. }
\item{choices}{Ordination axes. }
\item{knots}{Number of initial knots in \code{\link[mgcv]{gam}} (one
- more than degrees of freedom). }
+ more than degrees of freedom). If \code{knots = 0} or
+ \code{knots = 1} the function will fit a linear trend surface, and
+ if \code{knots = 2} the function will fit a quadratic trend surface
+ instead of a smooth surface. }
\item{family}{ Error distribution in \code{\link[mgcv]{gam}}. }
\item{col}{ Colour of contours. }
\item{thinplate}{Use thinplate splines in \code{\link[mgcv]{gam}}.}
Modified: branches/1.15/man/tsallis.Rd
===================================================================
--- branches/1.15/man/tsallis.Rd 2009-06-02 13:28:33 UTC (rev 851)
+++ branches/1.15/man/tsallis.Rd 2009-06-02 16:35:56 UTC (rev 852)
@@ -8,7 +8,7 @@
Function \code{tsallis} find Tsallis diversities with any scale or the corresponding evenness measures. Function \code{tsallisaccum} finds these statistics with accumulating sites.
}
\usage{
-tsallis(x, scales = seq(0, 2, 0.2), norm = FALSE)
+tsallis(x, scales = seq(0, 2, 0.2), norm = FALSE, hill = FALSE)
tsallisaccum(x, scales = seq(0, 2, 0.2), permutations = 100, raw = FALSE, ...)
\method{persp}{tsallisaccum}(x, theta = 220, phi = 15, col = heat.colors(100), zlim, ...)
}
@@ -17,6 +17,7 @@
\item{x}{Community data matrix or plotting object. }
\item{scales}{Scales of Tsallis diversity.}
\item{norm}{Logical, if \code{TRUE} diversity values are normalized by their maximum (diversity value at equiprobability conditions).}
+ \item{hill}{Calculate Hill numbers.}
\item{permutations}{Number of random permutations in accumulating sites.}
\item{raw}{If \code{FALSE} then return summary statistics of permutations, and if TRUE then returns the individual permutations.}
\item{theta, phi}{angles defining the viewing direction. \code{theta} gives the azimuthal direction and \code{phi} the colatitude.}
@@ -32,29 +33,36 @@
where \eqn{q} is a scale parameter, \eqn{S} the number of species in the sample (Tsallis 1988, Tothmeresz 1995). This diversity is concave for all \eqn{q>0}, but non-additive (Keylock 2005). For \eqn{q=0} it gives the number of species minus one, as \eqn{q} tends to 1 this gives Shannon diversity, for \eqn{q=2} this gives the Simpson index (see function \code{\link{diversity}}).
-When \code{norm = TRUE}, \code{tsallis} gives values normalized by the maximum:
+If \code{norm = TRUE}, \code{tsallis} gives values normalized by the maximum:
\deqn{H_q(max) = \frac{S^{1-q}-1}{1-q}}{H.q(max) = (S^(1-q)-1)/(1-q)}
where \eqn{S} is the number of species. As \eqn{q} tends to 1, maximum is defined as \eqn{ln(S)}.
+If \code{hill = TRUE}, \code{tsallis} gives Hill numbers (numbers equivalents, see Jost 2007):
+
+\deqn{D_q = (1-(q-1) H)^(1/(1-q))}{D.q = (1-(q-1)*H)^(1/(1-q))}
+
Details on plotting methods and accumulating values can be found on the help pages of the functions \code{\link{renyi}} and \code{\link{renyiaccum}}.
}
\value{
Function \code{tsallis} returns a data frame of selected indices. Function \code{tsallisaccum} with argument \code{raw = FALSE} returns a three-dimensional array, where the first dimension are the accumulated sites, second dimension are the diveristy scales, and third dimension are the summary statistics \code{mean}, \code{stdev}, \code{min}, \code{max}, \code{Qnt 0.025} and \code{Qnt 0.975}. With argument \code{raw = TRUE} the statistics on the third dimension are replaced with individual permutation results.
}
\references{
-Tsallis, C. (1988). Possible generalization of Boltzmann-Gibbs statistics.
+Tsallis, C. (1988) Possible generalization of Boltzmann-Gibbs statistics.
\emph{J. Stat. Phis.} 52, 479--487.
-Tothmeresz, B. (1995). Comparison of different methods for diversity
- ordering. \emph{Journal of Vegetation Science} 6, 283--290.
+Tothmeresz, B. (1995) Comparison of different methods for diversity
+ ordering. \emph{Journal of Vegetation Science} \bold{6}, 283--290.
-Patil, GP and Taillie, C. (1982). Diversity as a concep and its measurement.
- \emph{J. Am. Stat. Ass.} 77, 548--567.
+Patil, G. P. and Taillie, C. (1982) Diversity as a concep and its measurement.
+ \emph{J. Am. Stat. Ass.} \bold{77}, 548--567.
-Keylock, CJ (2005). Simpson diversity and the Shannon-Wiener index as special cases of a generalized entropy.
- \emph{Oikos} 109, 203--207.
+Keylock, C. J. (2005) Simpson diversity and the Shannon-Wiener index as special cases of a generalized entropy.
+ \emph{Oikos} \bold{109}, 203--207.
+
+Jost, L (2007) Partitioning diversity into independent alpha and beta components.
+ \emph{Ecology} \bold{88}, 2427--2439.
}
\author{\enc{P\'eter S\'olymos}{Peter Solymos}, \email{solymos at ualberta.ca}, based on the code of Roeland Kindt and Jari Oksanen written for \code{renyi}}
\seealso{
More information about the Vegan-commits
mailing list