[Vegan-commits] r950 - in branches/1.15: R inst man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Aug 25 20:41:00 CEST 2009
Author: jarioksa
Date: 2009-08-25 20:41:00 +0200 (Tue, 25 Aug 2009)
New Revision: 950
Added:
branches/1.15/R/mantel.correlog.R
branches/1.15/R/plot.mantel.correlog.R
branches/1.15/R/print.mantel.correlog.R
branches/1.15/man/mantel.correlog.Rd
Modified:
branches/1.15/inst/ChangeLog
Log:
copied mantel.correlog to branches/1.15
Copied: branches/1.15/R/mantel.correlog.R (from rev 949, pkg/vegan/R/mantel.correlog.R)
===================================================================
--- branches/1.15/R/mantel.correlog.R (rev 0)
+++ branches/1.15/R/mantel.correlog.R 2009-08-25 18:41:00 UTC (rev 950)
@@ -0,0 +1,153 @@
+`mantel.correlog` <-
+ function(D.eco, D.geo=NULL, XY=NULL, n.class=0, break.pts=NULL,
+ cutoff=TRUE, r.type="pearson", nperm=999, mult="holm",
+ progressive=TRUE)
+{
+ r.type <- match.arg(r.type, c("pearson", "spearman", "kendall"))
+ mult <- match.arg(mult, c("sidak", p.adjust.methods))
+
+ epsilon <- .Machine$double.eps
+ D.eco = as.matrix(D.eco)
+
+ ## Geographic distance matrix
+ if(!is.null(D.geo)) {
+ if(!is.null(XY))
+ stop("You provided both a geographic distance matrix and a list of site coordinates. Which one should the function use?")
+ D.geo <- as.matrix(D.geo)
+ } else {
+ if(is.null(XY)) {
+ stop("You did not provided a geographic distance matrix nor a list of site coordinates.")
+ } else {
+ D.geo <- as.matrix(dist(XY))
+ }
+ }
+
+ n <- nrow(D.geo)
+ if(n != nrow(D.eco))
+ stop("Numbers of objects in D.eco and D.geo are not equal.")
+ n.dist <- n*(n-1)/2
+ vec.D <- as.vector(as.dist(D.geo))
+ vec.DD <- as.vector(D.geo)
+
+ ## Number of classes
+ if(n.class > 0) {
+ if(!is.null(break.pts))
+ stop("You provided both a number of classes and a list of break points. Which one should the function use?")
+ } else {
+ if(is.null(break.pts)) {
+ ## Use Sturge's rule to determine the number of classes
+ n.class <- round(1 + log(n.dist, base=2))
+ start.pt <- min(vec.D)
+ end.pt <- max(vec.D)
+ width <- (end.pt - start.pt)/n.class
+ break.pts <- vector(length=n.class+1)
+ break.pts[n.class+1] <- end.pt
+ for(i in 1:n.class) {
+ break.pts[i] <- start.pt + width*(i-1) }
+ } else {
+ ## Use the list of break points
+ n.class <- length(break.pts) - 1
+ }
+ }
+ half.cl <- n.class %/% 2
+
+ ## Move the first breakpoint a little bit to the left
+ break.pts[1] = break.pts[1] - epsilon
+
+ ## Find the break points and the class indices
+ class.ind <- break.pts[1:n.class] +
+ (0.5*(break.pts[2:(n.class+1)]-break.pts[1:n.class]))
+
+ ## Create the matrix of distance classes
+ vec2 <- vector(length=n^2)
+ for(i in 1:n^2) {
+ vec2[i] <- min( which(break.pts >= vec.DD[i]) ) - 1
+ }
+
+ ## Start assembling the vectors of results
+ class.index <- NA
+ n.dist <- NA
+ mantel.r <- NA
+ mantel.p <- NA
+ ## check.sums = matrix(NA,n.class,1)
+
+ ## Create a model-matrix for each distance class, then compute a Mantel test
+ for(k in 1:n.class) {
+ class.index <- c(class.index, class.ind[k])
+ vec3 <- rep(0, n*n)
+ sel <- which(vec2 == k)
+ vec3[sel] <- 1
+ mat.D2 <- matrix(vec3,n,n)
+ diag(mat.D2) <- 0
+ n.dis <- sum(mat.D2)
+ n.dist <- c(n.dist, n.dis)
+ if(n.dis == 0) {
+ mantel.r <- c(mantel.r, NA)
+ mantel.p <- c(mantel.p, NA)
+ } else {
+ row.sums <- apply(mat.D2, 1, sum)
+ ## check.sums[k,1] = length(which(row.sums == 0))
+ if((cutoff==FALSE) ||
+ !(cutoff==TRUE && k > half.cl && any(row.sums == 0))) {
+ temp <- mantel(mat.D2, D.eco, method=r.type, permutations=nperm)
+ mantel.r <- c(mantel.r, -temp$statistic)
+ temp.p <- temp$signif
+ if(temp$statistic >= 0) {
+ temp.p <- ((temp.p*nperm)+1)/(nperm+1)
+ } else {
+ temp.p <- (((1-temp.p)*nperm)+1)/(nperm+1)
+ }
+ mantel.p <- c(mantel.p, temp.p)
+ } else {
+ mantel.r <- c(mantel.r, NA)
+ mantel.p <- c(mantel.p, NA)
+ }
+ }
+ }
+
+ mantel.res <- cbind(class.index, n.dist, mantel.r, mantel.p)
+ mantel.res <- mantel.res[-1,]
+
+ ## Note: vector 'mantel.p' starts with a NA value
+ mantel.p <- mantel.p[-1]
+ n.tests <- length(which(!is.na(mantel.p)))
+
+ if(mult == "none") {
+ colnames(mantel.res) <-
+ c("class.index", "n.dist", "Mantel.cor", "Pr(Mantel)")
+ } else {
+ ## Correct P-values for multiple testing
+ if(progressive) {
+ p.corr <- mantel.p[1]
+ if(mult == "sidak") {
+ for(j in 2:n.tests) {
+ p.corr <- c(p.corr, 1-(1-mantel.p[j])^j) }
+ } else {
+ for(j in 2:n.tests) {
+ temp <- p.adjust(mantel.p[1:j], method=mult)
+ p.corr <- c(p.corr, temp[j])
+ }
+ }
+ } else {
+ ## Correct all p-values for 'n.tests' simultaneous tests
+ if(mult == "sidak") {
+ p.corr <- 1 - (1 - mantel.p[1:n.tests])^n.tests
+ } else {
+ p.corr <- p.adjust(mantel.p[1:n.tests], method=mult)
+ }
+ }
+ temp <- c(p.corr, rep(NA,(n.class-n.tests)))
+ mantel.res <- cbind(mantel.res, temp)
+ colnames(mantel.res) <-
+ c("class.index", "n.dist", "Mantel.cor", "Pr(Mantel)", "Pr(corrected)")
+ }
+ rownames(mantel.res) <-
+ rownames(mantel.res,do.NULL = FALSE, prefix = "D.cl.")
+
+ ## Output the results
+ res <- list(mantel.res=mantel.res, n.class=n.class, break.pts=break.pts,
+ mult=mult, n.tests=n.tests, call=match.call() )
+ class(res) <- "mantel.correlog"
+ res
+}
+
Copied: branches/1.15/R/plot.mantel.correlog.R (from rev 949, pkg/vegan/R/plot.mantel.correlog.R)
===================================================================
--- branches/1.15/R/plot.mantel.correlog.R (rev 0)
+++ branches/1.15/R/plot.mantel.correlog.R 2009-08-25 18:41:00 UTC (rev 950)
@@ -0,0 +1,17 @@
+`plot.mantel.correlog` <-
+ function(x, alpha=0.05, ...)
+{
+ lim <- max(x$n.tests)
+ plot(x$mantel.res[1:lim,1],x$mantel.res[1:lim,3],
+ xlab="Distance class index", ylab="Mantel correlation", pch=22)
+ if(x$mult == "none") {
+ signif <- which((x$mantel.res[1:lim,4] <= alpha))
+ } else {
+ signif <- which((x$mantel.res[1:lim,5] <= alpha))
+ }
+ lines(x$mantel.res[1:lim,1], x$mantel.res[1:lim,3])
+ points(x$mantel.res[1:lim,1], x$mantel.res[1:lim,3], pch=22, bg="white")
+ points(x$mantel.res[signif,1], x$mantel.res[signif,3], pch=22, bg="black")
+ abline(a=0, b=0, col="red")
+ invisible()
+}
Copied: branches/1.15/R/print.mantel.correlog.R (from rev 949, pkg/vegan/R/print.mantel.correlog.R)
===================================================================
--- branches/1.15/R/print.mantel.correlog.R (rev 0)
+++ branches/1.15/R/print.mantel.correlog.R 2009-08-25 18:41:00 UTC (rev 950)
@@ -0,0 +1,9 @@
+'print.mantel.correlog' <- function(x, ...)
+{
+ cat('\nMantel Correlogram Analysis\n')
+ cat('\nCall:\n','\n')
+ cat(deparse(x$call),'\n')
+ cat('\n')
+ printCoefmat(x$mantel.res, P.values=TRUE, signif.stars=TRUE, Pvalues = TRUE)
+ invisible(x)
+}
\ No newline at end of file
Modified: branches/1.15/inst/ChangeLog
===================================================================
--- branches/1.15/inst/ChangeLog 2009-08-25 18:24:54 UTC (rev 949)
+++ branches/1.15/inst/ChangeLog 2009-08-25 18:41:00 UTC (rev 950)
@@ -43,6 +43,8 @@
* copied indpower (with rev948 for beals.Rd).
+ * copied mantel.correlog.
+
Version 1.15-3 (released June 17, 2009)
* merged revs 866 to 868: changed the way capscale displays
Copied: branches/1.15/man/mantel.correlog.Rd (from rev 949, pkg/vegan/man/mantel.correlog.Rd)
===================================================================
--- branches/1.15/man/mantel.correlog.Rd (rev 0)
+++ branches/1.15/man/mantel.correlog.Rd 2009-08-25 18:41:00 UTC (rev 950)
@@ -0,0 +1,159 @@
+\name{mantel.correlog}
+\alias{mantel.correlog}
+\alias{print.mantel.correlog}
+\alias{plot.mantel.correlog}
+\title{ Mantel Correlogram }
+
+\description{
+ Function \code{mantel.correlog} computes a multivariate
+ Mantel correlogram. Proposed by Sokal (1986) and Oden and Sokal
+ (1986), the method is also described in Legendre and Legendre (1998,
+ pp. 736-738).
+}
+
+\usage{
+mantel.correlog(D.eco, D.geo=NULL, XY=NULL, n.class=0, break.pts=NULL,
+cutoff=TRUE, r.type="pearson", nperm=999, mult="holm", progressive=TRUE)
+\method{plot}{mantel.correlog}(x, alpha=0.05, ...)
+}
+
+\arguments{
+ \item{D.eco}{ An ecological distance matrix, with class
+ either \code{dist} or \code{matrix}. }
+
+ \item{D.geo}{ A geographic distance matrix, with class either
+ \code{dist} or \code{matrix}. Provide either \code{D.geo} or
+ \code{XY}. Default: \code{D.geo=NULL}. }
+
+ \item{XY}{ A file of Cartesian geographic coordinates of the
+ points. Default: \code{XY=NULL}. }
+
+ \item{n.class}{ Number of classes. If \code{n.class=0}, the Sturge
+ equation will be used unless break points are provided. }
+
+ \item{break.pts}{ Vector containing the break points of the distance
+ distribution. Default: \code{break.pts=NULL}. }
+
+ \item{cutoff}{ For the second half of the distance classes,
+ \code{cutoff = TRUE} limits the correlogram to the distance classes
+ that include all points. If \code{cutoff = FALSE}, the correlogram
+ includes all distance classes. }
+
+ \item{r.type}{ Type of correlation in calculation of the Mantel
+ statistic. Default: \code{r.type="pearson"}. Other choices are
+ \code{r.type="spearman"} and \code{r.type="kendall"}, as in functions
+ \code{\link{cor}} and \code{\link{mantel}}. }
+
+ \item{nperm}{ Number of permutations for the tests of
+ significance. Default: \code{nperm=999}. For large data files,
+ permutation tests are rather slow. }
+
+ \item{mult}{ Correct P-values for multiple testing. The correction
+ methods are \code{"holm"} (default), \code{"hochberg"},
+ \code{"sidak"}, and other methods available in the
+ \code{\link{p.adjust}} function: \code{"bonferroni"} (best known, but
+ not recommended because it is overly conservative), \code{"hommel"},
+ \code{"BH"}, \code{"BY"}, \code{"fdr"}, and \code{"none"}. }
+
+ \item{progressive}{ Default: \code{progressive=TRUE} for progressive
+ correction of multiple-testing, as described in Legendre and Legendre
+ (1998, p. 721). Test of the first distance class: no correction;
+ second distance class: correct for 2 simultaneous tests; distance
+ class k: correct for k simultaneous tests. \code{progressive=FALSE}:
+ correct all tests for \code{n.class} simultaneous tests. }
+
+ \item{x}{ Output of \code{mantel.correlog}. }
+
+ \item{alpha}{ Significance level for the points drawn with black
+ symbols in the correlogram. Default: \code{alpha=0.05}. }
+
+ \item{...}{ Other parameters passed from other functions. }
+}
+
+\details{ A correlogram is a graph in which spatial correlation values
+ are plotted, on the ordinate, as a function of the geographic distance
+ classes among the study sites along the abscissa. In a Mantel
+ correlogram, a Mantel correlation (Mantel 1967) is computed between a
+ multivariate (e.g. multi-species) distance matrix of the user's choice
+ and a design matrix representing each of the geographic distance
+ classes in turn. The Mantel statistic is tested through a
+ permutational Mantel test performed by \code{vegan}'s
+ \code{\link{mantel}} function.
+
+ When a correction for multiple testing is applied, more permutations
+ are necessary than in the no-correction case, to obtain significant
+ p-values in the higher correlogram classes.
+
+ The \code{print.mantel.correlog} function prints out the
+ correlogram. See examples. }
+
+\value{
+
+ \item{mantel.res }{A table with the distance classes as rows and the
+ class indices, number of distances per class, Mantel statistics
+ (computed using Pearson's r, Spearman's r, or Kendall's tau), and
+ p-values as columns. A positive Mantel statistic indicates positive
+ spatial correlation. An additional column with p-values corrected for
+ multiple testing is added unless \code{mult="none"}. }
+
+ \item{n.class }{The n umber of distance classes. }
+
+ \item{break.pts }{The break points provided by the user or computed by
+ the program. }
+
+ \item{mult }{The name of the correction for multiple testing. No
+ correction: \code{mult="none"}. } #
+
+ \item{progressive }{A logical (\code{TRUE}, \code{FALSE}) value
+ indicating whether or not a progressive correction for multiple
+ testing was requested. } \item{n.tests }{The number of distance
+ classes for which Mantel tests have been computed and tested for
+ significance. }
+
+\item{call }{The function call. } }
+
+\author{ Pierre Legendre, Universite de Montreal }
+
+\references{
+
+ Legendre, P. and L. Legendre. 1998. Numerical ecology, 2nd English
+ edition. Elsevier Science BV, Amsterdam.
+
+ Mantel, N. 1967. The detection of disease clustering and a generalized
+ regression approach. Cancer Res. 27: 209-220.
+
+ Oden, N. L. and R. R. Sokal. 1986. Directional autocorrelation: an
+ extension of spatial correlograms to two dimensions. Syst. Zool. 35:
+ 608-617.
+
+ Sokal, R. R. 1986. Spatial data analysis and historical
+ processes. 29-43 in: E. Diday et al. [eds.] Data analysis and
+ informatics, IV. North-Holland, Amsterdam. }
+
+\examples{
+# Mite data from "vegan"
+data(mite)
+data(mite.xy)
+mite.hel <- decostand(mite, "hellinger")
+mite.hel. <- dist(mite.hel)
+
+mite.correlog <- mantel.correlog(mite.hel.D, XY=mite.xy, nperm=99)
+summary(mite.correlog)
+mite.correlog
+plot(mite.correlog)
+
+mite.correlog2 <- mantel.correlog(mite.hel.D, XY=mite.xy, cutoff=FALSE,
+r.type="spearman", nperm=99)
+summary(mite.correlog2)
+mite.correlog2
+plot(mite.correlog2)
+
+## Mite correlogram after spatially detrending the mite data
+mite.h.det <- resid(lm(as.matrix(mite.hel.D) ~ ., data=mite.xy))
+mite.correlog3 <- mantel.correlog(mite.h.det, XY=mite.xy, nperm=99)
+mite.correlog3
+plot(mite.correlog3)
+
+}
+
+\keyword{ multivariate }
\ No newline at end of file
More information about the Vegan-commits
mailing list