[Vegan-commits] r2529 - in pkg: . ordiconsensus ordiconsensus/R ordiconsensus/data ordiconsensus/man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Jun 17 15:44:09 CEST 2013


Author: gblanchet
Date: 2013-06-17 15:44:09 +0200 (Mon, 17 Jun 2013)
New Revision: 2529

Added:
   pkg/ordiconsensus/
   pkg/ordiconsensus/DESCRIPTION
   pkg/ordiconsensus/NAMESPACE
   pkg/ordiconsensus/R/
   pkg/ordiconsensus/R/RV.R
   pkg/ordiconsensus/R/coeffCompare.R
   pkg/ordiconsensus/R/consensusRDA.R
   pkg/ordiconsensus/data/
   pkg/ordiconsensus/data/beetle.expl.rda
   pkg/ordiconsensus/data/beetle.rda
   pkg/ordiconsensus/man/
   pkg/ordiconsensus/man/RV.Rd
   pkg/ordiconsensus/man/beetle.Rd
   pkg/ordiconsensus/man/coeffCompare.Rd
   pkg/ordiconsensus/man/consensusRDA.Rd
   pkg/ordiconsensus/man/ordiconsensus-package.Rd
Log:
added the package ordiconsensus

Added: pkg/ordiconsensus/DESCRIPTION
===================================================================
--- pkg/ordiconsensus/DESCRIPTION	                        (rev 0)
+++ pkg/ordiconsensus/DESCRIPTION	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,11 @@
+Package: ordiconsensus
+Type: Package
+Title: Consensus of canonical ordinations through the canonical redundancy analysis
+Version: 0.3-1
+Date: 2012-11-09
+Author: F. Guillaume Blanchet
+Maintainer: F. Guillaume Blanchet <guillaume.blanchet at helsinki.fi>
+Description: This package include function to calculate a consensus of canonical redundancy analyses performed using different association coefficients
+Depends: vegan
+Suggests: FactoMineR
+License: Unlimited

Added: pkg/ordiconsensus/NAMESPACE
===================================================================
--- pkg/ordiconsensus/NAMESPACE	                        (rev 0)
+++ pkg/ordiconsensus/NAMESPACE	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,8 @@
+### Export
+
+export(coeffCompare,consensusRDA,RV)
+
+### Import
+
+import(stats)
+import(vegan)

Added: pkg/ordiconsensus/R/RV.R
===================================================================
--- pkg/ordiconsensus/R/RV.R	                        (rev 0)
+++ pkg/ordiconsensus/R/RV.R	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,16 @@
+RV <-
+function(X, Y) {
+	#CC# General checks
+	if (nrow(X) != nrow(Y)) stop("'X' needs to have the same number of rows as 'Y'")
+	if (nrow(X) == 1) stop("Impossible to calculate RV using 1 object")
+	
+	Y <- scale(Y, scale = FALSE)
+	X <- scale(X, scale = FALSE)
+	
+	XXt <- tcrossprod(X)
+	YYt <- tcrossprod(Y)
+	
+	rv <- sum(diag(XXt %*% YYt))/(sum(diag(XXt %*% XXt)) * sum(diag(YYt %*% YYt)))^0.5
+	
+	return(rv)
+}

