主成分分析和MATLAB应用.doc

主成分分析
类型：一种处理高维数据的方法。
降维思想：在实际问题的研究中，往往会涉及众多有关的变量。但是，变量太多不但会增加计算的复杂性，而且也会给合理地分析问题和解释问题带来困难。一般说来，虽然每个变量都提供了一定的信息，但其重要性有所不同，而在很多情况下，变量间有一定的相关性，从而使得这些变量所提供的信息在一定程度上有所重叠。因而人们希望对这些变量加以“改造”，用为数极少的互补相关的新变量来反映原变量所提供的绝大部分信息，通过对新变量的分析达到解决问题的目的。
一、总体主成分
1.1 定义
设 X1，X2，…，Xp 为某实际问题所涉及的 p 个随机变量。记
X=(X1，X2，…,Xp)T，其协方差矩阵为

它是一个 p 阶非负定矩阵。设
（1）
则有
（2）
第 i 个主成分：
一般地，在约束条件

及

下，求 li 使 Var(Yi)达到最大，由此 li 所确定的

称为 X1，X2，…，Xp 的第 i 个主成分。
1.2 总体主成分的计算
设 是的协方差矩阵，的特征值及相应的正交单位化特征向量分别为

及

则 X 的第 i 个主成分为
（3）
此时

记 为主成分向量，则 Y=PTX，其中，且

定义第 i个主成分的贡献率：；
前m个主成分累计贡献率：，它表明前 m 个主成分Y1，Y2，…，Ym综合提供 X1，X2，…，Xp中信息的能力。
1.4 标准化变量的主成分
在实际问题中，不同的变量往往有不同的量纲，由于不同的量纲会引起各变量取值的分散程度差异较大，这时总体方差则主要受方差较大的变量的控制。为了消除由于量纲的不同可能带来的影响，常采用变量标准化的方法，即令
（5）
其中 这时

的协方差矩阵便是

的相关矩阵 ，其中
（6）
利用 X 的相关矩阵 作主成分分析，有如下结论：
设 为标准化的随机向量，其协方差矩阵（即 X 的相关矩阵）为 ，则 的第 i 个主成分为
（7）
并且
（8）
其中 为 的特征值，为相应于特征值 的正交单位特征向量。
第 i 个主成分的贡献率：；
前 m 个主成分的累计贡献率：；
与的相关系数为 。

二、样本主成分
前面讨论的是总体主成分，但在实际问题中，一般 （或）是未知的，需要通过样本来估计。设

为取自

的一个容量为n的简单随机样本，则样本协方差矩阵及样本相关矩阵分别为
（9）
其中

分别以 S 和 R 作为和的估计，然后按总体主成分分析的方法作样本主成分分析。

三、 例 题
某市为了全面分析机械类个企业的经济效益，选择了8个不同的利润指标，14企业关于这8个指标的统计数据如下表所示，试进行主成分分析。
表1 14家企业的利润指标的统计数据
变量
企业序号
净产值利润率（%）
xi1
固定资产利润率（%）
xi2
总产值利润率（%）
xi2
销售收入利润率（%）
xi3
产品成本利润率(%)
xi5
物耗利润率(%)
xi6
人均利润率
xi7(千元/人)
流动资金利润率(%)
xi8
1
40.4
24.7
7.2
6.1
8.3
8.7
2.442
20.0
2
25.0
12.7
11.2
11.0
12.9
20.2
3.542
9.1
3
13.2
3.3
3.9
4.3
4.4
5.5
0.578
3.6
4
22.3
6.7
5.6
3.7
6.0
7.4
0.176
7.3
5
34.3
11.8
7.1
7.1
8.0
8.9
1.726
27.5
6
35.6
12.5
16.4
16.7
22.8
29.3
3.017
26.6
7
22.0
7.8
9.9
10.2
12.6
17.6
0.847
10.6
8
48.4
13.4
10.9
9.9
10.9
13.9
1.772
17.8
9
40.6
19.1
19.8
19.0
29.7
39.6
2.449
35.8
10
24.8
8.0
9.8
8.9
11.9
16.2
0.789
13.7
11
12.5
9.7
4.2
4.2
4.6
6.5
0.874
3.9
12
1.8
0.6
0.7
0.7
0.8
1.1
0.056
1.0
13
32.3
13.9
9.4
8.3
9.8
13.3
2.126
17.1
14
38.5
9.1
11.3

9.5
12.2
16.4
1.327
11.6

解：样本均值向量为：
,
样本协方差矩阵为：

