多元统计分析-(4)

R语言与数据分析培训第4讲多元统计分析多元数据3NjPixxxxxxxxxxxxxXPNpppNNij,2,1;21}{32122322211131211，，，多元数据基本方法4聚类（cluster)排序(ordination)5Doubs鱼类数据集6该数据集包括3个矩阵，第一个矩阵是27种鱼类在每个样方的多度，第二个矩阵包括11个与河流的水文、地形和水体化学属性相关的环境变量，第三个矩阵是样方的地理坐标（笛卡尔坐标系，X和Y）。变量名称代码单位离源头距离daskm海拔Altma.s.l.坡度pen‰平均最小流量debm3s-1pH值pH-钙浓度（硬度）durmgL-1磷酸盐浓度phomgL-1硝酸盐浓度nitmgL-1铵浓度ammmgL-1氧含量oxymgL-1生物需氧量dbomgL-1法国和瑞士边境的Jura山脉的Doubs河7聚类分析8层次法（hierarchical）和非层次法（non-hierarchical）9距离（distance）（相异性）矩阵。例如欧氏距离、Jaccard相似系数排序的对象-距离矩阵212121(a,a)=()niiidistyy基于连接的层次聚类10平均聚合聚类11#加载所需的程序包library(ade4)library(vegan)#应该先加载ade4后加载vegan以避免冲突library(gclus)library(cluster)library(RColorBrewer)library(labdsv)library(mvpart)library(MVPARTwrap)#MVPARTwrap这个程序包必须从本地zip文件安装#导入CSV格式的数据spe-read.csv(DoubsSpe.csv,row.names=1)env-read.csv(DoubsEnv.csv,row.names=1)spa-read.csv(DoubsSpa.csv,row.names=1)#删除无物种数据的样方8spe-spe[-8,]env-env[-8,]spa-spa[-8,]12#物种多度数据：先计算样方之间的弦距离矩阵，然后进行单连#接聚合聚类spe.norm-decostand(spe,normalize)spe.ch-vegdist(spe.norm,euc)spe.ch.single-hclust(spe.ch,method=single)par(mfrow=c(2,2))#使用默认参数选项绘制聚类树plot(spe.ch.single)#计算完全连接聚合聚类#********************spe.ch.complete-hclust(spe.ch,method=complete)plot(spe.ch.complete)13#计算UPGMA聚合聚类#***********************spe.ch.UPGMA-hclust(spe.ch,method=average)plot(spe.ch.UPGMA)#这个UPGMA聚合聚类树看起来介于单连接聚类和完全连接聚类之间。这种#情况经常发生。#计算鱼类数据的形心聚类#***********************spe.ch.centroid-hclust(spe.ch,method=centroid)plot(spe.ch.centroid)14Ward最小方差聚类1516#计算Ward最小方差聚类#***********************par(mfrow=c(2,1))spe.ch.ward-hclust(spe.ch,method=ward)plot(spe.ch.ward)#使用距离平方造成此聚类树上半部分过于膨胀。为了使聚类树比例看起来#更协调而不影响结构，可以使用当前融合水平的平方根重新绘图spe.ch.ward$height-sqrt(spe.ch.ward$height)plot(spe.ch.ward)非层次聚类•非层次聚类（non-hierarchicalclustering）是对一组对象进行简单分组的方法，也可以表述为：在p维空间内有n个对象（点），将n个对象分为k组（或称为聚类簇），分组的依据是尽量使组内的对象之间比组间对象之间的相似度更高。此时用户需要自己决定分组的数量k。非层次聚类的算法首先需要有个初始的结构，即首先将所有对象任意分为k组，然后在初始结构的基础上进行不断替换迭代，以达到最优化的分组结果。初始结构的设定可以依据某种理论，但大多数情况下是随机分配。通常是设定不同的初始结构，然后通过大量的迭代以找到最佳的解决方案。1718#预转化后物种数据k-均值划分#****************************spe.kmeans-kmeans(spe.norm,centers=4,nstart=100)spe.kmeans#注意：即使给定的nstart相同，每次运行上述命令，所产生的结果也不一定#完全相同，因为每次运算设定的初始结构是随机的。#如果不知道多多少组，怎么办，下面SSI图帮助你19#k-均值划分，2组到10组#************************spe.KM.cascade-cascadeKM(spe.norm,inf.gr=2,sup.gr=10,iter=100,criterion=ssi)plot(spe.KM.cascade,sortg=TRUE)#该图显示每个对象在每种分类组数下的归属（图上每行代表一种组数）。图#内的表格有不同的颜色，每行两种颜色，代表分两组k=2，三种颜色代表k=3，#依此类推。右图代表不同k值条件下的中止标准的统计量。此系列中，到底#多少组数是最佳方案？如果倾向于较大的组数，哪个是最佳方案呢？2021物种环境因子多元回归树（MRT）空间时间MRT22海拔=340铵浓度0.11铵浓度=0.45海拔340铵浓度=0.11铵浓度0.453.08:n=130.95:n=30.356:n=31.27:n=10误差:0.37CV误差:0.7标准差:0.10123cpx变量相对误差0.20.40.60.81.01.2Inf0.210.0910.0650.0480.0440.0310.0210.0140.008201234567891011树的大小Min+1SE24#多元回归树#**********library(mvpart)#载入mvpart程序包spe.ch.mvpart-mvpart(data.matrix(spe.norm)~.,env,margin=0.08,cp=0,xv=pick,xval=nrow(spe),xvmult=100,which=4)#：参数xv=pickpick允许通过人机交互的方式选择所需的分类组数。如果想自动选择具#有最小CVRE值的回归树，可以设定xv=min。#此处可以点击想要的分组数（例如4组）25#MRT作为空间和时间系列约束聚类方法spe.ch.seq-mvpart(as.matrix(spe)~das,env,cp=0,xv=pick,margin=0.08,xval=nrow(spe),xvmult=100,which=4)#此时可以点击所期望的分组组数summary(spe.ch.seq)#组的组成（终端节点的标识）(gr-spe.ch.seq$where)#重新编排聚类簇的编号aa-1gr2-rep(1,length(gr))for(iin2:length(gr)){if(gr[i]!=gr[i-1])aa-aa+1gr2[i]-aa}26#在Doubs河地图上标注样方所属的聚类簇plot(spa,asp=1,type=n,main=MRTgroups,xlab=xcoordinate(km),ylab=ycoordinate(km))lines(spa,col=lightblue)text(50,10,Upstream,cex=1.2,col=red)text(25,115,Downstream,cex=1.2,col=red)k-length(levels(factor(gr2)))for(iin1:k){points(spa[gr2==i,1],spa[gr2==i,2],pch=i+20,cex=3,col=i+1,bg=i+1)}text(spa,row.names(spa),cex=0.8,col=white,font=2)legend(bottomright,paste(Group,1:k),pch=(1:k)+20,col=2:(k+1),pt.bg=2:(k+1),pt.cex=2,bty=n)27282930排序（Ordination)•非约束排序（UnconstrainedOrdination）•典范排序（CanonicalOrdination）31PetrTerBraak333435=37131&highlight=2个种的排序图3个种的排序图4个种的排序图？？？40个种排序图？？？Whyordination?37Whyordination?排序的目标：1.降低维数，减少坐标轴的数目；2.由降低维数引起的信息损失尽量少，即发生最小的畸变，也就是让新的坐标系前几轴排序轴包含大量的生态信息。排序类别•非约束排序UnconstrainedOrdinationPCA(Principalcomponentsanalysis)CA(Correspondenceanalysis)(DCA)PCoA(PrincipalCoordinateAnalysis)NMDS(NonmetricMultidimensionalScaling)非度量多维尺度分析•约束排序(典范排序)constrainedOrdinationRDA(Redundanceanalysis)CCA(Canonicalcorrespondenceanalysis)LDA(lineardiscriminantanalysis)(CCorA),canonicalcorrelationanalysis;(CoIA)co-inertiaanalysis,(MFA)multiplefactoranalysis（对称）非对称特征根方法非约束排序（UnconstrainedOrdination）•主成分分析（principalcomponentanalysis，PCA）：基于特征向量的主要排序方法。分析对象是原始的定量数据。排序图展示样方之间的欧氏距离。•对应分析（correspondenceanalysis，CA）：分析对象必须是频度或类频度、同量纲的非负数据。排序图展示行（对象）或列（变量）之间的卡方距离。在生态学研究中主要用于分析物种数据。•主坐标分析（principalcoordinateanalysis，PCoA）：分析对象为距离矩阵（大部分为Q模式），而非原始的样方-变量矩阵表格。因此，可以灵活选择关联测度。•非度量多维尺度分析（nonmetricmultidimensionalscaling，NMDS）：与前面三种排序方法不同，NMDS不是基于特征向量提取的排序方法。NMDS尝试在预先设定数量的排序轴去排序对象，目标是保持这些对象排位关系（orderingrelationships）不变。40PCA41协方差42#导入本章所需的程序包library(ade4)library(vegan)library(gclus)library(ape)#导入CSV文件数据spe-read.csv(DoubsSpe.csv,row.names=1)env-read.csv(DoubsEnv.csv,row.names=1)spa-read.csv(DoubsSpa.csv,row.names=1)#删除