Added: pkg/ordiconsensus/R/coeffCompare.R
===================================================================
--- pkg/ordiconsensus/R/coeffCompare.R	                        (rev 0)
+++ pkg/ordiconsensus/R/coeffCompare.R	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,122 @@
+coeffCompare <-
+function(ordires, ordisigniaxis,pval=0.05){
+###
+### Compare a series of association coefficient calculated through an RDA
+###
+### Arguments :
+###
+### ordires : A list of "cca rda" object from the vegan package gathering a series of RDA performed with different association coefficients on the same data.
+### ordisigniaxis : A list of anova.cca object where each axis was tested for each RDA in ordires or a vector defining the number of significant axes in the RDA. See details 
+### pval : Numeric. P-value threshold to select the number of axes to use. This argument is only active if a list of anova.cca object is given for the argument ordisigniaxis, otherwise it is not considered.
+###
+###
+### Details :
+### 
+### For the argument ordisigniaxis, if a vector of number of significant axes is given, for each RDA, it is assumed that the significant axes are selected in sequential order from the first axis.  
+### 
+### Value :
+###
+### RVmat : A matrix of RV coefficients calculated for all pairs of association coefficients
+### mst : minimum spanning tree calculated on (1-RVmat)
+###
+### F. Guillaume Blanchet - February 2012. (Modified November 2012, June 2013)
+################################################################################
+	#----------------
+	#CC# Basic object
+	#----------------
+	#CC# Number of sites
+	nsites<-nrow(scores(ordires[[1]],display="sites"))
+	#CC# Number of juges (association coefficients)
+	njuges<-length(ordires)
+	
+	#------------------
+	#CC# General checks
+	#------------------
+	#### Check if ordires contains only RDAs
+	allRDA<-sapply(ordires,function(x) any(class(x)=="rda"))
+	if(!all(allRDA)){
+		stop("One or more canonical ordination in 'ordires' is not an RDA")
+	}
+	
+	#### Check if ordires and ordisigniaxis have the same number of components
+	if(length(ordisigniaxis)!=length(ordires)){
+		stop("'ordires' is not the same length as 'ordisigniaxis'")
+	}
+	
+	ordisigniClass<-sapply(ordisigniaxis,function(x) class(x))
+	
+	#### Check if the capscale objects have the right number of species
+	anycapscale<-which(sapply(ordires,function(x) any(class(x)=="capscale")))
+	if(length(anycapscale) > 0){
+		nspcapscale<-numeric()
+		counter<-1
+		for(i in anycapscale){
+			nspcapscale[counter]<-nrow(scores(ordires[[i]],display="sp"))
+			counter<-counter+1
+		}
+		if(any(nspcapscale==nsites)){
+			stop("One or more of the analysis performed with capscale() did not include a site by species community matrix")
+		}
+	}
+	
+	#### If ordisigniaxis is a list
+	if(is.list(ordisigniaxis)){
+		#### Check P-values
+		if(pval < 0 | pval > 1){
+			stop("'pval' must range between 0 and 1")
+		}
+	
+		allsigni<-sapply(ordisigniaxis,function(x) any(class(x)=="anova.cca"))
+		if(!all(allsigni)){
+			stop("One or more canonical ordination test in 'ordisigniaxis' is not an 'anova.cca' object")
+		}
+		
+		anovaname<-unique(unlist(strsplit(sapply(ordisigniaxis,function(x) rownames(x)[1]),"1")))
+		if(!all(anovaname == "RDA" | anovaname == "CAP")){
+			stop("anova.cca by axis should be either 'RDA' or 'CAP'")
+		}
+		
+		#CC# Extract the number of signicant axes to use for each canonical ordinations
+		ordisigniaxis<-sapply(ordisigniaxis,function(x) length(which(x[,5]<=pval)))
+	#### If ordisigniaxis is a vector
+	}else{
+		if(is.vector(ordisigniaxis)){
+			if(!is.numeric(ordisigniaxis)){
+				stop("'ordisigniaxis' should be numeric")
+			}
+		}else{
+			stop("'ordisigniaxis' should either be a list of anova.cca objects or a vector of number of significant axes")
+		}
+	}
+	
+	#### If there are no significant axes
+	if(any(ordisigniaxis < 1)){
+		stop("One or more analysis does not have any significant axis remove it and start again")
+	}
+	
+	#-------------------------------------------------------------
+	#CC# Extract all the significant axes of matrix Z in scaling 1
+	#-------------------------------------------------------------
+	Zsigni<-vector("list",length=njuges)
+	
+	for(i in 1:njuges){
+		Zsigni[[i]]<-scores(ordires[[i]],display="lc",choices=1:ordisigniaxis[i],scaling=1)
+	}
+	
+	#CC# Coefficient RVf between the different association ceofficients
+	RVassoCoeff<-matrix(NA,nrow=njuges,ncol=njuges)
+	
+	for(i in 1:njuges){
+		for(j in 1:njuges){
+			RVassoCoeff[i,j]<-RV(Zsigni[[i]],Zsigni[[j]])
+		}
+	}
+	
+	#CC# Construct a minimum spanning tree
+	mst<-spantree(as.dist(1-RVassoCoeff))
+	
+	#CC# results
+	res<-list(RVmat=RVassoCoeff,mst=mst)
+	class(res)<-"coeffCompare"
+	return(res)
+}

