机器学习之PCA原理

PCA(principal Component Analysis)是一种常用的数据分析方法。PCA通过线性变换将原始数据变换为一组各维度线性无关的表示,可用于提取数据的主要特征分量,常用于高维数据的降维。这篇文章旨在分析PCA的数学原理,了解PCA的工作机制是什么。

数据的向量表示及降维问题

一般情况下,在数据挖掘和机器学习中,数据被表示为向量。例如某个淘宝店2018年全年的流量及交易情况可以看成一组记录的集合,其中每一天的数据是一条记录,格式如下:

(日期,浏览量,访客数,下单数,成交数,成交金额)

PS:其实上面是一个列向量

其中『日期』是一个记录标志而非度量值,而数据挖掘关心的大多是度量值,因此如果我们忽略日期这个字段后,我们得到一组记录,每条记录可以被表示为一个五维向量,其中一条看起来应该是这个样子:

( 500 , 240 , 25 , 13 , 2312.15 ) T (500,240,25,13,2312.15)^{T} (500,240,25,13,2312.15)T

注意这里我用了转置,因为习惯上使用列向量表示一条记录(后面会看到原因),本文后面也会遵循这个准则。不过为了方便有时我会省略转置符号,但我们说到向量默认为都是列向量

我们当然可以对这一组五维向量进行分析和挖掘,不过我们知道,很多机器学习算法的复杂度和数据的维数有着密切关系,甚至与维度呈指数级关联。当然,这里区区五维的数据,也许还无所谓,但是在实际机器学习中处理成千上万甚至几十万维的情况也并不罕见,在这种情况下,计算机的资源(如内存)消耗是不可接受的,因此我们必须对数据进行降维。

但降维当然意味着信息的丢失,不过鉴于实际数据本身常常存在的相关性,我们可以想办法在降维的同时将信息的损失尽量降低。

举个例子,假如某学籍数据有两列M和F,其中M列的取值是如果此学生为男性则M=1,为女性则M=0;而F列是学生为女性时F=1,为男性时F=0.此时如果我们统计全部学籍数据,会发现对于任何一条记录来说,当M为1时F必定为0,反之当M为0时F必定为1。在这种情况下,我们将M或F去掉实际上没有任何信息的损失,因为只要保留一列就可以完全还原另一列。

当然上面是一个极端的情况,在现实中也许不会出现,不过类似得情况还是很常见的。例如上面淘宝店铺的数据,从经验我们可以知道,『浏览量』和『访客量』往往具有较强的相关性,而『下单数』和『成交数』也具有较强的相关性。这里我们非正式的使用『相关性』这个词,可以直观理解为『当某一天这个店铺的浏览量较高(或较低)时,我们应该很大程度上认为这天的访客数也较高(或较低)』。后面的章节中我们会给出相关性的严格数学定义。

这种情况表明,如果我们删除『浏览量』或『访客数』其中一个指标,我们应该不会丢失太多的信息。因此我们可以删除一个,以降低机器学习算法的复杂度。

上面给出的是降维的朴素思想描述,可以有助于直观理解降维的动机和可行性,但并不具有操作指导意义。例如,我们到底删除哪一列损失的信息才最小?亦或是根本不是单纯删除几列,而是通过某些变换将原始数据变为更少的列但又使得丢失的信息最小?到底如何度量丢失信息的多少?如何根据原始数据决定具体的降维操作步骤?

要回答上面的问题,就要对降维问题进行数学化和形式化的讨论。而PCA是一种具有严格数据基础并且已被广泛采用的降维方法。下面我不会直接描述PCA,而是通过逐步分析问题,让我们一起重新『发明』一遍PCA。

向量的表示及其变换

下面先来看一个高中就学过的向量运算:内积。两个维数相同的向量的内积被定义为:

( a 1 , a 2 , . . . , a n ) T ⋅ ( b 1 , b 2 , . . . , b n ) = a 1 b 1 + a 2 b 2 + . . . + a n b n (a_1,a_2,...,a_n)^{T}\cdot (b_1,b_2,...,b_n)=a_{1}b_{1}+a_{2}b_{2}+...+a_{n}b_{n} (a1,a2,...,an)T(b1,b2,...,bn)=a1b1+a2b2+...+anbn

内积运算将两个向量映射为一个实数。其计算方式非常容易理解,但是其意义并不明显。下面我们分析内积的几何意义。假如A和B是两个n维向量,我们知道n维向量可以等价表示为n维空间中的一条从原点发射的有向线段,为了简单起见我们假设A和B均为二维向量,则 A = ( x 1 , y 1 ) , B = ( x 2 , y 2 ) A=(x_1,y_1),B=(x_2,y_2) A=(x1,y1),B=(x2,y2)。则在二维平面上A和B可以用两条发自原点的有向线段表示,见下图:

在这里插入图片描述

现在我们从A点向B所在直线引一条垂线。我们知道垂线与B的焦点叫做A在B上的投影,再设A与B的夹角是 α \alpha α,则投影的矢量长度为 ∣ A ∣ c o s α \left | A \right |cos\alpha Acosα,其中 ∣ A ∣ = x 2 + y 2 \left | A \right |=\sqrt{x^{2}+y^{2}} A=x2+y2 ,是向量A的模,也就是A线段的标量长度。

注意这里我们专门区分了矢量长度和标量长度,标量长度总是大于等于0,其值就是线段的长度;而矢量长度可能为负,其绝对值是线段长度,而符号取决于其方向与标准方向相同或相反。

到这里还是看不出内积和这东西有什么关系,不过如果我们将内积表示为另一种我们熟悉的形式:

A ⋅ B = ∣ A ∣ ∣ B ∣ c o s α A\cdot B=\left | A \right |\left | B \right |cos\alpha AB=ABcosα

现在事情似乎是有些眉目了:A与B的内积等于A到B的投影长度乘以B的模。再进一步,如果我们假设B的模为1,即让 ∣ B ∣ = 1 \left | B \right |=1 B=1,那么就变成了:

A ⋅ B = ∣ A ∣ c o s α A\cdot B=\left | A \right |cos\alpha AB=Acosα

也就是说,设向量B的模等于1,则A与B的内积值等于A向B所在直线投影的矢量长度,这就是内积的一种几何解释,也是我们得到的第一个重要结论。再后面的推导中,将反复使用这个结论。

下面我们继续在二维空间内讨论向量。上文说过,一个二维向量可以对应二维笛卡尔直角坐标系中从原点出发的一个有向线段。例如下面这个向量:

在这里插入图片描述

在代数表示方面,我们经常用线段终点的点坐标表示向量,例如上面的向量可以表示为(3,2),这是我们再熟悉不过的向量表示。

不过我们常常忽略,只有一个(3,2)本身是不能够精确表示一个向量的。我们仔细看一下,这里的3实际表示的是向量在x轴上的投影值是3,在y轴上的投影值是2.也就是说我们其实隐式引入了一个定义:以x轴和y轴上正方向长度为1的向量为基。那么一个向量(3,2)实际上是说在x轴投影为3,而在y轴的投影为2.注意投影是一个矢量,所以可以为负。

更正式的说,向量(x,y)实际上表示线性组合:

x ( 1 , 0 ) T + y ( 0 , 1 ) T x(1,0)^T+y(0,1)^T x(1,0)T+y(0,1)T

不难证明所有二维向量都可以表示为这样的线性组合。此处(1,0)和(0,1)叫做二维空间的一组正交单位基。

在这里插入图片描述

所以,要准确描述向量,首先要确定一组基,然后给出在基所在的各个直线上的投影值,就可以了。只不过我们经常省略第一步,而默认以(1,0)和(0,1)为基。