由于S中主对角线元素差异较大，因此我们样本相关矩阵R出发进行主成分分析。样本相关矩阵R为：

矩阵R的特征值及相应的特征向量分别为：
特征值
特征向量
6.1366
0.32113 0.29516 0.38912 0.38472 0.37955 0.37087 0.31996 0.35546
1.0421
-0.4151 -0.59766 0.22974 0.27869 0.31632 0.37151 -0.27814 -0.15684
0.43595
-0.45123 0.10303 -0.039895 0.053874 -0.037292 0.075186 0.77059 -0.42478
0.22037
-0.66817 0.36336 -0.22596 -0.11081 0.14874 0.069353 -0.13495 0.55949
0.15191
-0.038217 0.62435 0.12273 -0.036909 0.15928 0.21062 -0.43006 -0.58105
0.0088274
-0.10167 0.13584 -0.15811 0.86226 -0.25204 -0.34506 -0.13934 -0.026557
0.0029624
0.1596 -0.061134 -0.53966 0.046606 0.7609 -0.27809 0.06203 -0.13126
0.0012238
0.19295 -0.031987 -0.64176 0.11002 -0.25397 0.68791 -0.006045 -0.0054031

R的特征值及贡献率见下表
特征值
贡献率（%）
累计贡献率（%）
6.1366
0.76708
0.76708
1.0421
0.13027
0.89734
0.43595
0.054494
0.95184
0.22037
0.027547
0.97938
0.15191
0.018988
0.99837
0.0088274
0.0011034
0.99948
0.0029624
0.0003703
0.99985
0.0012238
0.00015297
1
前3个标准化样本主成分类及贡献率已达到95.184%，故只需取前三个主成分即可。
前3个标准化样本主成分中各标准化变量 前的系数即为对应特征向量，由此得到3个标准化样本主成分为
注意到，y1近似是8个标准化变量的等权重之和，是反映各企业总效应大小的综合指标，y1的值越大，则企业的效益越好。由于y1的贡献率高达76.708%，故若用y1的得分值对各企业进行排序，能从整体上反映企业之间的效应差别。将S中sii的值及中各的值以及各
企业关于xi的观测值代入y1的表达式中，可求得各企业y1的得分及其按其得分由大到小的排序结果。
企业序号
得分
12
-0.97354
4
-0.64856
3
-0.62743
11
-0.48558
10
-0.21949
7
-0.189
14
-0.004803
5
0.016879
8
0.17711
13

0.18925
1
0.29351
2
0.65315
6
0.85566
9
0.96285

所以，第9家企业的效益最好，
第12家企业的效益最差。
Matlab程序：[coeff,score,latent]=princomp(X)
注：该函数使用协方差阵作主成分分析。
主成分分析程序
a=[];
b=corrcoef(zscore(a))%计算相关系数矩阵
D=tril(b)%得到三角矩阵
[d,v]=eig(b)%计算特征值和特征向量
y1=zscore(a)*d(:,7)%计算第一主成分数值
[f1,i1]=sort(y1);
[f2,i2]=sort(i1);
[flipud(i1),flipud(f1),f2]%第一主成分得分排序
y2=zscore(a)*d(:,6)%计算第二主成分数值
[f1,i1]=sort(y2);
[f2,i2]=sort(i1);
[flipud(i1),flipud(f1),f2]%第二主成分得分排序
y3=zscore(a)*d(:,5)%计算第三主成分数值
[f1,i1]=sort(y3);
[f2,i2]=sort(i1);
[flipud(i1),flipud(f1),f2]%第三主成分得分排序
y4=zscore(a)*d(:,4)计算第三主成分数值
[f1,i1]=sort(y4);
[f2,i2]=sort(i1);
[flipud(i1),flipud(f1),f2]%第si主成分得分排序
y5=zscore(a)*d(:,3)%计算第一主成分数值
[f1,i1]=sort(y1);
[f2,i2]=sort(i1);
[flipud(i1),flipud(f1),f2]%第一主成分得分排序
y=y5*(0.64/(0.64+0.84+1.04+1.17+2.36))+y4*(0.84/(0.64+0.84+1.04+1.17+2.36))+y3*(1.04/(0.64+0.84+1.04+1.17+2.36))+y2*(1.17/(0.64+0.84+1.04+1.17+2.36))+y1*(2.36/(0.64+0.84+1.04+1.17+2.36))
[f1,i1]=sort(y);
[f2,i2]=sort(i1);
[flipud(i1),flipud(f1),f2]%综合主成分得分排序