Added: pkg/ordiconsensus/R/consensusRDA.R
===================================================================
--- pkg/ordiconsensus/R/consensusRDA.R	                        (rev 0)
+++ pkg/ordiconsensus/R/consensusRDA.R	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,162 @@
+consensusRDA <-
+function(ordires, ordisigniaxis,resp.var,expl.var,pval=0.05,scaling=2){
+###
+### Calculates a consensus results for a series of canonical ordinations performed with different association coefficients on the same data. This function is only available for ordinations performed with canonical redundancy analysis (RDA) or variant of RDA such as distance-based RDA and transformation-based RDA
+###
+### Arguments :
+###
+### ordires : A list of "cca rda" object from the vegan package gathering a series of RDA performed with different association coefficients on the same data.
+### ordisigniaxis : A list of anova.cca object where each axis was tested for each RDA in ordires or a vector defining the number of significant axes in the RDA. See details. 
+### resp.var : Matrix of response variables
+### expl.var : Matrix of explanatory variables
+### pval : Numeric. P-value threshold to select the number of axes to use. This argument is only active if a list of anova.cca object is given for the argument ordisigniaxis, otherwise it is not considered.
+### scaling : Type of scaling used to project the results. Default is 1 (distance).
+###
+###
+### Details
+### 
+### For the argument ordisigniaxis, if a vector of number of significant axes is given, for each RDA, it is assumed that the significant axes are selected in sequential order from the first axis.  
+### 
+### Although it is possible to apply a scaling 3 to the RDA (it is available in the vegan package), this scaling should only be used for canonical correspondence analysis (CCA), it does not make any sense to use in the RDA framework.
+###
+### F. Guillaume Blanchet - February 2012 (Modified November 2012, June 2013)
+################################################################################
+	#----------------
+	#CC# Basic object
+	#----------------
+	#CC# Number of sites
+	nsites<-nrow(resp.var)
+	#CC# Number of species
+	nsp<-ncol(resp.var)
+	#CC# Number of juges (association coefficients)
+	njuges<-length(ordires)
+	
+	#------------------
+	#CC# General checks
+	#------------------
+	#### Check if there are completely collinear explanatory variables and count the number of non-collinear variables
+	expl.var.tmp<-expl.var
+	deter<-det(cor(expl.var))
+	
+	#CC# threshold
+	tresh<-10^-8
+	
+	while(deter<= tresh){
+		expl.var.tmp<-expl.var.tmp[,-1]	
+		deter<-det(cor(expl.var.tmp))
+	}
+	nNoncolldesc<-ncol(expl.var.tmp)
+	
+	#### Check if ordires contains only RDAs
+	allRDA<-sapply(ordires,function(x) any(class(x)=="rda"))
+	if(!all(allRDA)){
+		stop("One or more canonical ordination in 'ordires' is not an RDA")
+	}
+	
+	#### Check if ordires and ordisigniaxis have the same number of components
+	if(length(ordisigniaxis)!=length(ordires)){
+		stop("'ordires' is not the same length as 'ordisigniaxis'")
+	}
+	
+	ordisigniClass<-sapply(ordisigniaxis,function(x) class(x))
+	
+	#### Check if the capscale objects have the right number of species
+	anycapscale<-which(sapply(ordires,function(x) any(class(x)=="capscale")))
+	if(length(anycapscale) > 0){
+		nspcapscale<-numeric()
+		counter<-1
+		for(i in anycapscale){
+			nspcapscale[counter]<-nrow(scores(ordires[[i]],display="sp"))
+			counter<-counter+1
+		}
+		if(any(nspcapscale==nsites)){
+			stop("One or more of the analysis performed with capscale() did not include a site by species community matrix")
+		}
+	}
+	
+	#### If ordisigniaxis is a list
+	if(is.list(ordisigniaxis)){
+		#### Check P-values
+		if(pval < 0 | pval > 1){
+			stop("'pval' must range between 0 and 1")
+		}
+
+		allsigni<-sapply(ordisigniaxis,function(x) any(class(x)=="anova.cca"))
+		if(!all(allsigni)){
+			stop("One or more canonical ordination test in 'ordisigniaxis' is not an 'anova.cca' object")
+		}
+		
+		anovaname<-unique(unlist(strsplit(sapply(ordisigniaxis,function(x) rownames(x)[1]),"1")))
+		if(!all(anovaname == "RDA" | anovaname == "CAP")){
+			stop("anova.cca by axis should be either 'RDA' or 'CAP'")
+		}
+		
+		#CC# Extract the number of signicant axes to use for each canonical ordinations
+		ordisigniaxis<-sapply(ordisigniaxis,function(x) length(which(x[,5]<=pval)))
+	#### If ordisigniaxis is a vector
+	}else{
+		if(is.vector(ordisigniaxis)){
+			if(!is.numeric(ordisigniaxis)){
+				stop("'ordisigniaxis' should be numeric")
+			}
+		}else{
+			stop("'ordisigniaxis' should either be a list of anova.cca objects or a vector of number of significant axes")
+		}
+	}
+	
+	#### If there are no significant axes
+	if(any(ordisigniaxis < 1)){
+		stop("One or more analysis does not have any significant axis remove it and start again")
+	}
+	
+	#### Make sure that X is a matrix 
+	expl.var<-as.matrix(expl.var)
+		
+	#-------------------------------------------------------------
+	#CC# Extract all the significant axes of matrix Z in scaling 1
+	#-------------------------------------------------------------
+	Zsigni<-vector("list",length=njuges)
+	Zsignimat<-matrix(NA,ncol=0,nrow=nsites)
+	
+	for(i in 1:njuges){
+		Zsigni[[i]]<-scores(ordires[[i]],display="lc",choices=1:ordisigniaxis[i],scaling=1)
+		Zsignimat<-cbind(Zsignimat,Zsigni[[i]])
+	}
+	
+	#-------------------------------------------------
+	#CC# Perform an RDA between Zsignimat and expl.var
+	#-------------------------------------------------
+	Zrda<-rda(Zsignimat,expl.var)
+	ZrdaEigen<-eigenvals(Zrda)[1:nNoncolldesc]
+	naxes<-length(ZrdaEigen)
+	
+	#CC# Extract the Z consensus results from RDA
+	Zconsensus<-scores(Zrda,choice=1:nNoncolldesc,display="lc",scaling=1)
+	
+	#CC# Extract the C consensus results from RDA
+	Cconsensus<-scores(Zrda,choice=1:nNoncolldesc,display="bp",scaling=1)
+	
+	#### This is the procedure proposed by Legendre and Legendre (2012, Subsection 9.3.3)
+	#### It is also the procedure proposed by Oksanen et al. in vegan and described in Blanchet et al. (submitted)
+	Uconsensus<-t(scale(resp.var,scale=FALSE))%*%Zconsensus%*%diag(ZrdaEigen^(-0.5))/sqrt(nsites-1)
+	
+	#----------------------
+	#CC# Ordination Scaling
+	#----------------------
+	if(scaling==2){
+		#-------------
+		#CC# Consensus
+		#-------------
+		#CC# Sites
+		Zconsensus<-sweep(Zconsensus,2,sqrt(ZrdaEigen/sum(ZrdaEigen)),"/")
+		
+		#CC# Species
+		Uconsensus<-sweep(Uconsensus,2,sqrt(ZrdaEigen/sum(ZrdaEigen)),"*")
+		
+		#### No scaling is performed on descriptors
+	}
+	
+	res<-list(values=ZrdaEigen,siteConsensus=Zconsensus,spConsensus=Uconsensus,descConsensus=Cconsensus)
+	
+	return(res)
+}