我们之所以默认选择以(1,0)和(0,1)为基,当然是比较方便,因为它们分别是x和y轴正方向上的单位向量,因此就使得二维平面上点坐标和向量一一对应,非常方便。但实际上任何两个线性无关的二维向量都可以成为一组基,所谓线性无关在二维平面内可以直观认为是两个不在一条直线上的向量(或者说不平行)。

例如,(1,1)和(-1,1)也可以成为一组基。一般来说,我们希望基的模为1,因为从内积的意义可以看到,如果基的模为1,那么就可以方便地用向量内积来直接获得其在新基上的坐标了!实际上,对应任何一个向量我们总可以找到其同方向上模为1的向量,只要让两个分量分别除以模就好了。例如,上面的基可以变为( 1 2 , 1 2 ) \frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}}) 2 1,2 1)和( − 1 2 , 1 2 ) -\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}}) 2 1,2 1)

现在,我们想获得(3,2)在新基上的坐标,即在两个方向上的投影矢量值,那么根据内积的几何意义,我们只要分别计算(3,2)和两个基的内积,不难得到新的坐标为( 5 2 , − 1 2 \frac{5}{\sqrt{2}},-\frac{1}{\sqrt{2}} 2 5,2 1)。下图给出了新的基以及(3,2)在新基上坐标值的示意图:
在这里插入图片描述

另外这里要注意的是,我们列举的例子中基是正交的(即内积为0,或直观说相互垂直),但可以成为一组基的唯一要求就是线性无关,非正交的基也是可以的。不过因为正交基有较好的性质(欧几里得空间中的正交变换只包含旋转、反射及它们的组合,不改变向量的长度或向量之间的夹角),所以一般使用的基都是正交的。

基变换的矩阵表示

下面我们找一种简便的方式来表示基变换。还是拿上面的例子,想一下,将(3,2)变换为新基上的坐标,就是用(3,2)与第一个基坐内积运算,作为一个新的坐标分量,然后用(3,2)与第二个基做内积运算,作为第二个新坐标的分量。实际上,我们可以用矩阵相乘的形式简洁的表示这个变换:

( 1 2 1 2 − 1 2 1 2 ) ( 3 2 ) = ( 5 2 − 1 2 ) \begin{pmatrix} \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} & \\ -\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}} & \end{pmatrix}\binom{3}{2}=\binom{\frac{5}{\sqrt{2}}}{-\frac{1}{\sqrt{2}}} (2 12 12 12 1)(23)=(2 12 5)

太漂亮了!其中矩阵的两行分别为两个基,乘以原向量,其结果刚好为新基的坐标。可以稍微推广一下,如果我们有m个二维向量,只要将二维向量按列排成一个两行m列矩阵,然后用“基矩阵”乘以这个矩阵,就得到了所有这些向量在新基下的值。例如(1,1),(2,2),(3,3),想变换到刚才那组基上,则可以这样表示:

( 1 2 1 2 − 1 2 1 2 ) ( 1 2 3 1 2 3 ) = ( 2 2 4 2 6 2 0 0 0 ) \begin{pmatrix} \frac{1}{\sqrt{2}} &\frac{1}{\sqrt{2}} \\ -\frac{1}{\sqrt{2}} &\frac{1}{\sqrt{2}} \end{pmatrix}\begin{pmatrix} 1 & 2& 3\\ 1 & 2& 3 \end{pmatrix}=\begin{pmatrix} \frac{2}{\sqrt{2}} & \frac{4}{\sqrt{2}} &\frac{6}{\sqrt{2}} \\ 0& 0 & 0 \end{pmatrix} (2 12 12 12 1)(112233)=(2 202 402 60)

于是一组向量的基变换被简洁的表示为矩阵的相乘。

一般的,如果我们有M个N维向量,想将其变换为由R个N维向量表示的新空间中,那么首先将R个基按行组成矩阵A,然后将向量按列组成矩阵B,那么两矩阵的乘积AB就是变换结果,其中AB的第m列为A中第m列变换后的结果。

在这里插入图片描述

