从矩阵(matrix)角度讨论PCA(Principal Component Analysis 主成分分析)、SVD(Singular Value Decomposition 奇异值分解)相关原理

0. 引言

本文主要的目的在于讨论PAC降维和SVD特征提取原理,围绕这一主题,在文章的开头从涉及的相关矩阵原理切入,逐步深入讨论,希望能够学习这一领域问题的读者朋友有帮助。

这里推荐Mit的Gilbert Strang教授的线性代数课程,讲的非常好,循循善诱,深入浅出。

Relevant Link:  

Gilbert Strang教授的MIT公开课:数据分析、信号处理和机器学习中的矩阵方法
https://mp.weixin.qq.com/s/gi0RppHB4UFo4Vh2Neonfw

 

1. 可逆矩阵

0x1:可逆矩阵的基本概念

对于数域K上的矩阵A,如果存在数域K上的矩阵B,使得:

则称A是可逆矩阵(或非奇异矩阵)。

如果A是可逆矩阵,则适合上式的矩阵B称为A的逆矩阵,记作

0x2:可逆矩阵的基本性质

1. 可逆矩阵是互为可逆的(可交换)

如果A是可逆矩阵,则它又逆矩阵使得:

从上式可看出,也是可逆矩阵,并且:

2. 可逆矩阵一定是方阵

从矩阵可逆的公式可以看出,A与B可交换,因此可逆矩阵一定都是方阵。

3. 矩阵的逆矩阵如果存在,一定是唯一的

利用反证法,假如还有矩阵B1也适合上式,则:

因此,B1 = B。

4. 如果n级矩阵A,B都可逆,则AB也可逆

并且:

5. 矩阵的逆和矩阵的转置满足交换律

如果A可逆,则A'也可逆,并且:

6. 可逆矩阵经过初等行变换化成的简化行阶梯矩阵一定是单位矩阵 I

7. 用一个可逆矩阵去左(右)乘矩阵A,不改变A的秩rank()

0x3:可逆矩阵的充分必要条件

1. 行列式判别法

数据K上n级矩阵A可逆的充分必要条件是: 。

值得注意的是,也是矩阵对应线性方程组有解的充要条件。

2. 矩阵秩判别法

n级矩阵A可逆的充分必要条件是 rank(A) = n,即A为满秩矩阵。

3. 向量组线性相关性判别法

数据K上n级矩阵A可逆的充分必要条件是A的行(列)向量组线性无关。 

4. 矩阵A可逆的充分必要条件是它可以表示成一些初等矩阵的乘积

充分性:设A可以表示成一些初等矩阵的乘积,由于初等矩阵都可逆,因为他们的乘积A也可逆。

必要性:设A可逆,则A经过初等行变换化成的简化行阶梯矩阵一定是单位矩阵I,因此有初等矩阵,使得

因此

由于初等矩阵的逆矩阵仍是初等矩阵,因此上式表明:矩阵A可以表示为一些初等矩阵的乘积。

0x4:逆矩阵求解方法

1. 伴随矩阵法

当矩阵A可逆时,

上式给出了求逆矩阵的一种方法,称之为伴随矩阵法

2. 初等变换化简法

有的时候,当矩阵的阶数比较高的时候,使用其行列式的值和伴随矩阵求解其逆矩阵会产生较大的计算量。这时,我们可以采用初等变换法进行逆矩阵求解。

设A是n级可逆矩阵,根据可逆矩阵性质,我们知道,有初等矩阵,使得:

所以又有:

从上式可以看出,如果有一系列初等行变换把A化成了单位矩阵I,那么同样的这些初等行变换就把I化成了

因此我们可以把A与I并排放在一起,组成一个nx2n级矩阵(A,I),对(A,I)作一系列初等行变换,把它的左半部分化成I,这时的右半部分就是,即:

这种求逆矩阵的方法称为初等变换法,这时最常用的方法。

0x5:线性回归模型参数求解的一种方式 - 逆矩阵求解法

在线性回归模型中,根据最小二乘求解公式,我们有:

在满秩情况下(满足数据点个数>特征数量),系数求解等价于存在一个简单的闭式矩阵解,使得最小二乘最小化。由下式给出:

其中,是矩阵的逆矩阵,令,则最终学得的多元线性回归模型为:

可以看到,参数求解的本质就是矩阵乘法,矩阵求逆等操作。

当然,需要指出的是,对于linear regression model来说,同样可以使用GD梯度下降的方式进行参数求解,感兴趣的读者朋友可以下载sklearn的内核源代码进行一对一的分析,可以获得更直观的体验。

Relevant Link:  

https://www.cnblogs.com/LittleHann/p/10498579.html#_label3_1_2_1 
https://wenku.baidu.com/view/1963472d9a6648d7c1c708a1284ac850ad0204ad.html

 

2. 正交矩阵

0x1:正交矩阵的基本概念

在平面上取一个直角坐标系Oxy,设向量a,β的坐标分别是(a1,a2),(b1,b2)。如果a,β都是单位向量,并且互相垂直,则他们的坐标满足:

上述4个等式可以写成一个矩阵形式:

设矩阵A的第1,2行分别是a,β的坐标,则上面矩阵形式可以重写为:

我们把满足这种特性的矩阵A称为正交矩阵

一般的,如果是实数域上的方阵A满足:,则称A是正交矩阵

0x2:正交矩阵的基本性质

1. 单位矩阵 I 是正交矩阵

2. 若A与B都是n级正交矩阵,则AB也是正交矩阵

若A,B都是正交矩阵,则有下式:

因此AB也是正交矩阵。

3. 若A是正交矩阵,则也是正交矩阵

4. 若A是正交矩阵,则 |A| =1 或 -1 

若A是正交矩阵,则,从而,因为转置不概念行列式,即

