主成分分析

前言

好久没再“动笔”来写点东西了，看了好多资料，但是总感觉不成体系，四处开花却没有比较深入的研究。所以，终究还是觉得，写下一些札记，这样至少能够在需要的时候快速回忆起来。
当然，先写PCA也是因为绕不开最近看的内容——降维。总想在深度学习这个领域中多研究点东西，那就从基本的开始吧。

主成分分析的背景来源

主成分分析(Principal components analysis，PCA)是一种分析、简化数据集的技术，来源于多元统计分析，说简单点主要就是做数据降维的。特别是在机器学习领域当中，经常面临维数灾难(Curse of Dimensionality)的问题。所以，往往无论在实际应用还是用于科学研究中，都会首先通过各种技术来对数据进行降维。而PCA就是一种成熟的，并且应用广泛的数据降维技术。

PCA的理论基础

两种解释的理论推导都走了不少弯路，主要参照了资料
[1] PRML读书会第十二章 Continuous Latent Variables

符号说明

设：样本数量为$N$，样本的维数为$D$，对应样本矩阵$X=(x_1,x_2,…,x_D)$。
假设样本均值满足中心化条件
$$\overline{x}=\frac{1}{N}\sum_{i=1}^N x_i=0$$

样本协方差
$$S=\frac{1}{N}\sum_{i=1}^N (x_i-\overline{x})(x_i-\overline{x})^T=\frac{1}{N}\sum_{i=1}^N x_i x_i^T$$

最大方差解释

这部分的解释源于最大方差的解释关键在于，投影之后样本之间方差的最大化。
假设投影矩阵$W=(w_1,…,w_D)$，那么$X$投影之后的样本矩阵为$Y=W^TX=(y_1,…,y_D)$。
就单个样本$x_i$的投影$y_i$
$$y_i=W^Tx_i$$

所有样本在某一方向$w_k$上的方差
\begin{align}
\frac{1}{N}\sum_{i=1}^N (w_k^Tx_i-w_k^T\overline{x})^2=\frac{1}{N}\sum_{i=1}^N (w_k^Tx_i)^2=w_k^T S w_k
\end{align}

事实上，根据
$$y_i=W^Tx_i
=\begin{pmatrix}
w_1^T \\
w_2^T \\
\vdots \\
w_D^T
\end{pmatrix}x_i
=\begin{pmatrix}
w_1^T x_i \\
w_2^T x_i \\
\vdots \\
w_D^T x_i
\end{pmatrix}$$

$w_k^Tx_i$的含义是投影后的样本$y_i$在第$k$维的值。
为不失一般性，我们设$w_k^Tw_k=1$，然后我们可以明确目标函数和条件，即构造正交变换$w_k$使得样本方差最大化：
\begin{align}
\max & \sum_{k=1}^D w_k^T S w_k \\
s.t. & w_k^Tw_k=1, k=1,…,D
\end{align}

这里直接用拉格朗日乘法求解
$$F=\sum_{k=1}^D w_k^T S w_k-\lambda\sum_{k=1}^D(w_k^Tw_k-1)$$

对$w_k$求偏导并等于0，得到
$$S w_k=\lambda w_k$$

最小平方误差解释

最小平方误差的意义在于，寻找到一个投影$W$，使得所有样本点到直线的距离之和最短。这部分内容也主要参考[1]。
假设$U=(u_1,…,u_D)$是$D$维空间的完备基向量，那么对于任意$D$维样本$x_k$，我们可以使用基向量$U$线性表示
$$x_k=\sum_{i=1}^D c_{ki} u_i$$
这里有个我想了很久都不知道怎么证明，却忽略了一个关键的事实——$U$中的向量彼此正交，即
$$u_i^T u_j=
\begin{cases}
1,&i=j \\
0,&i \not=j
\end{cases}$$

过程如下
$$u_i^T x_k =\sum_{j=1}^D c_{kj} (u_i^T u_j)=c_{ki}$$

然而我们并不想用$D$个$D$维向量来表示样本，而是希望用$d(<D)$个$D$维向量来表示。于是这样就有了对样本的估计
\begin{align}
\widehat{x}_k &=\sum_{i=1}^d (u_k^T x_k) u_i \\
\Rightarrow x_k-\widehat{x}_k &= \sum_{i=d+1}^D (u_k^T x_k) u_i
\end{align}

