聚类分析算法解析

1、聚类分析算法解析一、不相似矩阵计算1.加载数据data(iris)str(iris)> data(iris)> str(iris)1 data,frame1:150oba.of5vatiatoles:$Sepal.Length:num5,.55.44.6S4.44.9SSepal.Width:num3<533*23.13.63*93.43.12,93*1一$Petal.Length:nuio1,.$Petal.Width.:numCL20i.0.40»

2、.一$Species:Factorw/3levelssetosa-j-?verslcolor*,.21111111111-分类分析是无指导的分类，所以删除数据中的原分类变量。iris$Species<-NULL2 .不相似矩阵计算不相似矩阵计算，也就是距离矩阵计算，在R中采用dist()函数，或者cluster包中的daisy()函数。dist()函数的基本形式是dist(x,method="euclidean",diag=FALSE,upper=FALSE,p=2)其中x是数据框(数据集)，而方法可以指定为欧式距离"euclidean”,最大距离"

3、;maximum",绝对值距离"manhattan","canberra",二进制距离非对称"binary"和明氏距离"minkowski"。默认是计算欧式距离，所有的属性必须是相同的类型。比如都是连续类型，或者都是二值类型。dd<-dist(iris)str(dd)>str(dd)Class'disfatomic1:111750,5390.51口.64号0.110.616.一一attr(*z,rSizeF,J=int150.attr=logiFALSE.atvr"Upp已l