所以,|A| = 1,或者 -1。

5. 实数域上n级矩阵A是正交矩阵的充要条件为:A的行(列)向量组是欧几里得空间Rn的一个标准正交基

该定理告诉我们,构建正交矩阵等价于求标准正交基,许多实际问题需要构造正交矩阵,于是我们要设法求标准正交基。

0x3:正交矩阵由一组正交单位行(列)向量组成

设实数域上n级矩阵A的行向量组为,列向量组为,则:

(1)A为正交矩阵当且仅当A的行向量满足下式:

(2)A为正交矩阵当且仅当A的列向量组满足下式:

我们引用一个符号,它的含义是:

称为Kronecker(克罗内克)记号,采用这个符号,上面两式可以分别写成:

从上式可以看到,从行向量和列向量的角度,正交矩阵都满足上面的定理。

这两组式子的左端都是两个n元有序数组的对应元素的乘积之和。与几何空间中两个向量的内积在直角坐标系中的计算公式相似。因此,我们可以在实数域上的n维向量空间也引入内积的概念。

0x4:向量组标准内积

1. 标准内积公式

中,任给两个向量组,规定

这个二元实值函数(a,β)称为的一个内积,通常称这个内积为的标准内积

上式也可以简写为:

2. 标准内积基本性质

对一切,有下列几个性质

1)对称性:(a,β)= (β,a)

2)线性性之一:(a+γ,β)=(a,β)+(γ,β)

3)线性性之二:(ka,β)= k(a,β) 

4)正定性:(a,a)>=0,等号成立当且仅当 a=0

0x5:欧几里得空间 

n维向量空间有了标准内积后,就称为一个欧几里得空间

在欧几里得空间中,向量a的长度 |a| 规定为:

长度为1的向量称为单位向量,a是单位向量的充要条件是

因为,于是对于一定是单位向量。

把非零向量 a 乘以,称为把a单位化。

0x6:正交向量组

1. 正交向量组基本定义

在欧几里得空间中,如果,则称a与β是正交的,记作

在欧几里得空间中,由非零向量组成的向量组如果其中每两个不同的向量都彼此正交,则称它们是正交向量组

特别的,仅由一个非零向量组成的向量组也是正交向量组。

同时,如果正交向量组的每个向量都是单位向量,则称它为正交单位向量组

2. 正交向量组基本性质

欧几里得空间中,正交向量组一定是线性无关的。

证明:

是正交向量组,设

把上式两端的向量都与作内积,得:

由于,当 j != i 时,因此上式得:

由于,因此,从而得:

根据线性相关基本定理,线性无关。

3. 正交基

欧几里得空间中,n个向量组成的正交向量组一定是的一个基,称它为正交基。n个单位向量组成的正交向量组称为的一个标准正交基。 

Relevant Link:  

https://baike.baidu.com/item/%E6%96%BD%E5%AF%86%E7%89%B9%E6%AD%A3%E4%BA%A4%E5%8C%96/756386?fr=aladdin
《建明线性代数》丘维声

 

3. 矩阵的相抵与相似

0x1:矩阵的相抵

1. 矩阵相抵的基本定义

如果矩阵A可以经过一系列初等行变换与初等列变换变成矩阵B,则称A与B是相抵的(或等价),记作

由于矩阵的初等行(列)变换可以通过初等矩阵与矩阵的乘法来实现,并且一个矩阵可逆的充要条件是它能表示成一些初等矩阵的乘积,因此:

s x n 矩阵A与B相抵,等价于:

→ A经过初等行变换和初等列变换变成B;

→ 存在s级初等矩阵与n级初等矩阵,使得:

→ 存在s级可逆矩阵P与n级可逆矩阵Q,使得:

2. 矩阵相抵的基本性质

s x n 矩阵之间的相抵关系具有下述性质

1)反身性:任一矩阵A与自身相抵;

2)对称性:如果A与B相抵,则B与A相抵;

3)传递性:如果A与B相抵,B与C相抵,则A与C相抵;

0x2:矩阵的相似

1. 矩阵相似的基本定义

设A与B都是数域K上n级矩阵(注意,必须是n级方阵才存在相似矩阵),如果存在数域K上的一个n级可逆矩阵U,使得

则称A与B是相似的,或说B是A的相似矩阵,记作

2. 矩阵相似的基本性质

数域K上n级矩阵之间的相似关系具有下列性质:

1)反身性:任一 n 级矩阵A与自身相似;

2)对称性:如果,则

3)传递性:如果,则

4)相似的矩阵其行列式的值相同

5)相似的矩阵或者都可逆,或者都不可逆

并且当它们可逆时,它们的逆矩阵也相似

6)相似的矩阵有相同的秩

7)相似的矩阵又相同的迹

n级矩阵A的主对角线上元素的和称为A的迹,记作tr(A)

3. 矩阵的相似数学本质 

先抛出结论:同一个线性变换,不同基下的矩阵,称为相似矩阵

这个结论不容易简单阐述清楚,我们先从一些简单的概念开始,逐步切入概念。上面说到不同基下的矩阵,那什么是不同向量基呢?即什么是向量基变换呢?

1)什么是向量基变换

我们通过一个简单的例子来直观地说明什么是坐标(向量基)变换。

坐标转换是数学中的常用手段,目的是简化运算,提高可视化理解程度。比如常见的,把直角坐标系(xy坐标系)的圆方程,换元为极坐标(\rho\theta坐标系):

下图展示了换元前后的函数形式:

从直观上看,换元后的代数式和图像都变简单了。从线性方程角度,换元将原方程从非线性空间转换到了线性空间。

2)矩阵代表了一种基下的一种线性变换

线性函数y=x可以认为是把(a,0)点映射到(0,a)点,我们称为线性变换T,记作:

T:(a,0)\to(0,a),a\in\mathbb{R},b\in\mathbb{R}

该线性变换矩阵的形式如下:

进一步推广开来,只要替换(a,0)为平面内所有的点(a,b),我们就可以对整个平面做变换,该线性变换记作:

进而可以写作矩阵的形式:

我们记:

\vec{y_{}}=\begin{pmatrix}b\\a\end{pmatrix}\qquad\vec{x_{}}=\begin{pmatrix}a\\b\end{pmatrix}\qquad A=\begin{pmatrix}0&1\\1&0\end{pmatrix}

我们可以得到更简便的记法:

\vec{y_{}}=A\vec{x_{}}

下图用淡蓝色网格来表示这个线性变换(这个变换实际上镜面反转):

上面的讨论中,可能读者朋友觉得非常自然,这是因为我们在高中、大学的学习中,对直角坐标系已经非常熟悉了,头脑中也已经建立了非常直观感性的空间想象能力。

其实在前面的讨论中隐含了一个特征重要的信息,就是坐标系(向量基)

y=x是基于直角坐标系的(标准正交基)

标准正交基是\{\vec{i_{}}=\begin{pmatrix}1\\0\end{pmatrix},\vec{j_{}}=\begin{pmatrix}0\\1\end{pmatrix}\},它们所张成的线性空间如下:

通过这个转换:

y=x\to A=\begin{pmatrix}0&1\\1&0\end{pmatrix}

得到的A也是基于直角坐标系的

3)相似矩阵 - 不同基下的相同线性变换 

我们先来说明不同基下的矩阵的变换思路:

上图中:

  • V_1:\{\vec{i_{}},\vec{j_{}}\}和 V_2:\{\vec{i'},\vec{j'}\}:代表两个不同的向量基;
  • V1\to V2:代表了向量基V1通过某种线性变换,转换为向量基V1,可以通过P^{-1}转换;
  • V2\to V1:同上,可以通过P转换;
  • \vec{v'}V2下的点;
  • \vec{v'}通过P变为V1下的点,即P\vec{v'}
  • V1下,通过A矩阵完成线性变换,即AP\vec{v'}
  • 通过P^{-1}从变回V2下的点,即P^{-1}AP\vec{v'}

综上,我们可以有:

B\vec{v'}=P^{-1}AP\vec{v'}

即矩阵A和矩阵B是可以互相转化的,它们是同一个线性变换在不同向量基下的不同表现形式。

推广到一般情况,我们可以认为,:

B=P^{-1}AP

那么BA互为相似矩阵。

4)得到相似矩阵经历的相似变换本质上是对向量基的变换

继续来讨论相似矩阵公式: 

B=P^{-1}AP

前面说到,P代表了向量基转换V2\to V1,上个小节一笔带过了,这个小节我们来详细讨论下这个转换是如何完成的。

首先给出我们空间中的一点,例如m点:

不论有没有向量基,这个点都是客观存在的。

然后,我们定义V2的向量基,\vec{i'},\vec{j'},之后,我们给出m点在\vec{i'},\vec{j'}的坐标\vec{v'}

重写\vec{v'}在该向量基下的表示:

\vec{v'}=\begin{pmatrix}a\\b\end{pmatrix}=a\vec{i'}+b\vec{j'}

现在我们来看另一个向量基V1,假设我们知道了\vec{i'},\vec{j'}在V1,\vec{i},\vec{j}下的坐标(如下图所示): 

那么将V1下的坐标带入原式:


\begin{aligned}
    \vec{v'} & = a\vec{i'}+b\vec{j'}\\
             & = a(c\vec{i}+d\vec{j})+b(e\vec{i}+f\vec{j})
\end{aligned}

此时,实际上m点的坐标,已经从V2变到了V1向量基\vec{i},\vec{j}下的\vec{v}

转换为矩阵形式:

所以P其实就是:

P=\begin{pmatrix}\vec{i'}&\vec{j'}\end{pmatrix},式中,上\vec{i'},\vec{j'}是在\vec{i},\vec{j}下的坐标(向量基)。

由此可以看到,P实际上代表了一种对向量基的转换,或者说是向量基的映射。 

Relevant Link:  

https://www.matongxue.com/madocs/491.html

 

4. 矩阵特征值、特征向量

0x1:特征值、特征向量的定义

设A是数域K上的n级矩阵,如果中有非零列向量a,使得,则称是A的一个特征值,称a是A的属于特征值的一个特征向量

例如,设

由于

因此,a是A的一个特征值,a是A的属于特征值2的一个特征向量。

0x2:特征值、特征向量基本性质

如果a是A的属于特征值的一个特征向量,则

从而对于任意的,有

0x3:如何判断矩阵A是否存在特征值以及求全部特征值

1. 特征值存在性判断定理

首先,不是所有的矩阵都存在特征值,我们来看下面这个例子。

设σ是平面上绕原点O的转角为π/3的旋转,旋转是一种线性变换,可以用矩阵来表示线性变换,则σ可以用下述矩阵A来表示:

由于平面上任一非零向量在旋转σ下都不会变成它的倍数,因此在中不存在非零列向量a满足。从而A没有特征值,没有特征向量。

接下来讨论如何判断数域K上的n级矩阵A是否有特征值和特征向量,如果有的话,怎样求A的全部特征值和特征向量。

设:

是数域K上的矩阵,则是A的一个特征值,a是A的属于的一个特征向量,等价于下列几个公式:

→ 

→ 

→ a 是齐次线性方程组的一个非零解,

→ ,a是的一个非零解,

由于:

我们把多项式称为矩阵A的特征多项式

于是从上面推导得到,是A的一个特征值,a是A的属于的一个特征向量,等价于:

→ 是A的特征多项式在K中的一个根,a是齐次线性方程组的一个非零解。

上述推导过程可以推广为在任意n级矩阵下适用的定理:

(1)是A的一个特征值,当且仅当是A的特征多项式在K中的一个根;

(2)a是A的属于特征值的一个特征向量,当且仅当a是齐次线性方程组的一个非零解

2. 特征值求解方法

n级矩阵的特征多项式写出来就是:

于是利用上式可判断数域K上n级矩阵A有没有特征值和特征向量,如果有的话,求A的全部特征值和特征向量的方法如下:

(1)第一步:计算A的特征多项式

(2)第二步:判别多项式在数域K中有没有根,如果它在K中没有根,则A没有特征值,从而A也没有特征向量。反之,如果在K中有根,则它在K中的全部根就是A的全部特征值,此时接着做第三步;

(3)第三步:对于A的每一个特征值,求齐次线性方程组的一个基础解系。于是A的属于的全部特征向量组成的集合是

是A的一个特征值,我们把齐次线性方程组的解空间称为A的属于特征子空间。它的全部非零向量就是A的属于的全部特征向量。

0x4:特征值和特征向量的物理意义 

如果把矩阵看作是运动

  • 特征值就是运动的速度
  • 特征向量就是运动的方向

特征值、特征向量自然可以称为运动(即矩阵)的特征。 

0x5:特征值和特征向量的几何意义

在一个特定向量基\vec{i_{}},\vec{j_{}}下面有一个向量\vec{v_{}},我们先随便左乘一个矩阵A,下图所示:

我们知道,矩阵乘法的意义是在一个特定的向量基下对向量进行一个线性变换,上图中的变换结果看起来很普通没啥特殊的。

我们接下来调整下\vec{v_{}}的方向:

可以观察到,调整后的\vec{v_{}}A\vec{v_{}}在同一根直线上,只是A\vec{v_{}}的长度相对\vec{v_{}}的长度变长了。
此时,我们就称\vec{v_{}}A的特征向量,而A\vec{v_{}}的长度是\vec{v_{}}的长度的\lambda倍,\lambda就是特征值。
从而,特征值与特征向量的几何意义如下图所示:

实际上,这个矩阵A的特征值有2个,对应的特征向量也有2个,笔者用gif动图展示了对应的特征向量。 

 

笔者这里打开了Av的迹追踪(图中的红点轨迹),可以看到,矩阵A对应的特征向量的方向正好就是A投影离散度最大的方向,同时两个特征向量彼此垂直,关于这点我们后面讨论PAC降维原理的时候会再谈到。

原始的geogebra链接在这里,读者朋友可以自己上手进行修改和可视化观察。

0x6:矩阵特征值分解

对于矩阵A如果可以对角化的话(特征值/特征向量存在),可以通过相似矩阵进行下面这样的特征值分解:

A=P\Lambda P^{-1},其中\Lambda为对角阵,P的列向量是单位化的特征向量。

下图我们展示了一个矩阵的相似对角矩阵分解过程:

可以看到,矩阵的对角化分解将原矩阵分解为了特征向量特征值两个部分,这个公式笔者认为是线性代数里最美妙的公式之一了。

首先,我们可以看到,矩阵的本质是将原始向量带入一个新的向量基视角下(可能升维也可能降维),之后进行两种形式的运动:
  • 旋转
  • 拉伸

整体矩阵代表的最后的运动结果就是这两种的合成。

读者思考:原始矩阵A和与之相似的对角矩阵(由特征值组成),本质上是同一个线性变换在不同向量基下的不同表现形式

Relevant Link:   

http://mini.eastday.com/bdmip/180328092726628.html#
https://blog.csdn.net/special00/article/details/84033124
http://mini.eastday.com/bdmip/180328092726628.html# 
https://www.geogebra.org/m/F5dBUqAt
https://www.matongxue.com/madocs/228/

 

5. 矩阵奇异值分解

我们上一个章节讨论了矩阵的特征值分解,需要特别注意的是,特征值分解只适用于方阵(矩阵相似也只适用于方阵)。而在现实的世界中,我们看到的大部分矩阵都不是方阵,虽然在图像处理领域可以通过resize强行获得一个方阵matrix,但是这意味着还没开始处理就先失真了一部分了。

相比于矩阵特征值分解,矩阵奇异值分解,就是一种通用性更强的矩阵特征分解方法,奇异值分解可以适用于任意的矩阵,它本质上也矩阵特征值分解的目的是一样的,都是分解出矩阵的向量方向分量不同方向的分量强度

0x1:用翻绳游戏类比矩阵奇异值分解

我们知道,对于翻绳的这个游戏而言,每轮的花型是由四只手完成的:

我们可以认为这个花型是由两个方向的力合成的:

我们可以想象,如果其中一个力(相比另外一个力而言)比较小的话,那么绳子的形状基本上由大的那个力来决定:

或者换句话说,每个花型从力学上都可以分解为两个不同手方向上的力。不同的花型,分解的结果也不同。

将奇异值分解类比于翻绳,我们可以认为:
  • 奇异值分解,就是把矩阵分成多个“力分量”
  • 奇异值的大小,就是各个“力分量”的实值大小

0x2:奇异值分解数学公式

了解了奇异值分解的基本感性概念,接下来需要先定义一下矩阵奇异值分解的数学公式定义:

SVD也是对矩阵进行分解,但是和特征分解不同,SVD并不要求要分解的矩阵为方阵。假设M是一个 m×n 阶矩阵,其中的元素全部属于域 K,也就是实数域或复数域,那么我们定义矩阵A的SVD为:

U是一个M * M的方阵(里面的向量是正交的,U里面的向量称为左奇异向量);

Σ是一个M * N的实数对角矩阵(对角线以外的元素都是0,对角线上的元素称为奇异值,常见的做法是为了奇异值由大而小排列);

VT(V的转置)是一个N * N的矩阵,里面的向量也是正交的,V里面的向量称为右奇异向量);

U 和 V 都是酉矩阵,满足U^{T}U=I,V^{T}V=I 。下图可以很形象的看出上面SVD的定义:

接下来的问题是如何求出SVD分解后的U,Σ,V这三个矩阵呢?方法就是借助特征值分解的方法。

1. 奇异值分解的具体方法 

1)V矩阵的计算方法

如果我们将一个任意mxn的矩阵A的转置和A做矩阵乘法,那么会得到n×n的一个方阵 A^{T}A 。既然 A^{T}A 是方阵,那么我们就可以进行特征分解,得到的特征值和特征向量满足下式:

image   

这样我们就可以得到矩阵 A^{T}A 的n个特征值和对应的n个特征向量v了。将 A^{T}A 的所有特征向量张成一个n×n的矩阵V,就是我们SVD公式里面的V矩阵了。一般我们将V中的每个特征向量叫做A的右奇异向量

2)U矩阵的计算方法

如果我们将A和A的转置做矩阵乘法,那么会得到m×m的一个方阵 AA^{T} 。既然 AA^{T} 是方阵,那么我们就可以进行特征分解,得到的特征值和特征向量满足下式:

这样我们就可以得到矩阵 AA^{T} 的m个特征值和对应的m个特征向量u了。将 AA^{T} 的所有特征向量张成一个m×m的矩阵U,就是我们SVD公式里面的U矩阵了。一般我们将U中的每个特征向量叫做A的左奇异向量

3)奇异值矩阵Σ计算方法

现在U和V我们都求出来了,接下来求奇异值矩阵Σ,由于Σ除了对角线上是奇异值其他位置都是0,那我们只需要求出每个奇异值σ就可以了。

推导公式如下:

这样我们可以求出我们的每个奇异值,进而求出奇异值矩阵Σ。

4)奇异值分解和特征值分解的关联关系

上面我们说到,A^{T}A 的特征向量组成的就是我们SVD中的V矩阵,而AA^{T} 的特征向量组成的就是我们SVD中的U矩阵,这里我们来推导证明下,我们以V矩阵的证明为例:

上式证明中,使用了 U^{T}U=I,\Sigma^{T}\Sigma=\Sigma^{2} ,所以等式可以化简。

可以看出 A^{T}A 的特征向量组成的的确就是我们SVD中的V矩阵。类似的方法可以得到 AA^{T} 的特征向量组成的就是我们SVD中的U矩阵。

进一步我们还可以看出:特征值矩阵等于奇异值矩阵的平方,也就是说特征值和奇异值满足如下关系:

这样也就是说,我们可以不用 \sigma_{i}=\frac{Av_{i}}{u_{i}} 来计算奇异值,也可以通过求出 A^{T}A 的特征值取平方根来求奇异值。

0x3:SVD分解计算举例

这里我们用一个简单的例子来说明矩阵是如何进行奇异值分解的。我们的矩阵A定义为:

根据上一章节所述,我们先来求V(右奇异向量),我们通过将矩阵A的转置乘以矩阵A,得到一个2级方阵,从而可以借助特征值分解来计算奇异值分解。 

通过计算矩阵的特征多项式,并求根的方式来计算矩阵的特征值,得到了特征值就得到了特征向量:

,令方程为零,得到方程的根为,对特征值开方既得奇异值:

将奇异值放置在对角位置既得奇异值矩阵:

矩阵的特征向量为:

两个特征向量合起来就是V(右奇异向量)

继续求U(左奇异向量)

合起来既得:

综上得

整个矩阵奇异值分解完毕。

0x4:SVD分解的几何意义

还是继续上个小节的实际例子,我们从几何的视角来理解下SVD的概念:

上图展示了一个奇异值分解过程,奇异值分解实际上把矩阵的变换分为了三部分:

  • 旋转
  • 拉伸
  • 投影

方阵没有投影,这里忽略,我们来看旋转和拉伸:

单位圆先被旋转(90°),x和y坐标系互换,是没有形变的。

再进行拉伸,这里决定了单位圆的形状,奇异值分别是椭圆的长轴和短轴,σ1和σ2代表了在各自方向上的“拉伸力度”:

可以看到,通过两个方向上的拉伸力的综合,棕色的圆被拉伸成了一个绿色椭圆,注意,因为之前已经进行了坐标系的变换,所以此时是一个竖着的椭圆。

最后,再次旋转90°,被旋转到最终的位置,成为一个横着的椭圆,这一过程也没有发生形变,最终就是上图上看到的绿色椭圆了。

0x4:奇异值的一些有用的性质 

对于奇异值,它跟我们特征分解中的特征值类似,在奇异值矩阵中也是按照从大到小排列,而且奇异值的减少特别的快,在很多情况下,前10%甚至1%的奇异值的和就占了全部的奇异值之和的99%以上的比例。也就是说,我们也可以用最大的k个的奇异值和对应的左右奇异向量来近似描述矩阵。

也就是说:

其中k要比n小很多,也就是一个大的矩阵A可以用三个小的矩阵来表示。

如下图所示,现在我们的矩阵A只需要灰色的部分的三个小矩阵就可以近似描述了。

由于这个重要的性质,SVD可以用于PCA降维,来做数据压缩和去噪。也可以用于推荐算法,将用户和喜好对应的矩阵做特征分解,进而得到隐含的用户需求来做推荐。同时也可以用于NLP中的算法,比如潜在语义索引(LSI)。

应用总是五花八门,千变万化的,笔者这里希望传递的更多是SVD中包含的内核思想:不管是image grid数据,还是NLP序列向量数据,其抽象本质都是matrix矩阵,而matrix作为数据的一种载体,其内部是包含着一定的信息熵的,SVD做的事情就是,找到数据矩阵中混乱度最大的k个方向(混乱度越大,信息熵越大),这是信息熵的本质。更高层来说,这k个方向可以代表原始数据matrix的特征空间,也可以代表着不同层次的语义信息,取决于我们用什么视角来看它,以及上下游stacking的其他model,这是非常灵活的

Relevant Link:

https://www.cnblogs.com/fuleying/p/4466326.html
https://www.geogebra.org/m/kfyxcbee
https://www.geogebra.org/m/yJK7STqg 

 

6. 二次型、矩阵的合同

0x1:二次型的基本定义

系数在数域K中的n个变量的一个二次齐次多项式,称为数域K上的一个n元二次型,它的一般形式是:

上式也可以写成:

,其中

0x2:二次型矩阵

把上面二次型方程中的系数排成一个n级矩阵A:

把A称为二次型的矩阵,它是堆成矩阵。

显然,二次型的矩阵是唯一的,它的主对角元依次是的稀疏,它的(i,j)元是的系数的一半。

令:

则二次型可以写成:

其中A是二次型的矩阵。

用一个例子来说明:

更一般的:

写成矩阵相乘的形式:

0x3:为什么需要二次型

简单来说,二次型就是通过矩阵研究二次函数,借助矩阵运算的各种高效技巧,可以对各种复杂的二次函数具备更好的处理能力和直观的表现。

设二次曲面S在直角坐标系1中的方程为:

这是一个什么样的曲面呢?

解决这个问题的思路是:依赖坐标系线性变换的不变性,作直角坐标变换,使得在直角坐标系2中,S的方程不包含交叉项,只包含单纯的平方项,那么久可以看出S是什么二次曲面。

设直角坐标系的变换公式为:

,其中T一定是正交矩阵。

从而,上面二次曲面方程可以重写为:

将上式中,中间的3级矩阵记作A,得:

为了使新的方程中不出现交叉项,只要使矩阵为对角矩阵,又由于,因此也就是只要使为对角矩阵,这就是希望A能对角化,并且要找一个正交矩阵T,使A对角化。

所以二次型的主要研究问题就是:对于实数域上的对称矩阵A,能否找到正交矩阵T,使为对角矩阵?

我们下面来看一些具体的例子。

1. 对圆锥曲面的变形

下面展示了一个原的方程:

我们通过线性变换,改变一下二次型矩阵:

可以看到,椭圆和圆之间是线性关系,通过矩阵变换就可以从圆变为椭圆,我们继续扭曲:

我们发现双曲线和圆之间也是线性关系,准确的说是仿射的。

其实圆、椭圆、双曲线之间关系很紧密的,统称为圆锥曲线,都是圆锥体和平面的交线。上面这几个例子可能还是不太直观,我们来看下面的动图:

Relevant Link:

https://www.matongxue.com/madocs/271.html  

 

6. PCA(Principal Component Analysis 主成分分析)

有了上面对矩阵相似分解,特征值/特征向量的讨论,我们这节开始讨论PCA主成分分析,这是一种降维手段,了解其底层数学原理后读者朋友会发现其并不复杂,甚至很简单质朴。

0x1:关于冗余特征的一个简单例子 - 只有一个特征值的情况

先搁置PCA的具体原理概念在一边,这个小节我们通过一个生活中的常见例子,展示在只有一个特征值的情况下进行降维处理,这里所谓的只有一个特征值指的是降维到可以得到一个1维向量,我们后面会解释原理,读者朋友这里姑且忽略。

假设现在我们采集了一组房屋的两个方面的特征数据,【房价、面积】,可以看出两者完全正相关,有一列其实是多余的: 

把这个二维数据画在坐标轴上,横纵坐标分别为“房价”、“面积”,可以看出它们排列为一条直线:

我们现在旋转坐标系,让横坐标和这条直线重合:

旋转后的坐标系,横纵坐标不再代表“房价”、“面积”了(因为我们转换了向量基),新的坐标系而是两者的混合,这里把它们称作“主元1”、“主元2”,坐标值很容易用勾股定理计算出来,比如a在“主元1”的坐标值为下图,同时,很显然a 在“主元2”上的坐标为0:

我们把所有的数据点换算到新的坐标系上:

可以看到,因为“主元2”全都为0,完全是多余的,我们只需要“主元1”就够了,这样就把数据降为了一维,而且没有丢失任何信息:

0x2:稍微复杂一些的一个例子 - 特征值数量为二

接下来我们继续上面那个例子,但是我们稍微修改一下采集到的数据集,将问题变得更贴近现实业务场景一些,标题这里写特征值数量为二,读者朋友也可以先忽略,我们后面解释:

观察这组数据,很明显,第二列和第一列依然存在着某种正相关,但是这回已经不像上个小节那样是完全相等了,而是有了一些自己的变化。

把这个二维数据画在坐标轴上,横纵坐标分别为“房价”、“面积”,虽然数据看起来很接近一条直线,但是终究不在一条直线上:

现在新问题来了,这个情况要怎么降维呢?降维后的新的数据又是什么样的呢?

要回答这个问题,我们需要将思路往回倒退,回到向量在坐标系(向量基)中的表示,以及坐标系的变换这个问题上来。

先讨论只有一个点的情况,在某坐标系有一个点,\boldsymbol{a}=\begin{pmatrix}x\\y\end{pmatrix} ,它表示在该坐标系向量基\boldsymbol{e_1},\boldsymbol{e_2} 的线性组合:

当使用不同的线性变换对坐标系进行线性变换时,在不同坐标系中,x,y 的值会有所不同(不同的坐标系代表了不同的视角),如下图:

注意,\boldsymbol{a} 到原点的距离d 不会因为坐标系改变而改变: 

所以,在某坐标系下分配给x 较多,那么分配给y 的就必然较少,反之亦然。最极端的情况是,在某个坐标系下,全部分配给了x ,使得y=0 ,如下图:

上图对应的就是上个小节我们的降维结果,将二维的数据完全无损地降维到了一维。

但是遗憾的是,这个情况在本例的数据集中是不可能的。因为其实上个小节的数据集中隐含了一个重要条件我没有指出,即因为上个小节中的数据集的特征值只有1个,特征向量也只有一个,所以可以从2列直接降维到1列。但是这个小节中的数据特征向量有2个,特征向量也有2个,是不可能无损地降维到一维的。

为什么无法无损的降维到1维呢?我们将问题扩展到两个点\boldsymbol{a}=\begin{pmatrix}x_1\\y_1\end{pmatrix},\boldsymbol{b}=\begin{pmatrix}x_2\\y_2\end{pmatrix} ,就能说清楚原因了,两个点以上本质是类似的,下图给出一个图示:

可以直观地看到,我们无法找到一个完美的坐标系,使得a,b点都同时能落在同一个坐标系上,也就说,无法无损地降维到1维。

注意,这里是说”不能无损地降到1维“,如果可以接受数据损失,我们强行地丢弃另一个方向上的向量,也是可以降到1维的。

好了!现在我们已经认命了,即无法无损地降维到1维,这个时候,我们要改变策略,我们希望找到一组新的坐标系,这组坐标系能够满足以下几个特征:

  • 坐标系的维数尽可能地少,越少越好,例如将向量分量尽量多分配给x_1,x_2 ,少分配给y_1,y_2 的坐标系
  • 原始数据投影到新坐标系后,数据的失真尽可能地小,还能最大程度地保留原始数据中的信息熵

如何实现上述两个目标呢?这就是接下来要详述的PCA算法的原理。

好了,讨论到了这里,感性部分已经阐述完毕,总结一下就是:

  • 原始数据都出在一个具体的向量基组中
  • 降维的核心就是寻找一种线性变换方案(一个矩阵),来对原始的向量基进行转换,以求数据能够尽量集中地分配在少部分的向量基中
  • 降维所谓的降,就是指有选择的丢弃部分分配向量分量较少的向量基,保留分配向量分量较多的向量基,使得总体的信息丢失尽可能地少

以上部分都是一些感性化的讨论,不具备实际操作性。接下来我们要开始讨论PCA的具体理论,看看如何用数学方式来描述上面的理论并进行实际降维操作。

0x3:PCA公式推导

假设有如下数据:

上图中的2个数据向量,在初始坐标系下(也就是自然基\boldsymbol{e_1}=\begin{pmatrix}1\\0\end{pmatrix},\boldsymbol{e_2}=\begin{pmatrix}0\\1\end{pmatrix} )下坐标值为:

图示如下:

随着坐标系的不同,X_1,X_2 的值会不断变化(在不同坐标系下的不同表现)。

现在我们的目标是:想尽量多分配给X_1,X_2 ,尽可能少的分配给Y1、Y2,根据最小二乘法的思想,就是让下式成立:

X_1^2+X_2^2=\sum_{i=0}^2 X_i^2\ \ 最大\\

要实现上式,就需要找到一个新的坐标系,我们假设新的坐标系如下:

\boldsymbol{e_1}=\begin{pmatrix}e_{11}\\e_{12}\end{pmatrix}\quad \boldsymbol{e_2}=\begin{pmatrix}e_{21}\\e_{22}\end{pmatrix}\\

则原始数据在新坐标系下其中一个坐标系e1的投影为:

带上上面最小二乘公式有:

上式其实是一个二次型:

这里矩阵P 就是二次型,是一个对称矩阵,可以进行如下的奇异值分解:

其中,U 为正交矩阵,即UU^\mathrm{T}=I 。

\Sigma 是对角矩阵:

\Sigma=\begin{pmatrix}\sigma_1&0\\0&\sigma_2\end{pmatrix}\\

其中,\sigma_1,\sigma_2 是奇异值,\sigma_1 > \sigma_2 。

P 代回原方程:

因为U 是正交矩阵,所以令:

所得的\boldsymbol{n} 也是单位向量,即:

\boldsymbol{n}=\begin{pmatrix}n_1\\n_2\end{pmatrix}\implies n_1^2+n_2^2=1\\

继续代回原方程:

至此,最初求最大值的问题就转化为了:

可以用拉格朗日乘子法计算上述条件极值,结果是当n_1=1,n_2=0 时取到极值。

因此可以推出要寻找的主元1,即:

\boldsymbol{n}=\begin{pmatrix}1\\0\end{pmatrix}=U^\mathrm{T}\boldsymbol{e_1}\implies \boldsymbol{e_1}=U\begin{pmatrix}1\\0\end{pmatrix}\\

即:

同样的思路可以求出:

\boldsymbol{e_2}=\begin{cases}P=U\Sigma U^\mathrm{T}\\\\最小奇异值\sigma_2对应的奇异向量\end{cases}\\

0x4:基于协方差矩阵的优化目标 - 找到损失最低的变换基

我们上个章节通过一个具体的例子,推导了如何得到一个最优的降维坐标系。

但是请读者朋友注意的是,上一节对最优化目标的描述是:想尽量多分配给X_1,X_2 ,尽可能少的分配给Y1、Y2。这句话还是过于直观和感性,我们这节用协方差矩阵的方式来形式化地定义出PCA降维的目标是什么。

PCA降维的本质是基变换,如果我们必须使用一组新的基来表示原始数据,又希望尽量保留原始的信息(保留原始数据的概率分布),我们应该如何选择?

而不同的基变换对应了不同的投影方向,我们希望投影后的投影值尽可能分散。数据越分散,可分性就越强,可分性越强,概率分布保存的就越完整。 这种分散程度,可以用数学上的方差来表述

但是光有方差还是不够的,考虑三维降到二维问题。首先我们希望找到一个方向使得投影后方差最大,这样就完成了第一个方向的选择,继而我们选择第二个投影方向。如果我们还是单纯只选择方差最大的方向,很明显,这个方向与第一个方向应该是“几乎重合在一起”,显然这样的维度是没有用的,因此,应该有其他约束条件。

从直观上说,让两个新的向量基尽可能表示更多的原始信息,同时我们不希望它们之间存在(线性)相关性的,因为相关性意味着两个字段不是完全独立,必然存在重复表示的信息。
数学上可以用两个字段的协方差表示其相关性:

,这里Y表示第二个向量基。

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

用协方差矩阵来统一表示同一个向量基内数据的离散方差不同向量基间相关性协方差这两个数学概念:

根据上述推导,我们要达到的优化目标,用数学语言等价于:将协方差矩阵对角化,即除对角线(方差要尽可能大)外的其它元素化为0(协方差为0),并且在对角线上将元素按大小从上到下排列(这样便于我们筛选出方差最大的k个向量基),这样我们就达到了优化目的。

P,Q 都可以进行奇异值分解:

P=U\begin{pmatrix}\sigma_1&0\\0&\sigma_2\end{pmatrix} U^\mathrm{T}\quad Q=\frac{1}{n}P=U\begin{pmatrix}\frac{\sigma_1}{n}&0\\0&\frac{\sigma_2}{n}\end{pmatrix} U^\mathrm{T}\\

矩阵P上个章节解释过了,指的是原始数据的二次型矩阵。协方差矩阵Q 的奇异值分解和P 相差无几,只是奇异值缩小了n 倍,但是不妨碍奇异值之间的大小关系,所以在实际问题中,往往都是直接分解协方差矩阵Q 。

0x5:PCA陷阱 - 注意PCA的假设前提条件!

1. PCA只能去除数据中的线性相关性

考虑采样y=sin(x)半个周期所获得的数据: 

虽然上述数据显然是线性不相关的(平均而言,当x值增加时,y值增加的和它减少的一样多),因此对应于一个对角的协方差矩阵,但这两个变量之间仍然有明确的非线性关系。

一般情况下,PCA会减少数据中的线性相关性,但不会删除统计相关。如果我们的数据集是非线性相关的,那么可能不适合用PCA进行降维。 

2. 可区分信息不存在于最大方差方向上 

PCA的一个最大假设是,最能区别的信息通过特征空间中的最大方差捕获。因为最大方差方向编码的信息最多,这很可能是真的。

然而,有一些情况,其中的区别信息实际驻留在最小方差的方向,使得PCA可以大大损害分类性能。举个例子,考虑下图的两种情况(将二维特征空间降到1D表示): 

  • 在第一幅图中,如果保留方差最大的方向,那么基本上得不到任何可区分的信息。因为PCA是根据方差大小来选择保留哪些方向的特征的,如果因为数据集中存在噪声而导致某个方向的方差特别大,那么PCA很可能会错误的保留这个噪声方向,而割舍真正有价值的相关性维度
  • 在第二张图中PCA则完美发挥作用,这是PCA的典型场景。

如果最能区别的信息包含在较小的特征向量中,应用PCA实际上可能恶化维数的诅咒,因为现在需要更复杂的分类模型(例如非线性分类器)来分类低维问题。在这种情况下,其他降维的方法可能会感兴趣,如线性判别分析(LDA),其试图找到能够一个映射向量,该向量最优的分开两个类别。 

0x6:PCA降维的应用 

1. 用SVD代替特征分解的庞大运算量

PCA降维,需要找到样本协方差矩阵 X^{T}X 的最大的k个特征向量,然后用这最大的k个特征向量张成的矩阵来做低维投影降维。可以看出,在这个过程中需要先求出协方差矩阵 X^{T}X ,当样本数多样本特征数也多的时候,这个计算量是很大的。

注意到我们的SVD也可以得到协方差矩阵 X^{T}X 最大的k个特征向量张成的矩阵,而且同时SVD有个好处,有一些SVD的实现算法可以不用先求出协方差矩阵 X^{T}X ,也能求出我们的右奇异矩阵V。也就是说,我们的PCA算法可以不用做特征分解,而是做SVD来完成。这个方法在样本量很大的时候很有效。实际上,scikit-learn的PCA算法的背后真正的实现就是用的SVD,而不是我们我们认为的暴力特征分解。

2. 左奇异分解矩阵 - 行压缩

另一方面,注意到PCA仅仅使用了我们SVD的右奇异矩阵,没有使用左奇异矩阵,那么左奇异矩阵有什么用呢?

假设我们的样本是m×n的矩阵X,如果我们通过SVD找到了矩阵 XX^{T} 最大的k个特征向量张成的m×k维矩阵U,则我们如果进行如下处理:

可以得到一个k x n的矩阵X‘,这个矩阵和我们原来的m×n维样本矩阵X相比,行数从m减到了k,可见对行数进行了压缩,左奇异矩阵可以用于行数的压缩。

3. 右奇异矩阵 - 列维数压缩

右奇异矩阵可以用于列数即特征维度的压缩,也就是我们的PCA降维。  

Relevant Link:    

https://www.cnblogs.com/LittleHann/p/6558575.html#_lab2_1_7
https://zhuanlan.zhihu.com/p/31386807
https://www.matongxue.com/madocs/306.html
https://www.cnblogs.com/LittleHann/p/6558575.html#4208658
https://www.zhihu.com/question/41120789/answer/481966094

   

posted @ 2019-06-07 17:45 郑瀚Andrew.Hann 阅读(...) 评论(...) 编辑 收藏