基于R的群落学多元统计分析

1基于的群落多元统计分析——用vegan包进行排序分析赖江山（janson)[e-mail]的缩写,是群落分析的package作者：JariOksanen://cc.oulu.fi/~jarioksa/softhelp/vegan.htmllibrary(vegan)什么是排序(ordination)?排序的过程是将样方或植物种排列在一定的空间，使得排序轴能够反映一定的生态梯度，从而，能够解释植被或植物种的分布与环境因子间的关系，也就是说排序是为了揭示植被-环境间的生态关系。因此，排序也叫梯度分析（gradientanalysis）。间接梯度分析(Indirectgradientanalysis)直接梯度分析(directgradientanalysis)2个种的排序图3个种的排序图4个种的排序图？？？40个种排序图？？？排序的目标：1.降低维数，减少坐标轴的数目；2.由降低维数引起的信息损失尽量少，即发生最小的畸变，也就是让新的坐标系第1-3轴排序轴包含大量的生态信息。排序的目的：表示植被与环境之间的关系：所有排序方法都反映植物种和环境之间的关系以及在某一环境梯度上的种间关系。1.线形模型（linearmodel），短的梯度，主成分分析（Principlecomponentanalysis)，需要对数据进行非线性转换，如取对数；2.非线性模型（non-linearmodel）如高斯模型，长的梯度，对应分析（Correspondenceanalysis)群落数据输入gtsdata=read.table(gtsdata.txt,header=T)gtsdatadim(gtsdata)环境因子数据输入gtsenv=read.table(gtsenv.txt,header=T);gtsenvdim(gtsenv)数据的标准化1.decostand(x,method,MARGIN,…)total:除以行和或列和(defaultMARGIN=1是row);max:除以行或列的最大值(defaultMARGIN=2是列);freq:除以行或列的最大值,并乘以非零值的个数，非零值的平均值为1(defaultMARGIN=2);normalize:使行或列的平方和等于1(defaultMARGIN=1);range:标准化使行或列的值在0...1(defaultMARGIN=2).standardize:标准化使行或列的和为1且方差为1(defaultMARGIN=2)；pa:将数据转换为0、1数据；chi.square:除以行和及列和的平方根；hellinger:采用total标准化以后再取平方根；log：对数化，默认自然对数，logbase参数是自选的base2.wisconsin(x):除以列最大值，再除以行和。排序类别(inCANOCO)间接梯度分析（IndirectGradientAnalysis）：PCA(Principalcomponentsanalysis)CA(Correspondenceanalysis)DCA(DetrendedCorrespondenceAnalysis)直接梯度分析（DirectGradientAnalysis）：RDA(Redundanceanalysis)CCA(Canonicalcorrespondenceanalysis)DCCA(DetrendedCCA)PCARDACACCADCADCCA13决定排序的模型：单峰还是线性？decorana(gtsdata)Call:decorana(veg=gtsdata)Detrendedcorrespondenceanalysiswith26segments.Rescalingofaxeswith4iterations.DCA1DCA2DCA3DCA4Eigenvalues0.39390.22390.095550.06226Decoranavalues0.50250.17560.067120.03877Axislengths3.25952.51301.214451.00854如果这四个轴中梯度最长（最大值）超过4，选择单峰模型排序（CA、CCA、DCA）更合适。如果是小于3，选择线性模型（PCA、RDA）比较合理。如果介于3-4之间，单峰模型和线性模型结果差不多。间接梯度分析（IndirectGradientAnalysis）PCA(Principalcomponentsanalysis)CA(Correspondenceanalysis)DCA(DetrendedCorrespondenceAnalysis)主成分分析(Principlecomponentanalysis,PCA)主成分分析的主要原理是：使坐标旋转一定的角度后，使第一轴表示数据最大的方差，使第二轴表示数据第二的方差。而且轴与轴之间是正交的（orthogonal)。-4-2024-4-2024PCA和RDA都采用函数rda实现：在vegan包中，rda(formula,data,scale=FALSE,...)rda(X,Y,Z,scale=FALSE,...)scores(x,choices,display=c(sites,species),...)在stat包中：princomp(x,...)主成分分析princomp(formula,data=NULL,subset,na.action,...)gts.rda=rda(gtsdata)gts.rdaCall:rda(X=gtsdata)InertiaRankTotal352.1Unconstrained352.122InertiaisvarianceEigenvaluesforunconstrainedaxes:PC1PC2PC3PC4PC5PC6PC7PC8111.77973.58054.60732.95926.48118.06312.7637.637scores(gts.rda,choices=c(1:4),display=c(si,sp))summary(gts.rda)#类似Canoco的log文件和.sol文件的信息plot(gts.rda,choices=c(1,2),display=c(sp,si))biplot(gts.rda,choices=c(1,2),display=c(sp,si))plot(rda(gtsdata,scale=T))plot(rda(gtsdata))!!如果不对数据做标准化的话，丰富种的值就非常大，排序时就只能看清丰富种的位置，其它种就拥挤在一起。•如用x1,x2,x3,x4,x5,x6…分别表示原先的变量，而用y1,y2,y3,y4,y5,y6…表示新的主成分，那么，第一和第二主成分为11234562123456-0.806-0.674-0.6750.8930.8250.8360.3530.5310.5130.3060.4350.425yxxxxxxyxxxxxx•这些系数称为主成分载荷（loading），它表示主成分和相应的原先变量的相关系数。比如y1表示式中x1的系数为-0.806，这就是说第一主成分和的x1变量的相关系数为-0.806。相关系数(绝对值）越大，主成分对该变量的代表性也越大。负荷(loading)gts.pca=princomp(gtsdata)gts.pca$loadingsgtsenv.pca=princomp(gtsdenv)gtsenv.pca$loadingsbiplot(gtsenv.pca)第一主成分代表海拔高度，第二主成分代表坡向对应分析(Correspondenceanalysis,CA)1.PCA在迭代运算过程是采用线性模型2.CA在迭代运算过程采用单峰模型（加权平均法）CA在vegan中也是用cca函数来实现：gts.ca=cca(gtsdata)gts.casummary(gts.ca)Call:cca(X=gtsdata)InertiaRankTotal1.424Unconstrained1.42421InertiaismeansquaredcontingencycoefficientEigenvaluesforunconstrainedaxes:CA1CA2CA3CA4CA5CA6CA7CA80.502530.265640.140230.105020.091270.055400.050630.0420426plot(gts.ca)CANOCO里面scalingofordinationscores27plot(gts.ca,scaling=1)用物种数据对样方坐标进行加权平均,使样方坐标在物种数据的中心，因此对样方感兴趣的话，采用这种做图方法。plot(gts.ca,scaling=2)plot(gts.ca,scaling=3)用样方数据对物种坐标进行加权平均,使物种数据在样方数据的中心，因此对物种感兴趣的话，采用这种做图方法。如果一个物种靠近某个样方，表明该物种可能对该样方的位置起很大的作用。特别是对于二元数据的排序，这个样方可能就代表该物种。如图中，20号样方与短柄枹(QUESER)靠得比较近，表明：短柄枹表征了20号样方的特征，19号样方与20号样方距离近，生态关系也较近。•只在少数样方出现的物种通常在排序空间的边缘，表明它们只偶然发生，或它们只在稀有生境（如米槠CASCAR）。•在排序空间中心的物种，可能在取样区域是该物种最优分布区，如甜槠，或有两个或多个最优分布区，或与前两个轴不相关。除趋势对应分析(Detrendedcorrespondenceanalysis,DCA)：CA采用单峰曲线表示物种和环境关系CA产生的弓形效应CA的第二排序轴在许多情况下是第一轴的二次变形，即所谓的“弓形效应”（Archeffect）或者“马蹄形效应”（horse-shoeeffect）。（详见张金屯《群落生态学》168页）×和•分别代表除趋势前和除趋势后的样方排序（引自Hill和Gauch1980）DCA在R中的实现采用函数decorana。decorana(veg,iweigh=0,iresc=4,ira=0,mk=26,short=0,before=NULL,after=NULL)veg:群落数据;iweig:稀有物种的权重；(稀有物种影响比较大）iresc：纠正弓形效应的次数；ira:分析的类型(DCA:0,CA:1);mk：校正弓形效应轴的分段数;short:需要校正的最短梯度。plot(decorana(gtsdata))plot(cca(gtsdata))直接梯度分析（DirectGradientAnalysis）：RDA(Redundanceanalysis)CCA(Canonicalcorrespondenceanalysis)pRDA(partialRDA)pCCA(partialCCA)冗余分析(redundancyanalysis,RDA)及典范对应分析（Canonicalcorrespondenceanalysis,CCA)1.通常采用PCA处理环境数据，采用CA处理群落数据，但这些方法都只能处理一个数据表；2.RDA和CCA是多元分析（PCA，CA）和线性回归的结合，研究植被和环境之间的关系。38当我们在解释变量（环境因子数据）与响应变量（物种数据）之间建立预测模型的时候，经常会遇到这样的情况，往往我们仅仅考察解释变量中某几个环境因子的对物种数据的影响，但剩下的环境因子也会对物种产生影响，这些剩余环境因子我们经常称为协变量(Covariables)。在CANOCO中，协变量的