主成分分析|PCA算法大全

文章目录

主成分分析|PCA算法大全
- 1. PCA原理
- - 1.1 最大方差理论
  - 1.2 最小平方误差理论
  - 1.3 高维数据下的特征值分解
- 2. CCIPCA增量主元分析算法[1]
- 3. KPCA
- 4. 2DPCA[2]
- 5. 2D2DPCA
- 6. BDPCA
- 参考文献

主元分析提供了一种用低维数据来表示高维复杂数据最主要特征的途径。PCA的思想是将nnn维特征映射到kkk维空间上k<nk<nk<n，这kkk维特征是全新的正交特征，是新构造出来的kkk维特征，而不是简单地从nnn维特征中去除其余n−kn-kn−k维特征。在原始PCA的基础上，根据数据的特点，研究人员提出各种版本的PCA，但是核心思想一样。例如，针对非线性数据的核PCA，针对序列数据的CCIPCA，针对二维图片数据的2DPCA、2D2DPCA以及BDPCA等。下面首先详细介绍PCA的思想，然后简单介绍其他版本的PCA算法。

1. PCA原理

（本部分参考从零开始实现主成分分析(PCA)算法）
PCA可以被定义为数据在低维线性空间上的正交投影，这个线性空间被称为主子空间（principal subspace），使得投影数据的方差被最大化（Hotelling, 1933），即最大方差理论。等价地，它也可以被定义为使得平均投影代价最小的线性投影，即最小误差理论。平均投影代价是指数据点和它们的投影之间的平均平方距离（Pearson, 1901）。

1.1 最大方差理论

在信号处理中，认为信号具有较大的方差，噪声有较小的方差，信噪比就是信号与噪声的方差比，越大越好。因些我们认为，最好的kkk维特征是将nnn维样本点变换为kkk后，每一维上的样本方差都尽可能的大。

首先，考虑在一维空间(k=1)(k=1)(k=1)上的投影。我们可以使用nnn维向量uuu定义这个空间的方向。为了方便并不失一般性，我们假定选择一个单位向量，从而uTu=1u^T u=1uTu=1。

(假设数据是零均值化后的)
如上图所示，红色点表示原样本点x(i)x^{(i)}x(i)，uuu是蓝色直线的斜率也是直线的方向向量，而且是单位向量，直线上蓝色点表示原样本点x(i)x^{(i)}x(i)在uuu上的投影。容易知道投影点离原点的距离是x(i)Tux^{(i)T} ux(i)Tu,由于这些原始样本点的每一维特征均值都为0，因此投影到uuu上的样本点的均值仍然是0。

假设原始数据集为Xm×nX_{m\times n}Xm×n，我们的目标是找到最佳的投影空间Wn×k=(w1,w2,...,wk)W_{n\times k}=(w_1, w_2, ..., w_k)Wn×k=(w1,w2,...,wk)，其中wiw_iwi是单位向量且wiw_iwi是单位向量且wiw_iwi与wj(i≠j)w_j(i\neq j)wj(i=j)正交，何为最佳的WWW？就是原始样本点投影到WWW上之后，使得投影后的样本点方差最大。

由于投影后均值为0，因些投影后的总方差为：1m∑i=1m(x(i)Tw)2=1m∑i=1mwTx(i)x(i)Tw=∑i=1mwT(1mx(i)x(i)T)w\frac{1}{m}\sum_{i=1}^m (x^{(i)T}w)^2=\frac{1}{m}\sum_{i=1}^m w^Tx^{(i)} x^{(i)T}w=\sum_{i=1}^m w^T(\frac{1}{m}x^{(i)}x^{(i)T})wm1i=1∑m(x(i)Tw)2=m1i=1∑mwTx(i)x(i)Tw=i=1∑mwT(m1x(i)x(i)T)w

其中1mx(i)x(i)T\frac{1}{m}x^{(i)}x^{(i)T}m1x(i)x(i)T就是原始数据集XXX的协方差矩阵（因为x(i)x^{(i)}x(i)的均值为0，因为无偏估计的原因，一般协方差矩阵除以m−1m-1m−1，这是用m）。

