备忘录之主成分分析

备忘录之主成分分析

//这个暑期小学期，跟管理技术系的同学上营销数据分析。这门选修，主要是为凑学分了，偶尔也在课程bbs发些贴，权当备忘录了。这个主成分分析，主要来自高惠璇老师的两本教材，《应用多元统计分析》（北京大学出版社，2005）和《实用统计方法与SAS系统》（北京大学出版社，2001）。这个贴就发在学院bbs上，其他ip登不上，贴过来。在记事本上写的，像Z1这样的标记就是Z加一个下标1。

在实际问题中，如果变量个数太多：
1.变量彼此存在一定的相关性，因而观测到的数据在一定程度上反映的信息有所重叠（多重共线性）；
2.在高维空间中研究样本的分布规律比较复杂，变量太多增加分析问题的复杂性。

我们希望：
1.用较少的综合变量来代替原来较多的变量（降维，减少维度）；
2.这几个综合变量能够尽可能多地反映原来变量的信息；
3.这几个综合变量彼此互不相关。

主成分分析、因子分析、典型相关分析和偏最小二乘回归等统计方法就是以上降维思想的体现。

一、主成分分析

1.定义

这里说的主成分分析，又称主分量分析、主轴分析，是将多个指标（变量）化为少数几个综合指标的一种统计方法。
把p个变量X1，X2，…，Xp，记为一个p维的随机向量X=(X1，X2，…，Xp)，其协方差阵为D(X)。考虑X的线性变换：

Z1=A1*X
Z2=A2*X
……
Zp=Ap*X，

这里的X和A1、A2、…、Ap等都不妨看成向量形式。假如我们想用Z1来代替原来的p个变量，这就要求Z1尽可能多地反映原来p个变量的信息。这里“信息”可以用Z1的方差Var（Z1）来表示，方差Var（Z1）越大，表示Z1包含的信息越多。当然，这需要强加一些数学上的限制，否则Var（Z1）就可能是无限大了，这里的限制是向量A1和它自己的转置之积等于1，记为A1*Trans（A1）=1。就这样：

         若存在满足A1*Trans（A1）=1的A1，使得Var（Z1）最大，则称Z1为为第一主成分，或第一主分量，Z1=A1*X。

如果第一主成分不足以代表原来p个变量的绝大部分信息，我们就可以考虑X的第二个线性组合Z2=A2*X。此时，我们要求，已经体现在第一主成分Z1中的信息不要出现在Z2中，即Z1和Z2的协方差Cov（Z1，Z2）=0。就这样：

          在Cov（Z1，Z2）=0时，若存在满足A2*Trans（A2）=1的A2，使得Var（Z2）最大，则称Z2为为第二主成分，或第二主分量。

类似我们可以定义X的第三主成分，以致第p主成分（当然，对p维的随机向量X来说，第p主成分就没有必要了）。

2.直观解释

从代数上讲，主成分就是p个原始变量的一些特殊的线性组合。
从几何上讲，这些线性组合是把由X1，X2，…，Xp构成的坐标系通过旋转而产生的新坐标系。

3.求解

从上面的定义可知，求第一主成分Z1=A1*X的问题，就是求A1，使得在A1*Trans（A1）=1下，Var（Z1）达到最大，这是一个条件极值问题，可以用拉格朗日乘子法求解。

   Max Var(Z1)=Var(A1*X)
   s.t A1*Trans（A1）=1

运用拉格朗日乘子法,令L=Var（A1*X）-t(A1*Trans（A1）-1)，接下来就是一阶条件了，按下不表。
说说结果。求解以上条件极值问题，就是求X的协方差阵D(X)的特征值（就是拉格朗日乘子t）和特征向量（就是我们感兴趣的A1、A2、…、Ap）。一个定理：

设X=(X1，X2，…，Xp)是p维的随机向量，其协方差阵为D(X)，D(X)的特征值为t1>t2>…>tp,A1、A2、…、Ap为相应的单位正交特征向量，则X的第i主成分为Zi=Ai*X.