Added: pkg/ordiconsensus/data/beetle.expl.rda
===================================================================
(Binary files differ)


Property changes on: pkg/ordiconsensus/data/beetle.expl.rda
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/ordiconsensus/data/beetle.rda
===================================================================
(Binary files differ)


Property changes on: pkg/ordiconsensus/data/beetle.rda
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/ordiconsensus/man/RV.Rd
===================================================================
--- pkg/ordiconsensus/man/RV.Rd	                        (rev 0)
+++ pkg/ordiconsensus/man/RV.Rd	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,41 @@
+\name{RV}
+\alias{RV}
+
+\title{
+RV coefficient
+}
+\description{
+Calculates an RV coefficient.
+}
+\usage{
+RV(X, Y)
+}
+
+\arguments{
+  \item{X}{A matrix with the same number of rows as \code{Y}}
+  \item{Y}{A matrix with the same number of rows as \code{X}}
+}
+\details{
+ \code{RV} calculates the RV coefficient between two matrices.The RV coefficient is a multivariate generalization of the Pearson's correlation. The package \pkg{FactoMineR} offers a much more complete version of the RV coefficient (\code{\link[FactoMineR]{coeffRV}}). This function was designed to calculate only the RV coefficient making it much faster than \code{\link[FactoMineR]{coeffRV}}.
+}
+\value{
+The RV coefficient. A numeric value ranging between 0 and 1
+}
+\references{
+Escouffier, Y. (1973) \emph{Le traitement des variables vectorielles}. Biometrics \bold{29}:751--760.\cr
+}
+\author{
+F. Guillaume Blanchet
+}
+
+\seealso{
+\code{\link[FactoMineR]{coeffRV}}
+}
+\examples{
+data(beetle)
+data(beetle.expl)
+
+RV(beetle,beetle.expl)
+
+}
+\keyword{ multivariate }