其中 p i p_i pi是一个行向量,表示第i个基, a j a_j aj是一个列向量,表示第j个原始数据记录。

特别要注意的是,这里R可以小于N,而R决定了变换后数据的维数。也就是说,我们可以将N维数据变换到更低维度的空间中去,变换后的维度取决于基的数量。因此这种矩阵相乘的表示也可以表示降维变换。

最后,上述分析同时给矩阵相乘找到了一种物理解释:两个矩阵相乘的意义是将右边矩阵中的每一列列向量变换到左边矩阵种每一行行向量为基所表示的空间中去。更抽象的说,一个矩阵可以表示一种线性变换。

上面我们讨论了选择不同的基可以对同样一组数据给出不同的表示,而且如果基的数量少于向量本身的维数,则可以达到降维的效果。但是我们还没有回答一个最最关键的问题:如何选择基才是最优的或者说,如果我们有一组N维向量,现在要将其降到K维(K小于N),那么我们应该如何选择K个基才能最大程度保留原有的信息?

要完全数学化这个问题非常繁杂,这里我们用一种非形式化的直观方法看这个问题。

为了避免过于抽象的讨论,我们仍以一个具体的例子展开。假设我们的数据由五条记录组成,将它们表示成矩阵形式:

( 1 1 2 4 2 1 3 3 4 4 ) \begin{pmatrix} 1 &1 & 2 & 4 & 2\\ 1&3 &3 & 4& 4 \end{pmatrix} (1113234424)

其中每一列为一条数据记录,而一行为一个字段。为了后续处理方便,我们需要将数据中心化(将坐标原点移动至数据的中心,使每段的均值为0),方法是将每个字段内所有字都减去字段均值(这样做的道理和好处后面会看到)。

我们看上面的数据,第一个字段均值为2,第二个字段均值为3,所以变换后:

( − 1 − 1 0 2 0 − 2 0 0 1 1 ) \begin{pmatrix} -1 &-1 &0 &2 &0 \\ -2&0 &0 &1 &1 \end{pmatrix} (1210002101)

我们可以看下五条数据在平面直角坐标系内的样子:

在这里插入图片描述

现在问题来了:如果我们必须使用一维来表示这些数据,又希望尽量保留原始的信息,我们是如何选择的呢?

通过上面对基变换的讨论我们知道,这个问题实际上是要在二维平面中选择一个方向,将所有数据都投影到这个方向所在的直线上,用投影值表示原始记录。这是一个实际的二维降到一维的问题。

那么如何选择这个方向(或者说基)才能尽量保留最多的原始信息呢?一种直观的看法是:希望投影后的投影值尽可能分散

以上图为例,可以看出如果向x轴投影,那么最左边的两个点会重叠在一起,中间的两个点也会重叠在一起,于是本身四个各不相同的二维点投影后只剩下两个不同的值了,这是一种严重的信息丢失,同理,如果向y轴投影最上面的两个点和分布在x轴上的两个点也会重叠。所以看来x和y轴都不是最好的投影选择。我们直观目测,如果向通过第一象限和第三象限的斜线投影,则五个点在投影后还是可以区分的

下满我们用数学的方法表述这个问题最大方差法

方差

上文说到,我们希望投影后的投影值尽可能分散,而这种分散程度,可以用数学上的方差表述。此处,一个字段的方差可以看做是每个元素与字段均值的差的平方和的均值,即:

V a r ( a ) = 1 m ∑ i = 1 m ( a i − μ ) 2 Var(a)=\frac{1}{m}\sum _{i=1}^{m}(a_i-\mu )^{2} Var(a)=m1i=1m(aiμ)2

由于上面我们已经将每个字段的均值都化为0了,因此方差可以直接用每个元素的平方和除以元素个数表示:

V a r ( a ) = 1 m ∑ i = 1 m ( a i ) 2 Var(a)=\frac{1}{m}\sum _{i=1}^{m}(a_i)^{2} Var(a)=m1i=1m(ai)2