4.一些名词

原总体X的总惯量（总方差）：就是上面我们求出来的特征值之和，sum（t1，t2，…，tp）
贡献率：单个特征值ti与总方差之比，称为主成分Zi的贡献率。同样可以定义累计的贡献率。
因子负荷量（因子载荷量）：主成分Zi与原始变量Xi的相关系数。

5.计算机求解（SAS）

（在SAS系统中，主成分分析的程序步是princomp，因子分析的程序步是factor。所以可以知道主成分分析和因子分析毕竟不是一样的东西，具体见因子分析的部分。）

这里有一个数据，5个变量，id、x1身高、x2体重、x3胸围、x4坐高，30个观测值。我们想对x1身高、x2体重、x3胸围、x4坐高做主成分分析。

data princomp;
  input id x1 x2 x3 x4  ;
  cards;
1 148 41 72 78  
2 139 34 71 76
3 160 49 77 86
4 149 36 67 79
5 159 45 80 86
6 142 31 66 76
7 153 43 76 83  
8 150 43 77 79
9 151 42 77 80
10 139 31 68 74
11 140 29 64 74
12 161 47 78 84
13 158 49 78 83
14 140 33 67 77
15 137 31 66 73
16 152 35 73 79
17 149 47 82 79
18 145 35 70 77
19 160 47 74 87
20 156 44 78 85
21 151 42 73 82
22 147 38 73 78
23 157 39 68 80
24 147 30 65 75
25 157 48 80 88
26 151 36 74 80
27 144 36 68 76
28 141 30 67 76
29 139 32 68 73
30 148 38 70 78
;
proc princomp;
  var x1 x2 x3 x4  ;
run;

结果解释：

（1）结果首先有相关矩阵，可以看出这四个变量之间有较强的相关性。
                                          Correlation Matrix
                                      x1          x2          x3          x4
                          x1      1.0000      0.8632      0.7321      0.9205
                          x2      0.8632      1.0000      0.8965      0.8827
                          x3      0.7321      0.8965      1.0000      0.7829
                          x4      0.9205      0.8827      0.7829      1.0000

（2）然后是以上相关阵的特征值t（Eigenvalue，上面例子中的t）
                                         Eigenvalues of the Correlation Matrix
                                   Eigenvalue    Difference    Proportion    Cumulative
                        1    3.54109800    3.22771484        0.8853        0.8853
                        2    0.31338316    0.23397420        0.0783        0.9636
                        3    0.07940895    0.01329906        0.0199        0.9835
                        4    0.06610989                             0.0165        1.0000

从中可以看到，第一主成分的贡献率就达到 0.8853。单个贡献率就是某个特征值与总特征值的比，这里0.8853=3.54109800/(3.54109800 +0.31338316 +0.07940895+0.06610989) .还有，前两个主成分的累计是0.9636，所以用这两个主成分就可以很好地概括这些数据。

（3）还有特征向量（ Eigenvectors，我们例子中的A）
                                                    Eigenvectors
                                      Prin1         Prin2         Prin3         Prin4
                      x1      0.496966      -.543213      -.449627      0.505747
                      x2      0.514571      0.210246      -.462330      -.690844
                      x3      0.480901      0.724621      0.175177      0.461488
                      x4      0.506928      -.368294      0.743908      -.232343

我们可以由最大的两个特征值对应的特征向量，写出第一和第二主成分（结果中的Prin1和Prin2）：
Z1=Prin1= 0.496966*X1+0.514571 *X2+0.480901*X3+0.506928*X4
Z2=Prin2= -0.543213 *X1+0.210246 *X2+0.724621*X3 -0.368294*X4

这样我们就写出了开头那个线性组合。接下来说些因子分析，因子分析是主成分分析的推广，稍复杂些。（待续）

转载：http://johnthu.spaces.live.com/blog/cns!2053CD511E6D5B1E!270.entry