主成分分析程序matlab,主成分分析的matlab实现完整程序.doc

主成分分析的matlab实现完整程序.doc

《计量地理学》(徐建华，高等教育出版社，2005)配套实习指导58主成分分析1.概述Matlab语言是当今国际上科学界(尤其是自动控制领域)最具影响力、也是最有活力的软件。它起源于矩阵运算，并已经发展成一种高度集成的计算机语言。它提供了强大的科学运算、灵活的程序设计流程、高质量的图形可视化与界面设计、与其他程序和语言的便捷接口的功能。Matlab语言在各国高校与研究单位起着重大的作用。主成分分析是把原来多个变量划为少数几个综合指标的一种统计分析方法，从数学角度来看，这是一种降维处理技术。1.1主成分分析计算步骤①计算相关系数矩阵(1)ppprrrrR212112在(3.5.3)式中，rij(i，j=1，2，…，p)为原变量的xi与xj之间的相关系数，其计算公式为(2)nknkjjiikjkjiiijxxr1122)()(因为R是实对称矩阵(即rij=rji)，所以只需计算上三角元素或下三角元素即可。《计量地理学》(徐建华，高等教育出版社，2005)配套实习指导59②计算特征值与特征向量首先解特征方程，通常用雅可比法(Jacobi)求出特征值0RI，并使其按大小顺序排列，即；然后分别求),21(pi0,21p出对应于特征值的特征向量。这里要求=1，即，i),2(pieie12pjie其中表示向量的第j个分量。ijeie③计算主成分贡献率及累计贡献率主成分的贡献率为iz),21(1pipki累计贡献率为),21(1pipkik一般取累计贡献率达85—95%的特征值所对应的第一、第二，m,21…，第m(m≤p)个主成分。④计算主成分载荷其计算公式为(3)),21,(),(pjiexzpljiij《计量地理学》(徐建华，高等教育出版社，2005)配套实习指导60得到各主成分的载荷以后，还可以按照(3.5.2)式进一步计算，得到各主成分的得分(4)nmnmzzzzZ2121122.程序结构及函数作用在软件Matlab中实现主成分分析可以采取两种方式实现：一是通过编程来实现；二是直接调用Matlab种自带程序实现。下面主要主要介绍利用Matlab的矩阵计算功能编程实现主成分分析。2.1程序结构主函数子函数2.2函数作用Cwstd.m——用总和标准化法标准化矩阵Cwfac.m——计算相关系数矩阵；计算特征值和特征向量；对主成分进行排序；计算各特征值贡献率；挑选主成分(累计贡献率大于85%)，输出主成分个数；计算主成分载荷Cwprint.mCwstd.mCwfac.mCwscore.m《计量地理学》(徐建华，高等教育出版社，2005)配套实习指导61Cwscore.m——计算各主成分得分、综合得分并排序Cwprint.m——读入数据文件；调用以上三个函数并输出结果3.源程序3.1cwstd.m%cwstd.m,用总和标准化法标准化矩阵functionstd=cwstd(vector)cwsum=sum(vector,1);%对列求和[a,b]=size(vector);%矩阵大小,a为行数,b为列数fori=1:aforj=1:bstd(i,j)=vector(i,j)/cwsum(j);endend3.2cwfac.m%cwfac.mfunctionresult=cwfac(vector);fprintf( 相关系数矩阵:\n )std=CORRCOEF(vector)%计算相关系数矩阵fprintf( 特征向量(vec)及特征值(val)：\n )[vec,val]=eig(std)%求特征值(val)及特征向量(vec)newval=diag(val);[y,i]=sort(newval);%对特征根进行排序，y为排序结果，i为索引fprintf( 特征根排序：\n )forz=1:length(y)newy(z)=y(length(y)+1-z);endfprintf( %g\n ,newy)rate=y/sum(y);fprintf( \n贡献率：\n )《计量地理学》(徐建华，高等教育出版社，2005)配套实习指导62newrate=newy/sum(newy)sumrate=0;newi=[];fork=length(y):-1:1sumrate=sumrate+rate(k);newi(length(y)+1-k)=i(k);ifsumrate>0.85break;endend%记下累积贡献率大85%的特征值的序号放入newi中fprintf( 主成分数：%g\n\n ,length(newi));fprintf( 主成分载荷：\n )forp=1:length(newi)forq=1:length(y)result(q,p)=sqrt(newval(newi(p)))*vec(q,newi(p));endend%计算载荷disp(result)3.3cwscore.m%cwscore.m,计算得分functionscore=cwscore(vector1,vector2);sco=vector1*vector2;csum=sum(sco,2);[newcsum,i]=sort(-1*csum);[newi,j]=sort(i);fprintf( 计算得分：\n )score=[sco,csum,j]%得分矩阵：sco为各主成分得分；csum为综合得分；j为排序结果3.4cwprint.m%cwprint.mfunctionprint=cwprint(filename,a,b);%filename为文本文件文件名，a为矩阵行数(样本数)，b为矩阵列数(变量指标数)fid=(filename, r )vector=fscanf(fid, %g ,[ab]);fprintf( 标准化结果如下：\n )v1=cwstd(vector)result=cwfac(v1);《计量地理学》(徐建华，高等教育出版社，2005)配套实习指导63cwscore(v1,result);4.程序测试4.1原始数据中国大陆35个大城市某年的10项社会经济统计指标数据见下表。城市名称年底总人口(万人)非农业人口比(%)农业总产值(万元)工业总产值(万元)客运总量(万人)货运总量(万吨)地方财政预算内收入(万元)城乡居民年底储蓄余额(万元)在岗职工人数(万人)在岗职工工资总额(万元)北京1249.900.59781843427199997062032345562279086326806646410.805773301天津910.170.5809150113622645502325926317112807311301931202.682254343石家庄875.400.23322918680688576829291911352348709587595.60758877太原299.920.65632360382737750193711895203277394310088.65654023呼和浩特207.780.44123653438164522351