时间戳

首页 > 技术文章 > 简单介绍SVD和PCA

简单介绍SVD和PCA

HanDoMind 2020-05-20 11:39 原文

Research 我唯唯诺诺

划水摸鱼 Tutorial 我重拳出击。

结构如下：首先是大致叙述下SVD，然后介绍PCA，再对SVD给出一个大概的推导思路，最后是一些发散内容。本文仅对实数矩阵讨论，复数的情况均可以简单推广得到。数学推导相关的内容均参考于 Gilbert Strang 的 Linear Algebra and Its Application。

简单介绍

SVD 即奇异值分解(Singular Value Decomposition)，在本文不考虑复数矩阵的前提下即对于任意 \(m\times n\) 实矩阵 \(A_{m\times n}\)，我们都可以将其表示为 \(A_{m\times n} = U_{m\times m} \Sigma_{m\times n} V_{n\times n}^{\mathrm{T}}\)，其中，\(U_{m\times m}\)和\(V_{n\times n}\)都是正交矩阵。对于矩阵\(\Sigma_{m\times n}\)，我们将其对角线上第 \(i\) 个元素记作\(\sigma_i\) ，即矩阵中第 \(i\) 行第 \(i\) 列的元素，这个矩阵仅在这些位置可能不等于0。更具体的说，假设 \(r = rank(A_{m\times n})\)，那么 \(\Sigma_{m\times n}\) 有且仅有 \(r\) 个非零元。一般我们将这些元素在对角线上按从大到小的顺序排列，也即 \(\sigma_1 \ge \sigma_2 \ge \sigma_3 \ge \sigma_r > 0\)。

那么，这样的分解有什么好处呢？明明\(V_{n\times n}^{\mathrm{T}}\) 和 \(V_{n\times n}\) 都是正交矩阵，为什么偏偏要写成转置的样子呢？

为了让大家对SVD 有confidence，我们先来对这两个问题做一个简单的分析。

分析之前，我们先来引入一点记号。

关于实矩阵 \(A_{m\times n}\)，最朴素的观点就是将它看成\(m\times n\) 个实数的堆积。那么如果我们认为每一列元素看作一个 特殊的命运共同体 向量以后，我们就可以将这个矩阵看成 \(n\) 个 \(m\) 维的列向量横向排列的结果，假如我们把这些列向量按顺序记作 \(a_1,a_2,\dots,a_n\)，那么矩阵 \(A_{m\times n}\) 就可以写成 \(A_{m\times n} = [a_1,a_2,\dots,a_n]\)。

注：后文为了方便，在不引起歧义的情况下，我们便不对矩阵加下标了。

假设 \(D\) 为 \(r\times r\) 的对角矩阵，且 \(d_{ii} = \sigma_i\) ，则我们的奇异值分解就可以看做

\[\begin{aligned} A = & [u_1, u_2, \dots, u_m] \begin{pmatrix} D & O_{r\times (n-r)}\\ O_{(m-r)\times r} & O_{(m-r)\times (n-r)} \end{pmatrix} [v_1,v_2,\dots, v_n]^{\mathrm{T}}\\ = & [\sigma_1 u_1, \sigma_2 u_2, \dots, \sigma_r u_r, 0_{m\times 1},\dots, 0_{m\times 1}]_{m\times n} \begin{pmatrix} v_1^{\mathrm{T}}\\ v_2^{\mathrm{T}}\\ \vdots\\ v_n^{\mathrm{T}}\\ \end{pmatrix}_{n\times n}\\ = & \sigma_1 u_1 v_1^{\mathrm{T}} + \sigma_2 u_2 v_2^{\mathrm{T}} + \dots + \sigma_r u_r v_r^{\mathrm{T}}\\ \end{aligned} \]

对于实在纠结但是又不会的人提供此hint，第二个等式来自于名人名言“左乘矩阵就是行变换，右乘矩阵就是列变换”中的第二句，第三个等式严格说明需要把左边的矩阵拆成 n 个只有单列不一定为0的矩阵的和，右边的矩阵拆成m个只有单行不一定为0的矩阵的和，于是这个矩阵乘法变成了\(m\times n\) 个矩阵乘法的和，但是通过简单检查可以发现只有 \(u_iv_i^{\mathrm{T}}\) 的组合才可能产生非零矩阵，这个小结论对于一般的矩阵乘法也是成立的。

注意，此时 \(u_i v_i^{\mathrm{T}}\) 均为 \(m\times n\) 的矩阵，即我们将 \(A\) 写成了 \(r\) 个矩阵的线性组合。很明显凭直觉可以想到 \(\sigma_i\) 越大，\(\sigma_i u_i v_i^{\mathrm{T}}\) 包含的关于 \(A\) 的信息就越多，于是这就提供了一个很好的即给 \(A\) 降维又尽可能保持 \(A\) 的信息的工具。简单来讲我们可以利用上述表达对 \(A\) 做如下近似