于是上面的问题被形式化表述为:寻找一个一维基,使得所有数据变换为这个基上的坐标表示后,方差值最大

协方差

对于上面二维降成一维的问题来说,找到那个使得方差最大的方向就可以了。不过对于更高维,还有一个问题需要解决。考虑三维降到二维问题。与之前相同,首先我们希望找到一个方向使得投影后方差最大,这样就完成了第一个方向的选择,继而我们选择第二个投影方向。

如果我们还是单纯只选择方差最大的方向,很明显,这个方向与第一个方向应该是『几乎重合在一起』,显然这样的维度上没有用的,因此,应该有其他约束条件。从直观上说,让两个字段尽可能表示更多的原始信息,我们是不希望它们之间存在(线性)相关性的,因为相关性意味着两个字段不是完全独立,必然存在重复表示的信息。

数学上可以用两个字段的协方差表示其相关性,由于已经让每个字段均值为0,则:

C o v ( a , b ) = 1 m ∑ i = 1 m a i b i Cov(a,b)=\frac{1}{m}\sum _{i=1}^{m}a_ib_i Cov(a,b)=m1i=1maibi

可以看到,在字段均值为0的情况下,两个字段的协方差简洁的表示为其内积除以元素数m

当协方差为0时,表示两个字段完全独立。为了让协方差为0,我们选择第二个基时只能在与第一个基正交的方向上选择。因此最终选择的两个方向一定是正交的。

至此,我们得到了降维问题的优化目标:将一组N维向量降为K维(K大于0,小于N),其目标是选择K个单位(模为1)正交基,使得原始数据变换到这组基上后,不同维度(不同数据特征)两两间协方差为0,而相同维度的协方差(维度的方差)则尽可能大(在正交的约束下,取最大的K个方差)

协方差矩阵

上面我们导出了优化目标,但是这个目标似乎不能直接作为操作指南(或者说算法),因为它只说要什么,但根本没有说怎么做。所以我们要继续在数学上研究计算方案。

我们看到,最终要达到的目的与字段内方差及字段间协方差有密切关系。因此我们希望能将两者统一表示,仔细观察发现,两者均可以表示为内积的形式,而内积又与矩阵相乘密切相关。于是我们来了灵感:

假设我们只有a和b两个字段,那么我们将它们按行组成矩阵X:

在这里插入图片描述

然后我们用X乘以X的转置,并乘上系数1/m:

在这里插入图片描述

奇迹出现了!这个矩阵对角线上的两个元素分别是两个字段的方差,而其它元素是a和b的协方差。两者被统一到了一个矩阵的。

根据矩阵相乘的运算法则,这个结论很容易被推广到一般情况

设我们有m×n数据记录,将其按列排成n×m的矩阵X,令 C = 1 m X X T C=\frac{1}{m}XX^{T} C=m1XXT,则C是一个对称矩阵,其对角线分别是各个维度的方差,而第i行j列和j行i列元素相同,表示i和j两个维度(特征)的协方差

协方差矩阵对角化

根据上述推导,我们发现要达到优化目标,等价于将协方差矩阵对角化::即除对角线外的其它元素化为0,并且在对角线上将元素按大小从上到下排列,这样我们就达到了优化目的。这样说可能还不睡很明晰,我们进一步看下原矩阵与基变换后矩阵协方差矩阵的关系。

设原始数据矩阵X对应的协方差矩阵为C,而P是一组基按行组成的矩阵,设Y=PX,则Y为X使用P为基的变换后的数据。设Y的协方差矩阵为D,我们推导一下C与D的关系:

在这里插入图片描述

现在事情很明白了!我们要找的P不是别的,而是能让原始协方差矩阵对角化的P。换句话说,优化目标变成了寻找一个矩阵P,满足 P C P T PCP^T PCPT是对角阵,并且对角元素按从小到大依次排列,那么P的前K行就是要寻找的基,用P的前K行组成的矩阵乘以X就使得X从N维降到了K维并满足上述优化条件。