Added: pkg/ordiconsensus/man/beetle.Rd
===================================================================
--- pkg/ordiconsensus/man/beetle.Rd	                        (rev 0)
+++ pkg/ordiconsensus/man/beetle.Rd	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,88 @@
+\name{beetle}
+\alias{beetle}
+\alias{beetle.expl}
+\docType{data}
+\title{
+Carabid beetle data with explanatory variables
+}
+\description{
+Carabid beetle data collected at 192 sites by Colin Bergeron in the summer of 2003 in Northwestern Alberta, Canada. See Bergeron et al. (2011, 2012) and Blanchet et al. (2013) for details.
+}
+\usage{
+data(beetle)
+data(beetle.expl)
+}
+\format{
+  The \code{beetle} data set is a data frame that contains the data on the following 37 carabid species:
+  \describe{
+    \item{\code{Agongrat}}{\emph{Agonum gratiosum}}
+    \item{\code{Agonplac}}{\emph{Agonum placidum}}
+    \item{\code{Agonretr}}{\emph{Agonum retractum}}
+    \item{\code{Agonsord}}{\emph{Agonum sordens}}
+    \item{\code{Agonsupe}}{\emph{Agonum superioris}}
+    \item{\code{Amarlitt}}{\emph{Amara littoralis}}
+    \item{\code{Amarluni}}{\emph{Amara lunicollis}}
+    \item{\code{Badiobtu}}{\emph{Badister obtusus}}
+    \item{\code{Bembgrap}}{\emph{Bembidion grapii}}
+    \item{\code{Bembrupi}}{\emph{Bembidion rupicola}}
+    \item{\code{Calaadve}}{\emph{Calathus advena}}
+    \item{\code{Calaingr}}{\emph{Calathus ingratus}}
+    \item{\code{Calofrig}}{\emph{Calosoma frigidum}}
+    \item{\code{Caracham}}{\emph{Carabus chamissonis}}
+    \item{\code{Dichcogn}}{\emph{Dicheirotrichus cognatus}}
+    \item{\code{Elapamer}}{\emph{Elaphrus americanus}}
+    \item{\code{Elaplapp}}{\emph{Elaphrus lapponicus}}
+    \item{\code{Harpfulv}}{\emph{Harpalus fulvilabris}}
+    \item{\code{Loripili}}{\emph{Loricera pilicornis}}
+    \item{\code{Miscarct}}{\emph{Miscodera arctica}}
+    \item{\code{Nebrgyll}}{\emph{Nebria gyllenhali}}
+    \item{\code{Notibore}}{\emph{Notiophilus borealis}}
+    \item{\code{Notidire}}{\emph{Notiophilus directus}}
+    \item{\code{Patrfove}}{\emph{Patrobus foveocollis}}
+    \item{\code{Patrsept}}{\emph{Patrobus septentrionis}}
+    \item{\code{Platdece}}{\emph{Platynus decentis}}
+    \item{\code{Platmann}}{\emph{Platynus mannerheimii}}
+    \item{\code{Pteradst}}{\emph{Pterostichus adstrictus}}
+    \item{\code{Pterbrev}}{\emph{Pterostichus brevicornis}}
+    \item{\code{Pterpens}}{\emph{Pterostichus pensylvanicus}}
+    \item{\code{Pterpunc}}{\emph{Pterostichus punctatissimus}}
+    \item{\code{Pterripa}}{\emph{Pterostichus riparius}}
+    \item{\code{Seriquad}}{\emph{Sericoda quadripunctata}}
+    \item{\code{Sterhaem}}{\emph{Stereocerus haematopus}}
+    \item{\code{Synuimpu}}{\emph{Synuchus impunctatus}}
+    \item{\code{Trecapic}}{\emph{Trechus apicalis}}
+    \item{\code{Trecchal}}{\emph{Trechus chalybeus}}
+
+The \code{beetle.expl} data set is a data frame that contains the relative basal area of the 25 trees closest to the center of each sampling site. The relative basal area is presented by tree species. 
+
+    \item{\code{Pt}}{Aspen (\emph{Populus tremuloides})}
+    \item{\code{Bp}}{White birch (\emph{Betula papyrifera})}
+    \item{\code{Ab}}{Balsam fir (\emph{Abie balsamea})}
+    \item{\code{Ll}}{Tamarack (\emph{Larix laricina})}
+    \item{\code{Pb}}{Balsam poplar (\emph{Populus balsamifera})}
+    \item{\code{Pc}}{Lodgepole pine (\emph{Pinus contorta})}
+    \item{\code{Pm}}{Black spruce (\emph{Picea mariana})}
+    \item{\code{Pg}}{White spruce (\emph{Picea glauca})}
+  }
+}
+\details{
+The \code{beetle} data has been previously transformed so that the ecological illustration presented in Blanchet et al. (in press) can be easily reproduced. The abundance of carabids was divided by the number of days each trap (there were three traps per sites) was active in the field. 
+}
+\source{
+All questions about these data (including if there are interest for using this data in publications) should be adressed to Colin Bergeron (cb1[at]ualberta.ca).
+}
+\references{
+Bergeron J.A.C., J.R. Spence, and W.J.A. Volney. 2011. Landscape patterns of species-level associations between ground-beetles (Coleoptera: Carabidae) and overstory trees in boreal forests of western Canada (Coleoptera: Carabidae). In Erwin, TL (Ed), \emph{Proceedings of a Symposium honoring the careers of Ross and Joyce Bell and their contributions to scientific work}, Burlington, VT, 12-15 June 2010. \emph{ZooKeys} \strong{147}: 577--600.
+
+Bergeron J.A.C., F.G. Blanchet, J.R. Spence, and W.J.A. Volney. 2012. Ecosystem classification and inventory maps as surrogates for ground beetle assemblages in boreal forest. \emph{Journal of Plant Ecology} \strong{5}:97--108.
+
+Blanchet, F.G., J.A.C. Bergeron, J.R. Spence, and F. He. 2013. Landscape effects of disturbance, habitat heterogeneity and spatial autocorrelation for a ground beetle (Carabidae) assemblage in mature boreal forest. \emph{Ecography} \strong{36}:636--647.
+
+Blanchet, F.G., P. Legendre, J.A.C. Bergeron, J.R. Spence, and F. He. \emph{in press}. Consensus RDA across dissimilarity coefficients for canonical ordination of community composition data. \emph{Ecological Monographs}.
+
+}
+\examples{
+data(beetle)
+data(beetle.expl)
+}
+\keyword{datasets}