记λ=1m∑i=1m(x(i)Tw)2,Σ=1mx(i)x(i)T\lambda=\frac{1}{m}\sum_{i=1}^m (x^{(i)T}w)^2,\Sigma=\frac{1}{m}x^{(i)}x^{(i)T}λ=m1∑i=1m(x(i)Tw)2,Σ=m1x(i)x(i)T,则有λ=wTΣw\lambda=w^T\Sigma wλ=wTΣw

上式两边同时左乘www，注意到wwT=1ww^T=1wwT=1（单位向量），则有λw=Σw\lambda w=\Sigma wλw=Σw
所以www是矩阵Σ\SigmaΣ的特征值所对应的特征向量。

欲使投影后的总方差最大，即λ\lambdaλ最大，因此最佳的投影向量www是特征值λ\lambdaλ最大时对应的特征向量，因些，当我们将www设置为与具有最大的特征值向量相等时，方差会达到最大值。这个特征向量被称为第一主成分。

我们可以用一种增量的方式定义额外的主成分，方法为：在所有与那些已经考虑过的方向正交的所有可能的方向中，将新的方向选择为最大化投影方差的方向。如果我们考虑kkk维投影空间的一般情形，那么最大化投影数据方差的最优线性投影由数据协方差矩阵Σ\SigmaΣ的kkk个特征向量w1,...,wkw_1, ..., w_kw1,...,wk定义，对应于kkk个最大的特征值 λ1,...,λk\lambda_1,...,\lambda_kλ1,...,λk。可以通过归纳法很容易地证明出来。

因此，我们只需要对协方差矩阵进行特征值分解，得到的前kkk大特征值对应的特征向量就是最佳的kkk维新特征，而且这kkk维新特征是正交的。得到前kkk个uuu以后，原始数据集XXX通过变换可以得到新的样本。

1.2 最小平方误差理论

如上图所示，假设有这样的二维样本点（红色点），按照前文我们讲解的最大方差理论，我们的目标是是求一条直线，使得样本点投影到直线或者平面上的点的方差最大。本质是求直线或者平面，那么度量直线求的好不好，不仅仅只有方差最大化的方法。再回想我们最开始学习的线性回归等，目的也是求一个线性函数使得直线能够最佳拟合样本点，那么我们能不能认为最佳的直线就是回归后的直线呢？回归时我们的最小二乘法度量的是样本点到直线的坐标轴距离。比如这个问题中，特征是xxx，类标签是yyy。回归时最小二乘法度量的是距离ddd。如果使用回归方法来度量最佳直线，那么就是直接在原始样本上做回归了，跟特征选择就没什么关系了。

因此，我们打算选用另外一种评价直线好坏的方法，使用点到直线的距离d′d′d′来度量。