至此,我们离“发明”PCA仅一步之遥了!

现在所有焦点都聚焦在了协方差矩阵对角化问题上,由上文可知,协方差矩阵C是一个对称矩阵,在线性代数上,是对称矩阵由一系列非常好的性质:

  1. 实对称矩阵不同特征值所对应的特征向量正交。
  2. λ \lambda λ是实对称矩阵的r重根,则对应特征值 λ \lambda λ恰有r个线性无关的特征向量
  3. 设A为n阶对称矩阵,则必有正交矩阵p,使得 P − 1 A P = = P T A P = Λ P^{-1}AP==P^{T}AP=\Lambda P1AP==PTAP=Λ,其中 Λ \Lambda Λ是以A的n个特征值为对角元素的对角矩阵。

由上面两条可知,一个n行n列的实对称矩阵一定可以找到n个单位正交特征向量,设这n个特征向量为 e 1 , e 2 , . . . , e n e_1,e_2,...,e_n e1,e2,...,en,我们将其按列组成矩阵:
在这里插入图片描述

则对协方差矩阵C有如下结论:

在这里插入图片描述

其中 Λ \Lambda Λ为对角矩阵,其对角元素为各个特征向量对应的特征值。

对比上面公式,我们已经找到了需要的变换矩阵P:

P = E T P=E^{T} P=ET

P是协方差矩阵的特征向量单位化后按行排列出的矩阵,其中每一行都是C的一个特征向量。如果设P按照 Λ \Lambda Λ中特征值的从大到小,将特征向量按照相应排列下来,则用P的前K行组成的矩阵乘以原始数据矩阵X,就得到了我们需要降维后的数据矩阵Y。

我们会根据情况选择前K个比较大的特征值,因为前K个大的特征值所包含的信息一般占的比例高达99%。

至此我们完成了整个PCA的数学原理讨论。在下面的一阶,我们将给出PCA的一个实例。


事实上,通过计算数据矩阵的协方差矩阵,然后得到协方差矩阵的特征值特征向量,选择特征值最大(即方差最大)的k个特征所对应的特征向量组成的矩阵。这样就可以将数据矩阵转换到新的空间中,实现数据特征的降维。

由于得到协方差矩阵的特征值特征向量有两种方法:特征值分解协方差矩阵、奇异值分解协方差矩阵,所以PCA算法有两种实现方法:基于特征值分解协方差矩阵实现PCA算法、基于SVD分解协方差矩阵实现PCA算法


基于特征值分解协方差矩阵实现PCA算法

设有m条n维数据:

  1. 将原始数据按列组成n行m列矩阵X
  2. 将X的每一行(代表一个特征)进行零均值化(中心化),即减去这一行的均值
  3. 求出协方差矩阵 C = 1 m X X T C=\frac{1}{m}XX^T C=m1XXT
  4. 求出协方差矩阵的特征值及对应的特征向量
  5. 将特征向量按对应特征值大小从上到下按行排列称矩阵,取前k行组成矩阵P
  6. Y=PX即为降维后k维后的数据

实例

在这里插入图片描述

在这里插入图片描述

降维后的效果图为:

在这里插入图片描述


基于SVD分解协方差矩阵实现PCA算法

特征值分解矩阵的缺点

我们前面讲了很多特征值、特征向量和特征值分解,而且基于我们以前学习的线性代数知识,利用特征值分解提取特征矩阵是一个容易理解且便于实现的方法。但是为什么还存在SCD(奇异值分解)呢?特征值分解最大的问题上只能针对方阵,即n*n的矩阵。而在实际的应用中,我们分解的大部分都不是方阵。

举个例子

关系型数据库中的某一张表的数据存储结构就类似于一个二维矩阵,假设这个表有m行,有n个字段,那么这个表数据矩阵规模就是m*n。很明显,在绝大部分情况下,m与n是不相等的。如果这个时候要对这个矩阵进行特征提取,特征值分解的方法明显就不行了。此时,就可以用SVD对非方阵矩阵进行分解。