Added: pkg/ordiconsensus/man/coeffCompare.Rd
===================================================================
--- pkg/ordiconsensus/man/coeffCompare.Rd	                        (rev 0)
+++ pkg/ordiconsensus/man/coeffCompare.Rd	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,99 @@
+\name{coeffCompare}
+\alias{coeffCompare}
+
+\title{ Compare dissimilarities used within the RDA framework}
+\description{
+  This function compares association coefficients used through the RDA framework with a minimum spanning tree. It was designed to compare how information explained by one dissimilarity coefficient diverge from the information explained by another. The comparison is made simultaneously on the site scores, the species scores and the canonical coefficients. 
+}
+\usage{
+coeffCompare(ordires, ordisigniaxis, pval = 0.05)
+}
+
+\arguments{
+  \item{ordires}{
+ A list of \code{\link[vegan]{rda}} or \code{\link[vegan]{capscale}} result object that includes a series of RDA or distance-based RDA (dbRDA) performed with different association coefficients on the same data.
+}
+  \item{ordisigniaxis}{
+ A list of \code{\link[vegan]{anova.cca}} object where each axis was tested for each RDA (or dbRDA) in \code{ordires}. This argument can also be a vector defining the number of significant axes in the RDA. See details. 
+}
+  \item{pval}{
+ Numeric. P-value threshold to select the number of axes to use. This argument is only active if a list of \code{\link[vegan]{anova.cca}} object is given for the argument \code{ordisigniaxis}, otherwise it is not considered. Default is 0.05.
+}
+}
+\details{
+For the argument \code{ordisigniaxis}, if a vector of number of significant axes is given, it is assumed that the significant axes are selected in sequential order from the first axis.
+The comparison made here rely on the RV coefficient \code{\link{RV}}, a multivariate generalization of the Pearson's correlation where matrix with the same number of rows are compared.
+\code{coeffCompare} should be used prior to using \code{\link{consensusRDA}} because it informs the user about the different association coefficients considered interesting to perform a consensus RDA. An association coefficient presenting results too different from the others should not be included in the consensus RDA, it should be considered apart or discarded.
+}
+\value{
+\item{RVmat}{A matrix of RV coefficients calculated from the sites scores matrices of RDA for all pairs of association coefficients}
+\item{mst}{minimum spanning tree calculated on (1-siteRVmat)}
+}
+\author{
+F. Guillaume Blanchet
+}
+
+\seealso{
+\code{\link{RV}}, \code{\link{consensusRDA}}
+}
+\examples{
+###################################################################
+### This example reproduces Figure 7b of Blanchet et al. (in press)
+###################################################################
+data(beetle)
+data(beetle.expl)
+
+### Construct results object
+ndis<-10
+ordiRes<-vector("list",length=ndis)
+
+#---------------------------------------------
+### Perform the various constrained ordination
+#---------------------------------------------
+### RDA species profile
+sp<-beetle/apply(beetle,1,sum)
+ordiRes[[1]]<-rda(sp~.,data=beetle.expl)
+	
+### RDA chord
+chord<-beetle/sqrt(apply(beetle^2,1,sum))
+ordiRes[[2]]<-rda(chord~.,data=beetle.expl)
+	
+### RDA Hellinger
+hell<-decostand(beetle,method="hellinger")
+ordiRes[[3]]<-rda(hell~.,data=beetle.expl)
+	
+### RDA chi2
+chisq<-decostand(beetle,method="chi.square")
+ordiRes[[4]]<-rda(chisq~.,data=beetle.expl)
+
+### db-RDA Bray-Curtis
+ordiRes[[5]]<-capscale(sqrt(vegdist(beetle,method="bray"))~.,data=beetle.expl,comm=beetle)
+	
+### db-RDA square-root Bray-Curtis
+ordiRes[[6]]<-capscale(sqrt(vegdist(beetle^0.5,method="bray"))~.,data=beetle.expl,comm=beetle)
+	
+### db-RDA fourth-root Bray-Curtis
+ordiRes[[7]]<-capscale(sqrt(vegdist(beetle^0.25,method="bray"))~.,data=beetle.expl,comm=beetle)
+	
+### db-RDA modified Gower log 2
+ordiRes[[8]]<-capscale(vegdist(decostand(beetle, "log",logbase=2), "altGower")~.,data=beetle.expl,comm=beetle) ### Warning message stem from log transformation of 0
+	
+### db-RDA modified Gower log 5
+ordiRes[[9]]<-capscale(vegdist(decostand(beetle, "log",logbase=5), "altGower")~.,data=beetle.expl,comm=beetle) ### Warning message stem from log transformation of 0
+	
+### db-RDA modified Gower log 10
+ordiRes[[10]]<-capscale(vegdist(decostand(beetle, "log",logbase=10), "altGower")~.,data=beetle.expl,comm=beetle) ### Warning message stem from log transformation of 0
+
+### Compare association coefficients
+AssoComp<-coeffCompare(ordiRes,rep(7,ndis))
+
+#---------------------------------------------
+### Draw a graphic to visualize the comparison
+#---------------------------------------------
+### Name of association coefficient compared
+name<-c("Species profiles","Chord","Hellinger","Chi2","Bray-Curtis","(Bray-Curtis)^0.5","(Bray-Curtis)^0.25","mGowerlog2","mGowerlog5","mGowerlog10")
+
+plot(AssoComp$mst,type="t",labels=name,xlab="",ylab="",main="MST Sites scores")
+
+}
+\keyword{ multivariate }