根据最小均方误差得到
\begin{align}
J &=\sum_{k=1}^N | x_k-\widehat{x}_k |^2 \\
&=\sum_{k=1}^N (x_k-\widehat{x}_k)^T(x_k-\widehat{x}_k) \\
&=\sum_{k=1}^N \left( \sum_{i=d+1}^D c_{ki} u_i \right)^T\left( \sum_{j=d+1}^D c_{kj} u_j \right) \\
&=\sum_{k=1}^N \left( \sum_{i=d+1}^D \sum_{j=d+1}^D (c_{ki}c_{kj}) u_j^T u_i \right) \\
&=\sum_{k=1}^N \left( \sum_{i=d+1}^D c_{ki}^2 \right)
\end{align}

又且
$$c_{ki}^2=u_i^T x_k u_i^T x_k=u_i^T(x_k x_k^T)u_i$$

故得到
\begin{align}
J&=\sum_{k=1}^N \sum_{i=d+1}^D u_i^T(x_k x_k^T)u_i \\
&=\sum_{i=d+1}^D \sum_{k=1}^N u_i^T(x_k x_k^T)u_i \\
&=\sum_{i=d+1}^D u_i^T \left( \sum_{k=1}^N (x_k x_k^T) \right) u_i \\
&=\sum_{i=d+1}^D u_i^T S u_i \\
\end{align}

实际上依然可以看出这是一个优化问题，并且解法跟最大方差解释的解法是相似的：
\begin{align}
\min &\sum_{i=d+1}^Du_i^T S u_i \\
s.t. &u_i^T u_i=1,i=1,…,D
\end{align}

其他说明

对于PCA两种解释的几何意义可参照[2]和[3]，[4]给予了实际的例子来演示，过程相当的清晰，然而这三者的理论过程不如[1]严谨。
[2] 主成分分析（Principal components analysis）-最大方差解释
[3] 主成分分析（Principal components analysis）-最小平方误差解释
[4] PCA的数学原理
另外一些需要说明的是关于PCA的使用范围(很遗憾，对这个内容能够参考的文献比较少)。很多文献都说，PCA对于非高斯分布(也有说是对指数族分布)的数据处理起来效果相当不理想，至于原因为何就不知道了，不过对于这种情况可以考虑使用Kernel PCA。
另外有一些对于PCA降维的叙述：
[5] 浅议PCA降维

Kernel PCA

实际是上，对于Kernel PCA的引入是一种很自然的事情。PCA的线性变换过程，决定了PCA对很多数据的降维实际上时无能为力的，特别是一些高阶数据。这个时候就可以考虑引入核函数，通过变换将原本在特征空间难于线性分类的数据，映射到在高维空间中线性可分的数据。
事实上，以下叙述的内容主要来源于一下资料：
[5] Kernel PCA 原理和演示

Kernel PCA的理论基础

我们的基本想法是，引入线性函数$\phi$，将输入空间中的样本点$x_1,…,x_D$变换为高维空间中的样本点$\phi(x_1),…,\phi(x_D)$，然后在此基础上使用原始PCA过程。但是，这个过程面临两个问题：
1.在核函数的理论中，很难或者很少直接给出映射函数$\phi$，一般都是通过内积形式给出，即给定一个核函数$K(x_i,x_j)$
2.PCA的使用要求输入数据是中心化的(centred)
下面两节依次来解决这两个问题。

使用核函数来在高维空间做PCA

这里首先假定映射后的数据是满足中心化的条件，即：
$$\sum_{i=1}^N \phi(x_i)=0$$

那么，其对应的协方差矩阵为
$$C=\frac{1}{N}\sum_{i=1}^N \phi(x_i)\phi(x_i)^T$$

事实上，我们的本意是对矩阵$C$使用原始的PCA过程，但是我们并不知道矩阵$C$，所以我们希望通过内积形式来寻找一个矩阵来替代它。
假设矩阵$C$的特征值和特征向量为$\lambda_k$和$v_k$，满足
$$Cu_k=\lambda_k u_k$$

首先利用上式子将$C$展开
\begin{align}
\lambda_k u_k &= Cu_k \\
\lambda_k u_k &= \left( \frac{1}{N}\sum_{i=1}^N \phi(x_i)\phi(x_i)^T \right)u_k \\
\Rightarrow u_k &= \frac{1}{N\lambda_k}\sum_{i=1}^N \phi(x_i) \left( \phi(x_i)^T u_k \right) \\
\Rightarrow u_k &= \sum_{i=1}^N \left( \frac{\phi(x_i)^T u_k}{N\lambda_i} \right) \phi(x_i) \\
\Rightarrow u_k &= \sum_{i=1}^N \alpha_i^k \phi(x_i)
\end{align}

其中$\alpha_i^k=\frac{\phi(x_i)^T u_k}{N\lambda_k}$，此时$u_k$可以由$\phi(x_i)$表示