奇异值分解

奇异值分解是一个能适用于任意矩阵的一种分解方法,对于任意矩阵A总是存在一个奇异值分解:

A = U Σ V T A=U\Sigma V^T A=UΣVT

假设A是一个m×n的矩阵,那么得到的U是一个m×m的方阵,U里面的正交向量被称为左奇艺向量 Σ \Sigma Σ除了对角线其他原始都为0,对角线上的元素称为奇异值 V T V^T VT是V的转置矩阵,是一个n×n的矩阵,它里面的正交向量被称为右奇异矩阵。而且一般来讲,我们会将 Σ \Sigma Σ上的值按从大到小的顺序排列。上面矩阵的维度裱花可以参考下图4:

在这里插入图片描述

那么接下来的问题是:如何求得左奇异向量、右奇异向量和奇异值呢?

首先,我们用矩阵A的转置乘以A,得到一个方阵,用这样的方阵进行特征分解,得到的特征值和特征向量满足下面的等式:

( A T A ) v i = λ i v i (A^TA)v_i=\lambda_{i}v_i (ATA)vi=λivi

这里的 v i v_i vi就是我们要求的右奇异向量。

其次,我们将A和A的转置做矩阵的乘法,得到一个方阵,用这样的方阵进行特征分解,得到的特征和特征向量满足下满的等式:

( A A T ) u i = λ i u i (AA^T)u_i=\lambda_{i}u_i (AAT)ui=λiui

这里的 u i u_i ui就是左奇异向量。

下面我们来证明一下为什么 A T A A^TA ATA的特征向量组成的矩阵就是我们的SVD中的V矩阵,而 A A T AA^T AAT的特征向量组成的就是我们的SVD中的U矩阵:

A = U Σ V T ⇒ A T = V Σ T U T ⇒ A T A = V Σ T U T U Σ V T = V Σ 2 V T A=U\Sigma V^T\Rightarrow A^T=V\Sigma^TU^T\Rightarrow A^TA=V\Sigma^TU^TU\Sigma V^T=V\Sigma ^{2}V^T A=UΣVTAT=VΣTUTATA=VΣTUTUΣVT=VΣ2VT………………(*)

得:

V T A T A V = Σ 2 V^TA^TAV=\Sigma ^{2} VTATAV=Σ2

说明将 A T A A^TA ATA对角化。

上式证明使用了 U T U = I , Σ T Σ = Σ 2 U^TU=I,\Sigma^T\Sigma=\Sigma ^{2} UTU=I,ΣTΣ=Σ2.可是看出 A T A A^TA ATA的特征向量组成的矩阵就是我们SVD中的V矩阵,而 A A T AA^T AAT的特征向量组成了我们SVD中的U矩阵。

补充定义
若U是n阶方阵,满足 U T U = I U^TU=I UTU=I,则U是实正交矩阵。

此外,我们还可以得到奇异值,奇异值的求法有两种:

第一种求法

A = U Σ V T ⇒ A V = U Σ V T V ⇒ A V = U Σ ⇒ A v i = σ i μ i ⇒ σ i = A v i μ i A=U\Sigma V^T\Rightarrow AV=U\Sigma V^TV\Rightarrow AV=U\Sigma\Rightarrow Av_i=\sigma_i \mu_i\Rightarrow \sigma_i=\frac{Av_i}{\mu_i} A=UΣVTAV=UΣVTVAV=UΣAvi=σiμiσi=μiAvi

第二种求法

通过上面的(*)式的证明,我们还可以看出,特征值矩阵等于奇异值矩阵的平方,也就是特征值与奇异值满足如下关系:

σ i = λ i \sigma_i=\sqrt{\lambda _i} σi=λi

这里的 σ i \sigma_i σi就是奇异值,奇异值 σ i \sigma_i σi跟特征值类似,在矩阵 Σ \Sigma Σ也是从大到小排列。