Added: pkg/ordiconsensus/man/consensusRDA.Rd
===================================================================
--- pkg/ordiconsensus/man/consensusRDA.Rd	                        (rev 0)
+++ pkg/ordiconsensus/man/consensusRDA.Rd	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,115 @@
+\name{consensusRDA}
+\alias{consensusRDA}
+\title{
+Calculates a consensus RDA
+}
+\description{
+Calculates a consensus RDA. That is a consensus of a series of RDA performed on the same data using different association coefficients.
+}
+\usage{
+consensusRDA(ordires, ordisigniaxis, resp.var, expl.var, pval = 0.05, scaling = 2)
+}
+%- maybe also 'usage' for other objects documented here.
+\arguments{
+  \item{ordires}{
+ A list of \code{\link[vegan]{rda}} or \code{\link[vegan]{capscale}} result object that includes a series of RDA or distance-based RDA (dbRDA) performed with different association coefficients on the same data.
+}
+  \item{ordisigniaxis}{
+ A list of \code{\link[vegan]{anova.cca}} object where each axis was tested for each RDA (or dbRDA) in \code{ordires}. This argument can also be a vector defining the number of significant axes in the RDA. See details. 
+}
+  \item{resp.var}{
+ A matrix of the response variables (species) used to construct all the RDA (or dbRDA) in \code{ordires}.
+}
+  \item{expl.var}{
+ A matrix of the explanatory variables used to construct all the RDA (or dbRDA) in \code{ordires}.
+}
+  \item{pval}{
+ Numeric. P-value threshold to select the number of axes to use. This argument is only active if a list of \code{\link[vegan]{anova.cca}} object is given for the argument \code{ordisigniaxis}, otherwise it is not considered. Default is 0.05.
+}
+  \item{scaling}{
+ Type of scaling used to project the results. A distance scaling (scaling 1) or a correlation scaling (scaling 2) can be used. Default is scaling 1. See details for more information.
+}
+}
+\details{
+For the argument \code{ordisigniaxis}, if a vector of number of significant axes is given for each RDA, it is assumed that the significant axes are selected in sequential order from the first axis.
+Although it is possible for the scaling to be 3 (it is available in the vegan package), this scaling should only be used for canonical correspondence analysis (CCA), it does not make any sense to use in the RDA framework.
+}
+\value{
+ \item{value}{A vector of eigenvalues associated to the axes of the consensus RDA}
+ \item{siteConsensus}{A matrix of consensus sites scores}
+ \item{spConsensus}{A matrix of consensus species scores}
+ \item{descConsensus}{A matrix of consensus canonical coefficient}
+}
+\author{
+F. Guillaume Blanchet
+}
+
+\examples{
+###########################################################################
+### This example reproduces Figure 7c of Blanchet et al. (in press)
+### 
+### However, for illustration purposes RDA axes were tested using anova.cca
+### This only has minor influence on the triplot. 
+###########################################################################
+data(beetle)
+data(beetle.expl)
+
+### Construct results object
+ndis<-9
+ordiRes<-vector("list",length=ndis)
+
+### RDA species profile
+sp<-beetle/apply(beetle,1,sum)
+ordiRes[[1]]<-rda(sp~.,data=beetle.expl)
+	
+### RDA chord
+chord<-beetle/sqrt(apply(beetle^2,1,sum))
+ordiRes[[2]]<-rda(chord~.,data=beetle.expl)
+	
+### RDA Hellinger
+hell<-decostand(beetle,method="hellinger")
+ordiRes[[3]]<-rda(hell~.,data=beetle.expl)
+
+### db-RDA Bray-Curtis
+ordiRes[[4]]<-capscale(sqrt(vegdist(beetle,method="bray"))~.,data=beetle.expl,comm=beetle)
+	
+### db-RDA square-root Bray-Curtis
+ordiRes[[5]]<-capscale(sqrt(vegdist(beetle^0.5,method="bray"))~.,data=beetle.expl,comm=beetle^0.5)
+	
+### db-RDA fourth-root Bray-Curtis
+ordiRes[[6]]<-capscale(sqrt(vegdist(beetle^0.25,method="bray"))~.,data=beetle.expl,comm=beetle^0.25)
+	
+### db-RDA modified Gower log 2
+ordiRes[[7]]<-capscale(vegdist(decostand(beetle, "log",logbase=2), "altGower")~.,data=beetle.expl,comm=decostand(beetle, "log",logbase=2)) ### Warning message stem from log transformation of 0
+	
+### db-RDA modified Gower log 5
+ordiRes[[8]]<-capscale(vegdist(decostand(beetle, "log",logbase=5), "altGower")~.,data=beetle.expl,comm=decostand(beetle, "log",logbase=5)) ### Warning message stem from log transformation of 0
+	
+### db-RDA modified Gower log 10
+ordiRes[[9]]<-capscale(vegdist(decostand(beetle, "log",logbase=10), "altGower")~.,data=beetle.expl,comm=decostand(beetle, "log",logbase=10)) ### Warning message stem from log transformation of 0
+
+#----------------
+### Test RDA axis 
+#----------------
+ordiResTest<-vector("list",length=ndis)
+
+for(i in 1:ndis){
+	ordiResTest[[i]]<-anova.cca(ordiRes[[i]],by="axis",cutoff=0.1)
+}
+
+### Consensus RDA
+consRDA<-consensusRDA(ordiRes,ordiResTest,beetle,beetle.expl)
+summary(consRDA)
+
+axisLabels<-c(paste("Axis 1 - ",round(consRDA$values[1]/sum(consRDA$values),4)*100,"%",sep=""),paste("Axis 2 - ",round(consRDA$values[2]/sum(consRDA$values),4)*100,"%",sep=""))
+plot(consRDA$siteConsensus[,1:2],pch=19,xlab=axisLabels[1],ylab=axisLabels[2],las=1)
+abline(h=0,v=0,lty=3)
+
+arrows(0,0,consRDA$spConsensus[,1]*10,consRDA$spConsensus[,2]*10,col="red",length=0.1,angle=12)
+text(consRDA$spConsensus[,1:2]*10,labels=rownames(consRDA$spConsensus),col="red")
+
+arrows(0,0,consRDA$descConsensus[,1]*0.4,consRDA$descConsensus[,2]*0.4,col="blue",length=0.1,angle=12)
+text(consRDA$descConsensus[,1:2]*0.4,labels=rownames(consRDA$descConsensus),col="blue")
+}
+
+\keyword{ multivariate }

