一直都抽不出时间写第十七篇读书笔记,直到前两天去面PDD被问到主成分分析的推导,完全想不起来的同时也决定回来把PCA这个坑给填上,顺便加深一下记忆。主成分分析(principal component analysis, PCA)利用正交变换把由线性相关变量表示的观测数据转换为少数几个由线性无关变量表示的数据,这些线性无关的变量被称为主成分。
一、总体主成分分析
主成分分析中,首先对给定数据进行规范化,使得数据每一变量的平均值为0,方差为1,再对数据进行正交变换,把原来线性相关变量表示的数据,变成若干个线性无关的新变量表示的数据。这里认为新变量的方差表示新变量上信息的大小,因此为了尽量保存更多的信息,新变量选择所有可能的正交变换中变量的方差最大的若干个。
设x=(x1,x2,⋯ ,xm)Tx=(x_1,x_2,\cdots,x_m)^Tx=(x1,x2,⋯,xm)T是mmm维随机变量,其均值向量是μ\muμ,即
μ=E(x)=(μ1,μ2,⋯ ,μm)\mu=E(x)=(\mu_1,\mu_2,\cdots,\mu_m)μ=E(x)=(μ1,μ2,⋯,μm) 注意,上面的xix_ixi是随机变量,可以想象成xix_ixi是个列向量,它在不同的样本中有不同的取值,而μi\mu_iμi是一个值,可以想象成xix_ixi在所有样本中的平均值。接着计算xxx的协方差矩阵Σ\SigmaΣ
Σ=cov(x,x)=E[(x−E(x))(x−E(x))T]=E[(x−μ)(x−μ)T]\Sigma=cov(x,x)=E[(x-E(x))(x-E(x))^T]=E[(x-\mu)(x-\mu)^T]Σ=cov(x,x)=E[(x−E(x))(x−E(x))T]=E[(x−μ)(x−μ)T] 接着我们来考虑对mmm维随机变量xxx进行多个线性变换,第iii个yiy_iyi线性变换就可以表示为
yi=αiTx=α1ix1+α2ix2+⋯+αmixmy_i=\alpha_i^Tx=\alpha_{1i}x_1+\alpha_{2i}x_2+\cdots+\alpha_{mi}x_myi=αiTx=α1ix1+α2ix2+⋯+αmixm 如果我们有mmm个线性变换,就可以得到由xxx变换而来的新的mmm维随机变量y=(y1,y2,⋯ ,ym)Ty=(y_1,y_2,\cdots,y_m)^Ty=(y1,y2,⋯,ym)T。我们对新的mmm维随机变量进行研究,注意,这里是对单个yiy_iyi随机变量进行探讨,所以下面的E(yi),var(yi),cov(yi,yj)E(y_i),var(y_i),cov(y_i,y_j)E(yi),var(yi),cov(yi,yj)所指代的都是一个值,而非矩阵或者向量。
首先来看新随机变量yiy_iyi的均值E(yi)E(y_i)E(yi)
E(yi)=E(αiTx)=αiTE(x)=αiTμE(y_i)=E(\alpha_i^Tx)=\alpha_i^TE(x)=\alpha_i^T\muE(yi)=E(αiTx)=αiTE(x)=αiTμ 再来看yiy_iyi的方差var(yi)var(y_i)var(yi),由于
var(yi)=E[(yi−E(yi))(yi−E(yi))T]var(y_i)=E[(y_i-E(y_i))(y_i-E(y_i))^T]var(yi)=E[(yi−E(yi))(yi−E(yi))T] 代入yi,E(yi)y_i,E(y_i)yi,E(yi)的表达式可得
var(yi)=E[(αiTx−αiTμ)(αiTx−αiTμ)T]=αiTE[(x−μ)(x−μ)T]αivar(y_i)=E[(\alpha_i^Tx-\alpha_i^T\mu)(\alpha_i^Tx-\alpha_i^T\mu)^T]=\alpha_i^TE[(x-\mu)(x-\mu)^T]\alpha_ivar(yi)=E[(αiTx−αiTμ)(αiTx−αiTμ)T]=αiTE[(x−μ)(x−μ)T]αi 即
var(yi)=αiTΣαivar(y_i)=\alpha_i^T\Sigma\alpha_ivar(yi)=αiTΣαi 最后来看两个新随机变量yi,yjy_i,y_jyi,yj的协方差cov(yi,yj)cov(y_i,y_j)cov(yi,yj)
cov(yi,yj)=E[(αiTx−αiTμ)(αjTx−αjTμ)T]=αiTΣαjcov(y_i,y_j)=E[(\alpha_i^Tx-\alpha_i^T\mu)(\alpha_j^Tx-\alpha_j^T\mu)^T]=\alpha_i^T\Sigma\alpha_jcov(yi,yj)=E[(αiTx−αiTμ)(αjTx−αjTμ)T]=αiTΣαj 推导出这一系列公式后,下面给出总体主成分需要满足的三个条件。
(1)αi\alpha_iαi是单位向量,即αiTαi=1\alpha_i^T\alpha_i=1αiTαi=1。
(2)变量yi,yjy_i,y_jyi,yj互不相关,即cov(yi,yj)=0cov(y_i,y_j)=0cov(yi,yj)=0。
(3)变量y1y_1y1是所有线性变换中方差最大的,y2y_2y2是所有与y1y_1y1不相关的线性变换中方差最大的,依次类推。
因此,设Σ\SigmaΣ是mmm维随机变量xxx的协方差矩阵,且Σ\SigmaΣ的特征值分别是λ1≥λ2≥⋯≥0\lambda_1\ge\lambda_2\ge\cdots\ge0λ1≥λ2≥⋯≥0,特征值对应的特征向量分别是α1,α2,⋯ ,αm\alpha_1,\alpha_2,\cdots,\alpha_mα1,α2,⋯,αm,则xxx的第kkk主成分是
yk=αkTxy_k=\alpha_k^Txyk=αkTx xxx的第kkk主成分的方差是
var(yk)=αkTΣαk=λk, k=1,2,⋯ ,mvar(y_k)=\alpha_k^T\Sigma\alpha_k=\lambda_k, k=1,2,\cdots,mvar(yk)=αkTΣαk=λk, k=1,2,⋯,m 因此有
Σαk=λkαk\Sigma \alpha_k=\lambda_k\alpha_kΣαk=λkαk 用矩阵表示即为
ΣA=AΛ\Sigma A=A\LambdaΣA=AΛ 这里的A=[aij]m×mA=[a_{ij}]_{m\times m}A=[aij]m×m,Λ\LambdaΛ是对角矩阵,由上式可以得到两个公式
ATΣA=ΛA^T\Sigma A=\LambdaATΣA=ΛΣ=AΛAT\Sigma=A\Lambda A^TΣ=AΛAT 下面叙述总体主成分的性质
(1)总体主成分yyy的协方差矩阵是对角矩阵
cov(y)=cov(ATx)=ATcov(x)A=ATΣA=Λcov(y)=cov(A^Tx)=A^Tcov(x)A=A^T\Sigma A=\Lambdacov(y)=cov(ATx)=ATcov(x)A=ATΣA=Λ (2)总体主成分yyy的方差之和等于随机变量xxx的方差之和
∑i=1mvar(xi)=∑i=1mσii=tr(Σ)=tr(AΛAT)=tr(Λ)=∑i=1mλi=∑i=1mvar(yi)\sum_{i=1}^mvar(x_i)=\sum_{i=1}^m\sigma_{ii}=tr(\Sigma)=tr(A\Lambda A^T)=tr(\Lambda)=\sum_{i=1}^m\lambda_i=\sum_{i=1}^mvar(y_i)i=1∑mvar(xi)=i=1∑mσii=tr(Σ)=tr(AΛAT)=tr(Λ)=i=1∑mλi=i=1∑mvar(yi) 上面σii\sigma_{ii}σii是随机变量xix_ixi的方差,即协方差矩阵Σ\SigmaΣ的对角线元素。
(3)第kkk个主成分yky_kyk与变量xix_ixi的相关系数为
ρ(yk,xi)=λkαikσii\rho(y_k,x_i)=\frac{\sqrt{\lambda_k}\alpha_{ik}}{\sqrt{\sigma_{ii}}}ρ(yk,xi)=σiiλkαik (4)第kkk个主成分yky_kyk与mmm个变量的因子负荷量满足
∑i=1mσiiρ2(yk,xi)=λk\sum_{i=1}^m\sigma_{ii}\rho^2(y_k,x_i)=\lambda_ki=1∑mσiiρ2(yk,xi)=λk (5)mmm个主成分与第iii个变量xix_ixi的因子负荷量满足
∑k=1mρ2(yk,xi)=1\sum_{k=1}^m\rho^2(y_k,x_i)=1k=1∑mρ2(yk,xi)=1 主成分分析的主要目的是降维,所以一般选择k(k≤m)k(k\le m)k(k≤m)个主成分来使问题得以简化。具体选择kkk的方法,通常利用方差贡献率。
第kkk主成分yky_kyk的方差贡献率定义为yky_kyk的方差与所有方差之和的比,记作ηk\eta_kηk
ηk=λk∑i=1mλi\eta_k=\frac{\lambda_k}{\sum_{i=1}^m\lambda_i}ηk=∑i=1mλiλk 通常取kkk使得累计方差贡献率达到规定的百分比以上。
接下来要讨论的是规范化变量的总体主成分分析,设xxx是mmm维随机变量,令
xi∗=xi−E(xi)var(xi), i=1,2,⋯ ,mx_i^*=\frac{x_i-E(x_i)}{\sqrt{var(x_i)}}, i=1,2,\cdots,mxi∗=var(xi)xi−E(xi), i=1,2,⋯,m 显然,规范化随机变量x∗x^*x∗的协方差矩阵就是相关矩阵RRR,因此规范化随机变量的总体主成分y∗y^*y∗有以下性质:
(1)y∗y^*y∗的协方差矩阵是
Λ∗=diag(λ1∗,λ2∗,⋯ ,λm∗)\Lambda^*=diag(\lambda_1^*,\lambda_2^*,\cdots,\lambda_m^*)Λ∗=diag(λ1∗,λ2∗,⋯,λm∗) 其中λ1∗≥λ2∗≥⋯λm∗≥0\lambda_1^*\ge\lambda_2^*\ge\cdots\lambda_m^*\ge0λ1∗≥λ2∗≥⋯λm∗≥0是相关矩阵RRR的特征值。
(2)y∗y^*y∗的协方差矩阵的特征值之和为mmm
∑k=1mλk∗=m\sum_{k=1}^m\lambda_k^*=mk=1∑mλk∗=m
二、样本主成分分析
上面一节讲述的是总体主成分分析,针对的是xxx这个mmm维随机变量,但在实际问题中没有随机变量,只有一个变量在nnn个样本上的观测值,这需要用到样本主成分分析了。
设观测数据用样本矩阵XXX表示,记作
X=[x1,x2,⋯ ,xn]X=[x_1,x_2,\cdots,x_n]X=[x1,x2,⋯,xn] 给定样本矩阵后,可以估计样本均值向量xˉ\bar xxˉ
xˉ=1n∑j=1nxj\bar x=\frac{1}{n}\sum_{j=1}^nx_jxˉ=n1j=1∑nxj 也可以计算任意两个特征向量i,ji,ji,j的协方差与相关系数
sij=1n−1∑k=1n(xik−xˉi)(xjk−xˉj)s_{ij}=\frac{1}{n-1}\sum_{k=1}^n(x_{ik}-\bar x_i)(x_{jk}-\bar x_j)sij=n−11k=1∑n(xik−xˉi)(xjk−xˉj)rij=sijsiisjjr_{ij}=\frac{s_{ij}}{\sqrt{s_{ii}s_{jj}}}rij=siisjjsij 进而计算出样本协方差矩阵SSS与样本相关矩阵RRR
S=[sij]m×mS=[s_{ij}]_{m\times m}S=[sij]m×mR=[rij]m×mR=[r_{ij}]_{m\times m}R=[rij]m×m 这里的SSS与RRR都是对总体的无偏估计。
在使用样本主成分时,一般假设样本数据都是规范化的,即
xij∗=xij−xˉisiix_{ij}^*=\frac{x_{ij}-\bar x_i}{\sqrt{s_{ii}}}xij∗=siixij−xˉi 这时样本协方差矩阵SSS就是样本相关矩阵RRR
R=1n−1XXTR=\frac{1}{n-1}XX^TR=n−11XXT
接着叙述主成分分析的求解方法,传统的求解方法是通过数据的相关矩阵特征值分解进行的,具体步骤为:
(1)对数据矩阵XXX进行规范化处理,仍以XXX表示
(2)计算样本相关矩阵
R=1n−1XXTR=\frac{1}{n-1}XX^TR=n−11XXT (3)求RRR的kkk个特征值和对应的单位特征向量
除了传统的特征值分解求法,现在常用的是通过数据矩阵的奇异值分解进行,首先定义一个新的n×mn\times mn×m矩阵X′X'X′
X′=1n−1XTX'=\frac{1}{\sqrt{n-1}}X^TX′=n−11XT 由于
X′TX′=1n−1XXT=RX'^TX'=\frac{1}{n-1}XX^T=RX′TX′=n−11XXT=R 因此对X′X'X′求截断奇异值分解X′=UΣVTX'=U\Sigma V^TX′=UΣVT后,VVV的列向量就是RRR的单位特征向量。对X′X'X′进行截断奇异值分解,得到
X′=UΣVTX'=U\Sigma V^TX′=UΣVT 矩阵VVV的前kkk列构成kkk个样本主成分,最终求得k×nk\times nk×n样本主成分矩阵
Y=VTXY=V^TXY=VTX