\[A \approx \widetilde{A} = \sigma_1 u_1 v_1^{\mathrm{T}} + \sigma_2 u_2 v_2^{\mathrm{T}} + \dots + \sigma_k u_k v_k^{\mathrm{T}}, \quad k<r \]

基于这样的想法，各类降维的应用便可以基于SVD 展开了。

从SVD到PCA

PCA，即主成分分析，是一种数据科学中常用的数据降维手法。这里我们按照数据科学中常用的记号把数据记为\(X = [x_1, x_2, \cdots, x_n]^\top \in \mathbb{R}^{n\times d}\)，即\(X\)的第\(i\)行记录了样本数据\(x_i\)的转置，为叙述方便我们假设\(d \ll n\)。

在研究PCA前，我们从线性代数的角度看一个显然的等式\(X = X I_d\)，这一个等式告诉我们，\(x_{ij}\)其实记录了第\(i\)个样本关于第\(j\)个特征的权重，其中第\(j\)个特征在这里便是\(I_d\)的第\(j\)行，于是这也意味着，\(x_{ij}\)就是第\(j\)个特征的值，与常识相符合。

那么PCA想做什么呢，对于高维海量数据，其内部结构往往是冗余的，意味着我们可以用少量的具有"代表性"的向量，便可以将其中的主要信息抽出，数学上，就是找一组尽可能小的基，去提取尽可能多的信息。根据应用需求，我们会对基的数量或者保留信息的百分比做提出要求去进行计算，这里为了叙述方便，我们假设我们选取\(c\)个基，使得保留尽可能多的信息，描绘成数学语言便是在\(\mathbb{R}^d\)中找到\(c\)个向量 \(y_1, y_2, \cdots, y_c\) 使得 \(x_1, x_2, \cdots, x_n\) 在其张成的子空间上的投影损失信息最少，即最小化下述问题

\[\min_{y_1, y_2,\cdots,y_c} \sum_{i=1}^n \| x_i - P_Y(x_i) \|^2 \]

注意到投影本身是一个线性变换，于是上述问题也可以改写作【注意这里对\(X\)中的行向量投影，所以是右乘投影矩阵】

\[\min_Y \|X - XP_Y\|_F^2 \]

同时注意到，子空间\(Y\)和其对应的投影矩阵是相互决定的，结合投影矩阵的幂等性和对称性，我们得到问题有如下形式

\[\begin{align*} \min_P & \|X - XP\|_F^2 = \operatorname{tr}\big( (X - XP)(X-XP)^\top \big) = \operatorname{tr}(XX^\top + XPP^\top X^\top - XP^\top X^\top - XPX^\top)\\ \text{s.t. } & P^2 = P = P^\top \end{align*} \]

由投影矩阵的性质，我们有

\[\begin{aligned} \mathop{\arg\min}_P \|X - XP\|_F^2 = & \mathop{\arg\min}_P \operatorname{tr}(XX^\top + XPP^\top X^\top - XP^\top X^\top - XPX^\top) \\ = & \mathop{\arg\min}_P \operatorname{tr}(XPP^\top X^\top - XPP X^\top - XPPX^\top) \\ = & \mathop{\arg\min}_P -\operatorname{tr} \big(XP(XP)^\top\big) \end{aligned} \]

下面我们处理约束，由\(P^2 = P\) 知，\(P\) 可对角化，且特征值为0和1，于是存在正交矩阵\(Q\)使得

\[P = Q^\top \begin{pmatrix} E_c & O_{c\times(d - c)} \\ O_{(d-c)\times c} & O_{(d-c)\times (d-c)} \end{pmatrix} Q = Q_1^\top Q_1, \quad \text{其中 } Q = \begin{pmatrix} Q_1 \\ Q_2 \end{pmatrix}, Q_1 \in \mathbb{R}^{c \times d}, Q_2 \in \mathbb{R}^{(d-c)\times d} \]

同时，由正交性我们有\(Q_1 Q_1^\top = E_c\)，于是问题转化为

\[\begin{aligned} \min_{Q_1} \ & -\operatorname{tr} \big(XQ_1^\top Q_1(XQ_1^\top Q_1)^\top\big) = - \operatorname{tr} (X Q_1^\top Q_1 X^\top)\\ \text{s.t. }\ & Q_1 Q_1^\top = E_c \end{aligned} \]