我们已经知道如何用奇异值分解任何矩阵了,那么问题又来了,一个m×n的矩阵A,你把它分解成m×m的矩阵U、m×n的矩阵Σ和n×n的矩阵V的转置。这三个矩阵中任何一个的维度似乎一点也不比A的维度小,而且还要做两次矩阵的乘法,这不是没事找事干嘛!把简单的事情搞复杂了么!并且我们知道矩阵乘法的时间复杂度为 O ( n 3 ) O(n^3) O(n3)

其实,在奇异值分解矩阵中 Σ \Sigma Σ里面的奇异值按从大到小的顺序排列,奇异值 σ i \sigma_i σi从大到小的顺序减小的特别快。在很多情况下,前10%甚至1%的奇异值的和就占了全部的奇异值之和的99%以上。也就是说,剩下的90%甚至99%的奇异值几乎没有什么作用。因此,我们可以用前面r个大的奇异值来近似描述矩阵,于是奇异值分解公式可以写成如下:

A m × n ≈ U m × r Σ r × r V r × n T A_{m\times n}\approx U_{m\times r}\Sigma_{r\times r}V_{r\times n}^{T} Am×nUm×rΣr×rVr×nT

其中r是一个远远小于m和n的数,右边的三个矩阵相乘的结果将会使一个接近A的矩阵。如果r越接近于n,则相乘的结果越接近于A。如果r的取值远远小于n,从计算机内存的角度来说,右边三个矩阵的存储内存要远远小于矩阵A的。所以在奇异值分解中r的取值很重要,就是在计算精度和时间空间之间做选择


SVD例子

在这里插入图片描述
在这里插入图片描述


基于SVD分解协方差矩阵实现PCA算法

输入:数据集 X = { x 1 , x 2 , x 3 , . . . , x n } X=\left \{ x_1,x_2,x_3,...,x_n \right \} X={x1,x2,x3,...,xn},需要降到k维。

  1. 中心化,即每一位特征减去各自的平均值。
  2. 计算协方差矩阵。
  3. 通过SVD计算协方差矩阵的特征值与特征向量。
  4. 对特征值从大到小排序,选择其中最大的k个。然后将其对应的k个特征向量分别作为列向量组成特征向量矩阵。
  5. 将数据转到k个特征向量构建的新空间中。

在PCA降维中,我们需要找到样本协方差矩阵 X X T XX^T XXT的最大k个特征向量,然后用这最大的k个特征向量组成的矩阵来做低维投影降维。可以看出,在这个过程中需要先求出协方差矩阵 X X T XX^T XXT,当样本数多、样本特征数也多的时候,这个计算还是很大的,当我们用到SVD分解协方差矩阵的时候,SVD有两个好处:

  1. 有一些SVD的实现算法可以先不求出协方差矩阵 X X T XX^T XXT也能求出我们的右奇异矩阵V。也就是说,我们的PCA算法可以不用做特征分解而是通过SVD来完成,这个方法在样本量很大的时候很有效。实际上,scikit-learn的PCA算法的背后真正的实现就是用的SVD,而不是特征值分解.
  2. 假设我们的样本是m×n的矩阵X,如果我们通过SVD找到了矩阵 X X T XX^T XXT最大的k个特征向量组成的k×n的矩阵 V T V^T VT ,则我们可以做如下处理:

X m × k ′ = X m × n V k × n T {X}'_{m\times k}=X_{m\times n}V_{k\times n}^T Xm×k=Xm×nVk×nT

可以得到一个m×k的矩阵 X ′ {X}' X,这个矩阵和我们原来的m×n矩阵X相比,列数从n减到了k,可知对列数进行了压缩。

左奇异矩阵可以用于对行数(即样本个数)进行压缩;右奇异矩阵可以用于对列数(即特征维度)进行压缩。

这就是我们用SVD分解协方差矩阵实现PCA可以得到两个方向的PCA降维(即行和列两个方向)。

  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值