Added: pkg/ordiconsensus/man/ordiconsensus-package.Rd
===================================================================
--- pkg/ordiconsensus/man/ordiconsensus-package.Rd	                        (rev 0)
+++ pkg/ordiconsensus/man/ordiconsensus-package.Rd	2013-06-17 13:44:09 UTC (rev 2529)
@@ -0,0 +1,27 @@
+\name{ordiconsensus-package}
+\alias{ordiconsensus-package}
+\alias{ordiconsensus}
+\docType{package}
+\title{
+Consensus RDA package
+}
+\description{
+Package to perform a consensus RDA, a method that finds the consensus among different RDA (or distance-based RDA) performed using different dissimilarity coefficients (e.g. Hellinger, Chord, Bray-Curtis)
+}
+\details{
+\tabular{ll}{
+Package: \tab ordiconsensus\cr
+Type: \tab Package\cr
+Version: \tab 0.3-1\cr
+Date: \tab 2012-11-12\cr
+License: \tab Unlimited\cr
+}
+This package includes different functions to perform a consensus RDA. It follows the work of Blanchet et al. (In press)
+}
+\author{
+F. Guillaume Blanchet
+
+Maintainer: F. Guillaume Blanchet <guillaume.blanchet at helsinki.fi>
+
+}
+\keyword{ multivariate }



More information about the Vegan-commits mailing list