记\(Q_1^\top = [q_1, q_2, \cdots, q_c]\)，即\(q_i^\top\)为\(Q_1\)的第\(i\)行，则我们有

\[\begin{aligned} - \operatorname{tr} (X Q_1^\top Q_1 X^\top) = - \operatorname{tr} (Q_1 X^\top X Q_1^\top) = -\sum_{i=1}^c q_i^\top (X^\top X) q_i \end{aligned} \]

假设\(X\)的奇异值分解为\(X = U \Sigma V^\top\)，则\(X^\top X = V \Sigma^2 V^\top\)，于是

\[- \operatorname{tr} (X Q_1^\top Q_1 X^\top) = -\sum_{i=1}^c q_i^\top V \Sigma^2 V^\top q_i = -\sum_{i=1}^c (\Sigma V^\top q_i)^\top (\Sigma V^\top q_i) \]

由于\(\|q_i\|=1\)，\(q_i\) 之间相互正交，我们有\(- \operatorname{tr}(XQ_1^\top Q_1 X^\top) \ge -\sum_{i=1}^{c} \sigma_i^2\)，其中\(\sigma_i\)为\(X\)第\(i\)大的奇异值，等式成立当且仅当\(q_i = v_i\)，于是 \(Q_1 = V_c^\top\)，其中\(V_c\) 是\(V\) 的前\(c\)列。于是PCA的结果\(XP = XQ_1^\top Q_1 = XV_c V_c^\top = X V_c(V_c^\top V_c)V_c^\top\) 相当于，将\(X\) 中的行向量向\(V_c\) 记录的子空间作投影，也即\(V_c\)便是我们需要找的最优子空间\(Y\)。

数学推导

我们按如下顺序来 mathematically 地推导一下 SVD。

  • 实对称矩阵可以做谱分解，即对于实对称矩阵 \(A\)，存在正交矩阵 \(Q\)，使得\(A = Q\Lambda Q^{\mathrm{T}}\)

    • 实对称矩阵的特征值都是实数
    • 实对称矩阵显然是正规矩阵（\(A^{\mathrm{T}}A = AA^{\mathrm{T}}\)）
    • 任意实数矩阵都能被正交矩阵上三角化，即存在正交矩阵\(Q\)，使得\(A = Q \Delta Q^{\mathrm{T}}\) （思想类似于Gram-Schmidt正交化，过程上类似于LU分解）
    • 对于\(A = Q \Delta Q^{\mathrm{T}}\)，如果\(A\)是正规矩阵，那么\(\Delta\) 也是正规矩阵
    • 如果一个矩阵既是正规矩阵，又是上三角阵，那么它一定是对角阵

  • 实对称正定矩阵的特征值大于0，实对称半正定矩阵的特征值大于等于0

  • 对于任意矩阵\(A\) ，\(A^{\mathrm{T}}A\) 均为实对称半正定矩阵

  • 因此，利用谱分解得到 \(A^{\mathrm{T}}A = V\Lambda V^{\mathrm{T}}\)，并取 \(\sigma_i = \sqrt{\lambda_i}\) ，我们便得到了SVD中的 \(\Sigma\) 矩阵和 \(V\) 矩阵，最后利用 \(AV = U\Sigma\) 即可推出 \(U\)

    回顾SVD的含义，注意到过程中 \(V\) 和 \(U\) 实际上都是存在冗余信息的，这种冗余信息的产生可以在SVD 的实践中清晰的看到原因，在此不予赘述。

发散内容

实际上，SVD 对矩阵信息进行了非常深刻的提纯，两个正交矩阵分别可以解读为\(U = [\text{列空间} | \text{左零空间}]\) ，\(V = [\text{行空间} | \text{零空间}]\)。

同时，SVD一个非常显著的优点是它的普适性，这种普适而有效的分解手段给许多分析都带来了很多便利。

最后，再给个 SVD 的应用：求最小二乘

有的时候最小二乘解不唯一，这个时候我们希望找到其中模最小的解，一种简单有效的方法就是利用SVD求出伪逆，再利用伪逆去计算最优最小二乘，方法如下

假设问题 \(\min \|y - Ax\|\) 的最优解不唯一，找出\(x^{*}=\arg\min \|x\|\)

  • \(A = U\Sigma V^{\mathrm{T}}\)
  • 按照\((\Sigma^+)_{ji} = (\Sigma)_{ij},\ i\ne j,(\Sigma^+)_{ii} = 1/(\Sigma)_{ii}\)定义\(\Sigma^+\)
  • 计算\(A^+ = V \Sigma^+ U^{\mathrm{T}}\)
  • 则 \(x^* = A^+y\)

这个做法有个重要的点在于正交变换具有不改变模长的性质，正交变换nb！

推荐阅读