现在有mmm个样本点x(1),...,x(m)x^{(1)},...,x^{(m)}x(1),...,x(m)，每个样本点为nnn维。将样本点x(i)x^{(i)}x(i)在直线上的投影记为x(1)′′x^{(1)''}x(1)′′，那么我们就是要最小化∑i=1m(x(i)′−x(i))2\sum_{i=1}^m (x^{(i)'}-x^{(i)})^2i=1∑m(x(i)′−x(i))2

这个公式称作最小平方误差(Least Squared Error)。

首先，我们确定直线经过的点，假设要在空间中找一点x0x_0x0来代表这mmm个样本点，“代表”这个词不是量化的，因此要量化的话，我们就是要找一个nnn维的点x0x_0x0，使得J0(x0)=∑i=1m(x(i)′−x(i))2J_0(x_0)=\sum_{i=1}^m (x^{(i)'}-x^{(i)})^2J0(x0)=i=1∑m(x(i)′−x(i))2最小。其中J0(x0)J_0(x_0)J0(x0)是平方错误评价函数。假设x‾\overline{x}x为mmm个样本点的均值，即x‾=1m∑i=1mx(i)\overline{x}=\frac{1}{m}\sum_{i=1}^m x^{(i)}x=m1i=1∑mx(i)
则J0(x0)=∑i=1m(x0−x(i))2=∑i=1m((x0−x‾)−(x(i)−x‾))2=∑i=1m(x0−x‾)2−2∑i=1m(x0−x‾)T(x(i)−x‾)+∑i=1m(x(i)−x‾)2=∑i=1m(x0−x‾)2−2(x0−x‾)T∑i=1m(x(i)−x‾)+∑i=1m(x(i)−x‾)2=−∑i=1m(x0−x‾)2+∑i=1m(x(i)−x‾)2J_0(x_0)=\sum_{i=1}^m(x_0-x^{(i)})^2\\=\sum_{i=1}^m((x_0-\overline{x})-(x^{(i)}-\overline{x}))^2\\=\sum_{i=1}^m(x_0-\overline{x})^2-2\sum_{i=1}^m(x_0-\overline{x})^T(x^{(i)}-\overline{x})+\sum_{i=1}^m(x^{(i)}-\overline{x})^2\\=\sum_{i=1}^m(x_0-\overline{x})^2-2(x_0-\overline{x})^T\sum_{i=1}^m(x^{(i)}-\overline{x})+\sum_{i=1}^m(x^{(i)}-\overline{x})^2\\=-\sum_{i=1}^m(x_0-\overline{x})^2+\sum_{i=1}^m(x^{(i)}-\overline{x})^2J0(x0)=i=1∑m(x0−x(i))2=i=1∑m((x0−x)−(x(i)−x))2=i=1∑m(x0−x)2−2i=1∑m(x0−x)T(x(i)−x)+i=1∑m(x(i)−x)2=i=1∑m(x0−x)2−2(x0−x)Ti=1∑m(x(i)−x)+i=1∑m(x(i)−x)2=−i=1∑m(x0−x)2+i=1∑m(x(i)−x)2

显然，上式的第二项与x0x_0x0无关，因此，J0(x0)J_0(x_0)J0(x0)在x‾\overline{x}x处有最小值。

接下来，我们确定直线的方向向量。我们已经知道直线经过点x‾\overline{x}x，假设直线的方向是单位向量e→\overrightarrow{e}e。那么直线上任意一点x(i)′x^{(i)'}x(i)′有：x(i)′=x‾+aie→x^{(i)'}=\overline{x}+a_i \overrightarrow{e}x(i)′=x+aie,其中，aia_iai是x(i)′x_{(i)'}x(i)′到点x‾\overline{x}x的距离。

我们重新定义最小平方误差：

这时候点都聚集在新的坐标轴周围，因为我们使用的最小平方误差的意义就在此。另外，PRML书上从线性子空间的角度进行了详细的阐述，有兴趣的读者可以看看。

优点：

它是无监督学习，完全无参数限制的。在PCA的计算过程中完全不需要人为的设定参数或是根据任何经验模型对计算进行干预，最后的结果只与数据相关，与用户是独立的。

用PCA技术可以对数据进行降维，同时对新求出的“主元”向量的重要性进行排序，根据需要取前面最重要的部分，将后面的维数省去，可以达到降维从而简化模型或是对数据进行压缩的效果。同时最大程度的保持了原有数据的信息。

各主成分之间正交，可消除原始数据成分间的相互影响。

计算方法简单，易于在计算机上实现。

缺点：

如果用户对观测对象有一定的先验知识，掌握了数据的一些特征，却无法通过参数化等方法对处理过程进行干预，可能会得不到预期的效果，效率也不高。

贡献率小的主成分往往可能含有对样本差异的重要信息。

特征值矩阵的正交向量空间是否唯一有待讨论。

在非高斯分布的情况下，PCA方法得出的主元可能并不是最优的，此时在寻找主元时不能将方差作为衡量重要性的标准。

1.3 高维数据下的特征值分解

在现实生活中，由于样本的高维特征（图像识别每一个像素点都是一个维度），将会导致协方差矩阵Σ\SigmaΣ非常大，计算机难以存储与计算，那如何针对高维数据做特征值分解呢？我们知道Σ=XXT\Sigma=XX^TΣ=XXT，考虑替代矩阵P=XTXP=X^T XP=XTX，如果有100个样本，10000个维度，那么Σ\SigmaΣ就是一个10000维的方阵，而PPP只是一个100维的方阵。我们做如下推导：Pv=λvXTXv=λvXXTXv=XλvΣ(Xv)=λ(Xv)Pv=\lambda v\\X^TXv=\lambda v\\XX^TXv=X\lambda v\\\Sigma(Xv)=\lambda(Xv)Pv=λvXTXv=λvXXTXv=XλvΣ(Xv)=λ(Xv)
这样通过求PPP的特征值和特征向量，其特征值即Σ\SigmaΣ的特征值，其特征向量右乘随机向量XXX即为Σ\SigmaΣ的特征向量，从而完成高维数据的特征值分解（只能得到部分前几的主成分）。

求证：ATAA^TAATA与AATAA^TAAT的特征值相等？
补充：一般情况下，ABABAB与BABABA的特征值只差000特征值的个数。

若存在可逆矩阵PPP使得，P−1AP=BP^{-1}AP=BP−1AP=B，则矩阵AAA与矩阵BBB相似。若两矩阵相似，则它们的秩相等，行列式相等，迹（所有特征值的和，或主对角线元素的总和）相等，两者拥有同样的特征值（可能会多出一些为零的特征值）。

2. CCIPCA增量主元分析算法[1]

candid covariance-free incremental principle component analysis直接协方差无关的增量式主成分分析。

假设输入向量序列为u′(t),t=1,2,...u'(t),t=1,2,...u′(t),t=1,2,...;第nnn幅图像输入时均值为m(n)=1n∑t=1nu′(t)m(n)=\frac{1}{n}\sum_{t=1}^n u'(t)m(n)=n1∑t=1nu′(t),它的协方差矩阵为A(n)=1n∑t=1n[u′(t)−m(n)][u′(t)−m(n)]T=1n∑i=1nu(t)u(t)T(1)A(n)=\frac{1}{n}\sum_{t=1}^n[u'(t)-m(n)][u'(t)-m(n)]^T=\frac{1}{n}\sum_{i=1}^n u(t)u(t)^T \tag{1}A(n)=n1t=1∑n[u′(t)−m(n)][u′(t)−m(n)]T=n1i=1∑nu(t)u(t)T(1)这里u(t)=u′(t)−m(n)u(t)=u'(t)-m(n)u(t)=u′(t)−m(n)

u(n)u(n)u(n)的第iii个特征值和特征值和特征向量的计算公式为λixi(n)=A(n)xi(n)\lambda_ix_i(n)=A(n)x_i(n)λixi(n)=A(n)xi(n),其中xi(n)x_i(n)xi(n)为第nnn输入时待求的第iii个特征向量，λi\lambda_iλi为对应的特征值。CCIPCA算法为了加快迭代的速度，整个迭代是对特征值和特征向量的乘积λixi\lambda_ix_iλixi进行的，设第nnn个输入时有vi(n)=λixi(n)=A(n)xi(n)(2)v_i(n)=\lambda_ix_i(n)=A(n)x_i(n) \tag{2}vi(n)=λixi(n)=A(n)xi(n)(2)
把（1）式代入（2）式，可得vi(n)=1n∑t=1nu(t)u(t)Txi(n)(3)v_i(n)=\frac{1}{n}\sum_{t=1}^n u(t)u(t)^T x_i(n)\tag{3}vi(n)=n1t=1∑nu(t)u(t)Txi(n)(3)
若通过迭代获得特征值和特征向量的乘积viv_ivi，因特征向量是归一的，只要对（2）式求模（内积，开根），可求得λi=∣∣vi∣∣,xi=vi∣∣vi∣∣\lambda_i=||v_i||,x_i=\frac{v_i}{||v_i||}λi=∣∣vi∣∣,xi=∣∣vi∣∣vi。迭代采用（3）式，把vi(n−1)∣∣vi(n−1)∣∣\frac{v_i(n-1)}{||v_i(n-1)||}∣∣vi(n−1)∣∣vi(n−1)近似为xi(n)x_i(n)xi(n)代入（3）式，经变换可得CCIPCA的基本迭代公式，vi(n)=n−1nvi(n−1)+1nu(n)u(n)Tvi(n−1)∣∣vi(n−1)∣∣(4)v_i(n)=\frac{n-1}{n}v_i(n-1)+\frac{1}{n}u(n)u(n)^T\frac{v_i(n-1)}{||v_i(n-1)||} \tag{4}vi(n)=nn−1vi(n−1)+n1u(n)u(n)T∣∣vi(n−1)∣∣vi(n−1)(4)

其中n−1n\frac{n-1}{n}nn−1为上一步的迭代值 vi(n−1)v_i(n-1)vi(n−1)的权值，第2项的1n\frac{1}{n}n1相当于迭代的调整步长。u(n)u(n)u(n)作为第nnn幅新输入数据对迭代向量vi(n)v_i(n)vi(n)的调整，在迭代中vi(n)v_i(n)vi(n)逐步收敛到所求的第iii个特征向量。对不同序号的特征向量，都可以用（4）式迭代，只是输入的向量u(n)u(n)u(n)不同。求最大的特征值对应的特征向量时，u(n)u(n)u(n)为机器人眼睛直接采到的第nnn个数据。在求第2，第3乃至更高维特征向量时，须作以下处理，已经通过迭代得到第1个特征向量，先设u1(n)=u(n)u_1(n)=u(n)u1(n)=u(n)，并把u1(n)u_1(n)u1(n)投影到上一个已经求到的特征向量上（现为第1个特征向量），求出残差图像u2(n)u_2(n)u2(n)，如下式表示，u2(n)=u1(n)−u1T(n)v1(n)∣∣v1(n)∣∣v1(n)∣∣v1(n)∣∣u_2(n)=u_1(n)-u_1^T(n)\frac{v_1(n)}{||v_1(n)||}\frac{v_1(n)}{||v_1(n)||}u2(n)=u1(n)−u1T(n)∣∣v1(n)∣∣v1(n)∣∣v1(n)∣∣v1(n)
u2(n)u_2(n)u2(n)作为求第2个特征向量的输入，类似的可求出第3，第4，…个特征向量。因残差图像和上1个特征向量所恢复的图像正交，从而可求出所有相互正交的特征向量。另外，每输入1幅新的数据时，均值也要更新，对（1）式，输入第nnn幅图像时的均值采用如下迭代式，m(n)=n−1nm(n−1)+1nu′(n)m(n)=\frac{n-1}{n}m(n-1)+\frac{1}{n}u'(n)m(n)=nn−1m(n−1)+n1u′(n)

3. KPCA

一篇关于KPCA的超级全的博客

原始PCA
XXTwi=λiwiXX^T w_i=\lambda_i w_iXXTwi=λiwi
利用函数Φ\PhiΦ将数据XXX映射到高维空间Φ(X)\Phi(X)Φ(X)
Φ(X)Φ(X)Twi=λiwi\Phi(X)\Phi(X)^T w_i=\lambda_i w_iΦ(X)Φ(X)Twi=λiwi
令wi=∑k=1NαiΦ(xi)=Φ(X)αw_i = \sum_{k=1}^N \alpha_i\Phi(x_i) = \Phi(X)\alphawi=∑k=1NαiΦ(xi)=Φ(X)α得
Φ(X)Φ(X)TΦ(X)α=λiΦ(X)α\Phi(X)\Phi(X)^T \Phi(X) \alpha=\lambda_i \Phi(X)\alphaΦ(X)Φ(X)TΦ(X)α=λiΦ(X)α
两边左乘Φ(X)T\Phi(X)^TΦ(X)T得
Φ(X)TΦ(X)Φ(X)TΦ(X)α=λiΦ(X)TΦ(X)α\Phi(X)^T\Phi(X)\Phi(X)^T \Phi(X) \alpha=\lambda_i \Phi(X)^T \Phi(X)\alphaΦ(X)TΦ(X)Φ(X)TΦ(X)α=λiΦ(X)TΦ(X)α
由核函数的性质可得K=Φ(X)TΦ(X)K=\Phi(X)^T \Phi(X)K=Φ(X)TΦ(X)
KKα=λiKα⇒Kα=λiαKK \alpha=\lambda_i K\alpha \\ \Rightarrow K\alpha=\lambda_i \alphaKKα=λiKα⇒Kα=λiα

此处注意，特征向量wiw_iwi应该进一步由Φ(X)α\Phi(X)\alphaΦ(X)α计算得出。但是，Φ(X)\Phi(X)Φ(X)是无法计算的，从而也就没办法求得特征向量。

此处，我们分析是否可以直接利用α\alphaα来对数据进行降维呢，答案是肯定的。

假如新来一个数据xxx，需要利用以上结果来得到降维的数据x^\hat{x}x^：

x^=wiTx=(Φ(x)α)TΦ(x)=αTΦ(X)TΦ(x)=[α1,...,αN][k(x1,x),...,k(xN,x)]T\hat{x}=w_i^T x\\=(\Phi(x)\alpha)^T \Phi(x)\\=\alpha^T \Phi(X)^T\Phi(x)\\=[\alpha_1,...,\alpha_N][k(x_1,x),..., k(x_N, x)]^Tx^=wiTx=(Φ(x)α)TΦ(x)=αTΦ(X)TΦ(x)=[α1,...,αN][k(x1,x),...,k(xN,x)]T

PCA与KPCA的区别

PCA 仍然不失为一种好分析方法，数据呈非线性流形分布，或者说是各指标呈非线性关系时，对于线性分析分析方法来说可能效果不是特别好，但同时应注意的是它也是一种统计分析方法，实际经济指标中都存在线性相关性（信息冗余），这是符合统计规律的，完全不相关的经济数据是极其少见，也就是说要求的数据只要大致呈线性分布，而且有PCA计算简单，无需先验知识，无需参数设置等优点。

PCA与线性核KPCA不完全一样，对于有P个指标的N个数据样本。PCA计算协方差阵为PXP维矩阵，它可以提取的主成分数为P，而KPCA是从核矩阵出发计算的，最大可以提取的主成分为N，为满足在特征空间中的样本均值为零，还要对核矩阵K进行特殊处理，这也是导致与线性核与原样本内积的不一致的原因。

KPCA核函数与核参数难于选择。PCA的协方差矩阵的特征向量对应于各指标的主成分比重，从而能用原指标云解释主成分，而KPCA的是基于核矩阵的特征向量，与原指标没有对应关系，从而核主成分解释困难，其次KPCA将指标投影到高维特征空间后，而其实际的数据又是在原空间处理的，其数值在原空间中是否均有排序意义也值得进一步研究。

4. 2DPCA[2]

设Ai,i=1,2,...,NA_i, i=1,2,...,NAi,i=1,2,...,N是NNN个样本图像，Ai∈Rm×nA_i\in R^{m\times n}Ai∈Rm×n。首先计算协方差矩阵，也就是总体散布矩阵GtG_tGt Gt=1N∑i=1N(Ai−A^)T(Ai−A^)G_t=\frac{1}{N}\sum_{i=1}^N(A_i-\hat{A})^T(A_i-\hat{A})Gt=N1i=1∑N(Ai−A^)T(Ai−A^)其中A^=1N∑i=1NAi\hat{A}=\frac{1}{N}\sum_{i=1}^N A_iA^=N1∑i=1NAi为全体样本的均值。

计算GtG_tGt的特征值和特征向量，取特征值累积项献率α=0.90.99\alpha=0.9~0.99α=0.9 0.99所对应的特征向量组成投影矩阵U=[u1,u2,...,uk]∈Rn×kU=[u_1, u_2, ..., u_k]\in R^{n\times k}U=[u1,u2,...,uk]∈Rn×k 。则，Fi=AiU∈Rm×kF_i=A_iU\in R^{m\times k}Fi=AiU∈Rm×k就是AjA_jAj的特征。可知原来二维的图像大小m×nm\times nm×n，现在降维为m×km\times km×k，kkk是根据α\alphaα来确定的。也就是说，实行特征提取后，只是压缩了图像矩阵列向量的位数，行向量位数不变。

5. 2D2DPCA

同于2DPCA只在列方向降了维数，降维的效果不理想。为了更好的降维，D.Q.Zhang和Z.H.Zhou提出了双向的二维主成分分析方法(2D2DPCA)，也就是在行和列两个方向都进行2DPCA处理。

对所有训练样本利用上述2DPCA处理之后得到新的训练样本Fi,i=1,2,...,NF_i, i=1,2,...,NFi,i=1,2,...,N，其中Yi∈Rm×kY_i\in R^{m\times k}Yi∈Rm×k。在新样本上构造协方差矩阵Gt∗G_t^*Gt∗ Gt∗=1N∑i=1N(Fi−F^)(Fi−F^)TG_t^*=\frac{1}{N}\sum_{i=1}^N (F_i-\hat{F})(F_i-\hat{F})^TGt∗=N1i=1∑N(Fi−F^)(Fi−F^)T

同理，求Gt∗G_t^*Gt∗的特征向值与特征向量，取特征值累计项献率为α\alphaα的特征向量，得到行方向的投影矩阵V=[v1,v2,...,vd]∈Rm×dV=[v_1,v_2,...,v_d] \in R^{m\times d}V=[v1,v2,...,vd]∈Rm×d。FiF_iFi的投影为VTFi∈Rd×kV^TF_i\in R^{d\times k}VTFi∈Rd×k。

至此，两个投影方向的最优投影矩阵UUU和VVV都求得了，对于图像Ai∈Rm×nA_i\in R^{m\times n}Ai∈Rm×n最终降维矩阵为Yi=VTAiU∈Rd×kY_i=V^TA_iU\in R^{d\times k}Yi=VTAiU∈Rd×k，其重构图像为Aiapprox=VYiUTA_i^{approx}=VY_iU^TAiapprox=VYiUT

6. BDPCA

BDPCA分别对行和列方向进行数据降维。
将图像样本矩阵Xi′X_{i}^{\prime}Xi′分解成ppp个1×q1\times q1×q的行向量，则行方向总体散度矩阵为:

Sr=1np∑i=1n[Xi′−X‾(n)]T[Xi′−X‾(n)]S_{r}=\frac{1}{n p} \sum_{i=1}^{n}\left[X_{i}^{\prime}-\overline{X}(n)\right]^{\mathrm{T}}\left[X_{i}^{\prime}-\overline{X}(n)\right] Sr=np1i=1∑n[Xi′−X(n)]T[Xi′−X(n)]

其中，X‾(n)\overline{X}(n)X(n)为图像样本矩阵的均值。取前krk_rkr个最大特征值所对应的特征向量组成行方向投影矩阵为
Wr=[wr1,wr2,⋯,wrkr]\boldsymbol{W}_{r}=\left[w_{r 1}, w_{r 2}, \cdots, w_{r k_{r}}\right] Wr=[wr1,wr2,⋯,wrkr]

同理，列方向总体散度矩阵及前kck_ckc个最大特征值所对应的特征向量组成列方向投影矩阵分别为
Sc=1nq∑i=1n[Xi′−X‾(n)][Xi′−X‾(n)]TS_{c}=\frac{1}{n q} \sum_{i=1}^{n}\left[X_{i}^{\prime}-\overline{X}(n)\right]\left[X_{i}^{\prime}-\overline{X}(n)\right]^{\mathrm{T}} Sc=nq1i=1∑n[Xi′−X(n)][Xi′−X(n)]T
Wc=[wc1,wc2,⋯,wckc]\boldsymbol{W}_{c}=\left[w_{c 1}, w_{c 2}, \cdots, w_{c k_{c}}\right] Wc=[wc1,wc2,⋯,wckc]
图像样本矩阵Xi′X_{i}^{\prime}Xi′所对应的特征矩阵为
Y=WcTXi′WcY=W_{c}^{\mathrm{T}} X_{i}^{\prime} W_{c} Y=WcTXi′Wc

BDPCA的特征矩阵维数仅为kc×krk_c\times k_rkc×kr，因此其运算量要远小于CCIPCA。但BDPCA是2维的批处理计算。

参考文献

[1]Juyang Weng, Yilu Zhang, and Wey-Shiuan Hwang, “Candid covariance-free incremental principal component analysis,” IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 25, no. 8, pp. 1034–1040, Aug. 2003.
[2] Jian Y, David Z, Frangi A F, et al. Two-dimensional PCA: a new approach to appearance-based face representation and recognition.[J]. IEEE Transactions on Pattern Analysis & Machine Intelligence, 2004, 26(1):131-137.