其次在两边同时乘以向量$\phi(x_k)^T$
\begin{align}
\lambda_k u_k &=Cu_k \\
\lambda_k u_k &= \left( \frac{1}{N}\sum_{i=1}^N \phi(x_i)\phi(x_i)^T \right)u_k \\
\Rightarrow\lambda_k \phi(x_k) u_k &=\phi(x_k) \left( \frac{1}{N}\sum_{i=1}^N \phi(x_i)\phi(x_i)^T \right)u_k
\end{align}

最后再将得到的$u_k=\sum_{j=1}^N \alpha_j^k \phi(x_j)$展开(注意下标)
\begin{align}
\lambda_k \phi(x_k)^T \sum_{i=1}^N \alpha_i^k \phi(x_i) &=\phi(x_k)^T \left( \frac{1}{N}\sum_{i=1}^N \phi(x_i)\phi(x_i)^T \right)\sum_{j=1}^N \alpha_j^k \phi(x_j) \\
\Rightarrow \lambda_k N \sum_{i=1}^N \alpha_i^k \phi(x_k)^T\phi(x_i) &=\sum_{i=1}^N \phi(x_k)^T\phi(x_i) \sum_{j=1}^N \alpha_j^k \phi(x_i)^T\phi(x_j) \\
\Rightarrow \lambda_k N \sum_{i=1}^N \alpha_i^k K_{ij} &=\sum_{i=1}^N K_{ki} \sum_{j=1}^N \alpha_j^k K_{ij}
\end{align}

其中$K_{ij}=\phi(x_i)^T\phi(x_j)$。
由于这里$\forall k \in {1,…,n}$都成立，故得到
\begin{align}
\lambda_k N K \alpha^k &=K^2 \alpha^k \\
\Rightarrow K\alpha^k &=\widehat{\lambda}_k\alpha^k
\end{align}

其中$\widehat{\lambda}_k=\lambda_k N$(这里我推导了很久，要仔细注意等式$\forall k \in {1,…,n}$都成立就明白了)。
显然我们似乎已经获得了我们想要的结论：矩阵$K$的元素$K_{ij}$是由映射后的数据$\phi(x_i)$和$\phi(x_j)$的内积得到的，而核函数作用恰好体现在这里，所以紧接着的需要做的事情就是如何求矩阵$K$了。

Centering高维空间的数据

假设$\phi_i = \phi(x_i)$，$\phi_i^C = \phi_i – \frac{1}{N}\sum_{k=1}^N \phi_k$，得：
\begin{align}
K_{ij}^C
&= \phi_i^C \cdot \phi_j^C \\
&= (\phi_i – \frac{1}{N}\sum_{k=1}^N \phi_k)^T (\phi_j – \frac{1}{N}\sum_{l=1}^N \phi_l ） \\
&=\phi_i^T\phi_j – \frac{1}{N}\sum_{l=1}^N \phi_i^T \phi_l – \frac{1}{N}\sum_{k=1}^N \phi_k^T \phi_j + \frac{1}{N^2}\sum_{k=1}^N \sum_{l=1}^N \phi_k^T \phi_l \\
&=K_{ij} –\frac{1}{N}\sum_{l=1}^N K_{il} – \frac{1}{N}\sum_{k=1}^N K_{kj} + \frac{1}{N^2}\sum_{k=1}^N \sum_{l=1}^N K_{kl}
\end{align}

这样，就完成了高纬空间数据的中心化。

数据投影

与原始PCA相似，利用特征方程$K\alpha^k =\widehat{\lambda}_k\alpha^k$求得特征值和特征向量，并选择其中的一部分即可做降维。但是有些许不同是在最后数据投影时依然使用了内积，具体过程如下：
矩阵$K$的特征值$\widehat{\lambda}_k$和特征向量$\alpha^k$相对应，选取其中$t$个组成投影矩阵$\alpha$
$$\alpha=(\alpha^1,…,\alpha^t)^T$$

原始数据$x_i$降维后的结果等于$\phi(x_i)$在$\alpha$上的投影
\begin{align}
\alpha\phi(x_i)^T &=(\alpha^1,…,\alpha^t)^T \phi(x_i)\\
&=\left(\sum_{j=1}^N \alpha_j^1 \phi(x_j)^T \phi(x_i),…,\sum_{l=1}^N \alpha_l^t \phi(x_l)^T \phi(x_i)\right) \\
&=\left(\sum_{j=1}^N \alpha_j^k K_{ji},…,\sum_{l=1}^N \alpha_l^t K_{li}\right)
\end{align}

实验

占座用，以后补充。

参考资料

[6] Max Wellings’s notes for Kernel Principal Component Analaysis
[7] 核PCA特征提取方法及其应用研究_高绪伟