4、ogiFALSE.一attr(*,frmethcdr,)=chrfFeuclideanF,*一attr(*t*rcall=Languagedist(x=iris)距离矩阵可以使用as.matrix()函数转化了矩阵的形式，方便显示。Iris数据共150例样本间距离矩阵为150行列的方阵。下面显示了15号样本间的欧式距离。dd<-as.matrix(dd)>str(d可num1:150,1:150100,5390.510,6480.141.-attr("diimnaiiies"=Listof2.$schr1:150Hl,frr2H"4".：ch

5、r1:150Ttr2nft>ddl:5,1:51io.oonoooo2O.5365L6530,509902040+64日口74150.141421420.S38S15S0.00000000,30000000.331662504603276330.5099020,3000000.0000000.2449490.50990240.64807410,33166250,244949口0.00000000.648074150.14142140.60027630.50990200.6430741040000000二、用hclust()进行谱系聚类法(层次聚类)1 .聚类函数R中自带的聚类函数是hcl

6、ust(),为谱系聚类法。基本的函数指令是结果对象<-hclust(距离对象,method=方法)hclust()可以使用的类间距离计算方法包含离差法"ward",最短距离法"single",最大距离法"complete",平均距离法"average","mcquitty",中位数法"median"和重心法"centroid"。下面采用平均距离法聚类。hc<-hclust(dist(iris),method="ave")&g

7、t;hcCall:hclust(d-disc(itis)rmethod=rrave")Clustermethod:averageDistance:euclideanNumberofobjects：ISO2 .聚类函数的结果聚类结果对象包含很多聚类分析的结果，可以使用数据分量的方法列出相应的计算结果。str(hc)>str(he)Listof7$merge:int1:149,1:2J-102-9-1-10-129-11-5-20-30-58$height：mini1:1490,1,10.10.1.$order:intL:ISO1421516333437213224.$labels

8、:NULLSmethod：chr"aver日gie”5call:Languagehclust(ddist(iris),fmethod",fave*'$disc.method:chrrFeuclidean,r-attrpF，cla33r,)=chr"hclust"下面列出了聚类结果对象hc包含的merge和height结果值的前6个。其行编号表示聚类过程的步骤，X1,X2表示在该步合并的两类，该编号为负代表原始的样本序号，编号为正代表新合成的类；变量height表示合并时两类类间距离。比如第1步，合并的是样本102和143,其样本间距离是0.0,合

9、并后的类则使用该步的步数编号代表，即样本-102和-143合并为1类。再如第6行表示样本11和49合并，该两个样本的类间距离是0.1,合并后的类称为6类。head(hc$merge,hc$height)>head)XIX2he.height1-102-1430.020-400T13-1-180.14rlCl-350.1-129-1330*1-11490.1卜面显示的聚类结果对象hc包含的merge和height结果值的5055步的结果。第50步结果表明样本43与13类(即第13步的聚类合类结果)合并。所有的类编号负数j表示原数据的样本编号，正数i表示聚类过程的第i步形成的新类。再如54步

10、，表示聚类过程第7步和第37步形成的类合并为新类，新类成为第54类，其类间距离是0.2641715。data.frame(hc$merge,hc$height)50:55>data,ftame(hege,he?height)50i55,XIX2he.height50-4313D.25S1139519310.258680252-71180,261803453-6B160.263B95B547370.264171555-69一日日0,26457513 .绘制聚类图聚类完成后可以使用plot()绘制出聚类的树图。plot(hc,hang=-1,labels=iris$Species)Clust

11、erDendrogramdistfiris)hclustC,"average"4 .指定分类和类中心聚类通过观察树形图，可见由于数据例很多，使图形很乱，不容易确定合理的分类，为简化图形，使用cutree()来确定最初的分类结果，先初步确定各个样本的最初分类的类数后，然后用hclust()再次聚类,重新聚类后，得出最后的分类结果。memb<-cutree(hc,k=10)#确定10个分类table(memb)#各类中的样本数>>menii<-cutree(hcfk=ICQ>界各娄中的样本数>table(me描)menib1234567891

12、045413722424913为了进行重新聚类，需要计算各类的类中心的类中心，这里用各类的均向量来代表各类的类中心，计算程序如下。cent<-NULLfor(kin1:10)函数，从10类起重新开始聚类。hclust()函cent<-rbind(cent,colMeans(irisSamplememb=k,drop=FALSE)4-cent<-rbind(cent,colMeans(irisineirto=k,drop+>centSepaliLengthSepal*UidthPetal+LengthPetal.Width1,144968889343866671.4711

13、110.24666672J5,5500004,1750001,4000000.22500003z4.50000D2.3OD0OO1.3000000.30000006.2351352.851351447405411.5945946S.6318182.69545S4,0590911.24090915.000a0D2.3000003.2750001.0250000Q6.5416673.050000545291672.12OB3337.4333332,9222226.2666671.98888899z4.9000002.5000004.5000001.700000010J7.6000003.73333

14、3644000002.2333333FALSE)有了各类的类中心后,cent<-MULLlorI氏in1：!)再次使用hclust()数这时需要指定各类的类中心，并采用中心法来聚类。hc1<-hclust(dist(cent),method="centroid",members=table(memb)plot(hc1)hc1ClusterDendrogram4IO)迂工ir>oLTJCMdistfcent)hclust"centroid")再次聚类的树形图就很精简了。从树形图上看，确定为三类是比较合适的。4,输出最终分类结果使用cutr

15、ee()来确定输出各个样本的最后分类。memb<-cutree(hc,k=3#确定3个分类table(memb)#各类中的样本数> mert±i<-cutree(he,k:=3)> #各类中的样本数> merttocutree(he,k=3)#各类中的样本数> tab工已(merrib)uneirib123506436下面将聚类结果和原始的分类比较可见聚类结果是不错的。错误率约为14/150=9.33%。> data(itis)> table(iris$Speciesfmeinb)menib312setosaSO00versicoLor

16、0500virefinica014365,直接数据距离阵聚类可以直接输入距离矩阵，然后用as.dist()函数转为hclust()可以使用的距离阵对象,然后再用hclust()再进行聚类分析。x<-read.table(text="idBAFIMIVORMTOBA0662877255412996FI6620295468268400MI8772950754564138VO2554687540219869RM4122685642190669TO9964001388696690",header=T)s(x)<-x$idx$id<-NULLxBAF

17、IHIVORMTOBA0662877255412996FI£62295468268400MI377295075456413BVO25546875402198S9RM4122685642190669TO996400138369669转换为距离阵。x<-as.dist(x)>x<-a3.ciist(xJBAFIMIVORMFI662MIS77295VO255468754RM412268564219TO996400138369669hc<-hclust(x)显示聚类过程。data.frame(hc$merge,hc$height)> data.frameEme

18、rgeThc2Sheight)XIX2hc2*height1-3-613S2-4-52193-214004-12412534996绘制聚类图。plot(hc2)ClusterDendrogramLLwXhclust(*."complete")从聚类图可见聚为两类比较合适。cutree(hc2,2)>cutree(he2f2)BAFIMIVORMTO122112所以BAVQRM聚为1类而TQFI、MI聚为另一类。四、使用cluster包进行聚类分析1.计算距离矩阵除了可以使用dist()函数外，加载了cluster包后还可以使用daisy()函数，该函数可以使用的数据类

19、型更多一些。当原始数据表中包含的数据是混合型即属性变量既有有连续又有分类的时候，dist()计算就不够完美，而采用daisy()函数就可以毫无问题。使用daisy()必须加载cluster包。第一次使用还必须安装。install.packages('cluster',repos='htt/25/rpkg)library(cluster)daisy()函数的形式和dist()相似，可以采用的距离包含欧式距离"euclidean",绝对值距离和"manhattan","gower"。dd2&l

20、t;-daisy(iris)dd2<-as.matrix(dd1)> dd2<-daisy(iris)> str(dd2)Classes1dissimilarity',1dist1atomicr1:111750*5390.510,6480.1410.61,-attr(*,FSizer,)-int1.50T.-attrr，HetricM)-chrauclidean>dd2<-asTnieitri：w(dd2J>dd21:5,1:510.ooooooo20,538516530,509902040.64807415CL14142142O.538516

21、S,OOOOOOO0.30000000.33156250.60927630.50A902a.3000000.000000.2449490.50990240.6480741.33166250.24494900.OOOOOOO.64配74150,1414214口.60827630.5099020口.6480741.oaoooo2.使用agens()聚类cluster包提供了两种层次聚类的方法，一种是abottom-up自底而上，采用先将数据样本先各自看成一类，然后通过合并法的聚类过程，最终形成一个大类，包含全部样本，agnes()函数完成这个功能，这个称为合并层次聚类。而diana()则采用相反的

22、方法，即atop-down自顶而下，先将所有样本看成一类，然后通过分裂类的过程将样品最终分为各自样本类，这个称为分裂法。agens()函数的形式是agnes(x,diss=TRUE|FALSE,metric="euclidean",stand=FALSE|TRUE,method="average",par.method,keep.diss=n<100,keep.data=!diss)其中的diss指定x对象是否是不相似矩阵对象，与hclust()函数不同，agnes()既可以使用原始数据也可以采用距离阵。stand是要求对数据进行标准化，而可以采用

23、的类间距离计算方法有“average”平均法，"UPGMA,“single”最短距离法，"complete”最大距离法，"ward"离差法，"weighted"权重平均距离法和"flexible"弹性法。keep.diss和keep.data指定在结果中是否保留不相似矩阵和数据，保留这些结果需要更多的内存。ag<-agnes(iris,diss=F,metric="euclidean",stand=F,method="single")>ag<-agnes(d

24、iss=Fmetric=ppeuclideanrr.fstand=Fmethod=TFsinle'rp)>agCall:agnes(xiris,dissF,metrioeuclidean.*1*,standmethodppsingleA5Slomerativecoefficient:0-8493364Order1ofobjects;38B40SO2B293611492427442101IB41S193546133031263448939431272022471437J253721321761933344516152342515387525557SS5966767598787772

25、S4907D81826BB39373899596971009162566785649279748086607191J12813912412714715010214311412.27384134104117138105129ri0913311114B112142146113110121141145125116137149103126127J13065101120108131115636988106123119136135109110561459461991D7ll：z：132Height(sunimaxy):Min.1stQu.MedianMean3rdQu.Max.Oi.OOOD0.1732口

26、.26460.29210,3464la400Availablecomponents:1border-*rtheightw*ac1*pWrgepf?piMthodp?>data,freutie(ag$rrierge,agheight)50：55fXIX250-71175141-745242-44534965453525S5448ag.height1.64012190.26457510,3828427.31622780426457510.3162278plot(ag,ask=FALSE,which.plots=NULL)Bannerofagnes(x=iris,diss=Fsm00.20.4

27、0608HeightAgglomerativeCoefficient=0.85sineDendrogramofagnes(x=iris,diss=F,metric="euclidean",stand=F,method=irisAgglomerativeCoefficient二。85cutree(ag,k=3)>cutreek=3)1111111111111111111111111111111111111111111114S111111222222222222222222222222222222222222228922222222Z222Z222Z2

28、22ZZ222Z2223222Z2222Z222Z3133J222.222222222222222load'adultuci.rdata')str(AdultUCI)>str(AdultUCI)1data*frame1:$age$trotfeelassobs.39ofSOeducation:int7751663311215646237213334092S453216018720964inari-tai.status:Ord.factortj/ISlevels,Preschcol,f<t,13t;-4th,f<.:14149occupation:int13139

29、71314591413relationship:FactorW7levelsied-AF-3pouseF.：53race:Factorw/14levelsl，rAdw-clei:ical,f.:1466104341SEX：Factorw/6levelsUxiskjand",rrHct-in-famiLy,|r,.+:2121capItai.gain:Factortr/5levelsAnien-IrLdiarL-EsJciinc.:5553353capital.loss:Factorw/ZlevelsF，Female,r,MalE,r:2222111212,hours.per,week:int21740000000140G45178.一native.country:int0000000000.income;int4013404040401645500.1JA:Faecer时/41levelsrrCanribodiaFf/f,Can.ada,f,.:3g3939395nA：Ord.factorWZlevels"smallF，<，rlarge":111111122444644:Factorw/8levelsF，Federal-go