绿色健康小清新

耐得住寂寞,守得住繁华

LDA(Linear Discriminant Analysis)的原理详解

什么是LDA

LDA是一种监督学习的降维技术,也就是说它的数据集的每个样本是有类别输出的。这点和PCA不同。PCA是不考虑样本类别输出的无监督降维技术。LDA的思想可以用一句话概括,就是“投影后类内方差最小,类间方差最大”。什么意思呢? 我们要将数据在低维度上进行投影,投影后希望每一种类别数据的投影点尽可能的接近,而不同类别的数据的类别中心之间的距离尽可能的大。

我们来看一个例子,假设我们有两类数据 分别为红色和蓝色,如下图所示,这些数据特征是二维的,我们希望将这些数据投影到一维的一条直线,让每一种类别数据的投影点尽可能的接近,而红色和蓝色数据中心之间的距离尽可能的大。

上图提供了两种投影方式,哪一种能更好的满足我们的标准呢?从直观上可以看出,右图要比左图的投影效果好,因为右图的红色数据和蓝色数据各个较为集中,且类别之间的距离明显。左图则在边界处数据混杂。以上就是LDA的主要思想了,当然在实际应用中,我们的数据是多个类别的,我们的原始数据一般也是超过二维的,投影后的也一般不是直线,而是一个低维的超平面。

再看一个例子

左边的投影效果显然不如右边图的投影效果。假如我们知道了两类数据的中心点u1u_1u2u_2,以及两类数据的离散度S1S_1以及S2S_2,那我们现在要做的就是使u1u22S12+S22\frac{|u_1 - u_2|^2}{S_1^2 + S_2^2}尽量大。即最大化两类数据的距离,最小化两类数据各自的离散度。

那我们应该如何进行计算呢???

在我们将上面直观的内容转化为可以度量的问题之前,我们先了解些必要的数学基础知识,这些在后面讲解具体LDA原理时会用到。

瑞利商(Rayleigh quotient)与广义瑞利商(genralized Rayleigh quotient)

我们首先来看看瑞利商的定义。瑞利商是指这样的函数R(A,x):

R(A,x)=xHAxxHxR(A,x) = \frac{x^HAx}{x^Hx}

其中x为非零向量,而A为n×n的Hermitan矩阵。所谓的Hermitan矩阵就是满足共轭转置矩阵和自己相等的矩阵,即AH=A。如果我们的矩阵A是实矩阵,则满足AT=A的矩阵即为Hermitan矩阵。

瑞利商R(A,x)有一个非常重要的性质,即它的最大值等于矩阵A最大的特征值,而最小值等于矩阵A的最小的特征值,也就是满足

λminxHAHxHxλmaxλ_{min}≤\frac{x^HAH}{x^Hx}≤λ_{max}

具体的证明这里就不给出了,感兴趣的可以看瑞利定理的定理证明

当向量x是标准正交基时,即满足xHx=1时,瑞利商退化为:R(A,x)=xHAx,这个形式在谱聚类和PCA中都有出现。

以上就是瑞利商的内容,现在我们再看看广义瑞利商。广义瑞利商是指这样的函数R(A,B,x):

R(A,X)=xHAHxHBxR(A,X) = \frac{x^HAH}{x^HBx}

其中x为非零向量,而A,B为n×n的Hermitan矩阵。B为正定矩阵。它的最大值和最小值是什么呢?其实我们只要通过将其通过标准化就可以转化为瑞利商的格式。我们令x=B−1/2x′,则分母转化为:

xHBx=xH(B12)HBB12x=xHB12BB12x=xHxx^HBx = x'^H(B^{-\frac{1}{2}})^HBB^{-\frac{1}{2}}x' = x'^HB^{-\frac{1}{2}}BB^{-\frac{1}{2}}x = x'^Hx'

而分子转化为:

xHAx=xHB12AB12xx^HAx = x'^HB^{-\frac{1}{2}}AB^{-\frac{1}{2}}x'

此时我们的R(A,B,x)转化为R(A,B,x′):

R(A,B,X)=xHB12AB12xxHxR(A,B,X') = \frac{ x'^HB^{-\frac{1}{2}}AB^{-\frac{1}{2}}x'}{x'^Hx'}

利用前面的瑞利商的性质,我们可以很快的知道,R(A,B,x′)的最大值为矩阵B12AB12B^{-\frac{1}{2}}AB^{-\frac{1}{2}}的最大特征值,或者说矩阵B1AB^{-1}A的最大特征值,而最小值为矩阵B1AB^{-1}A的最小特征值。

二类LDA原理

现在我们回到LDA的原理上,我们在第一节说讲到了LDA希望投影后希望同一种类别数据的投影点尽可能的接近,而不同类别的数据的类别中心之间的距离尽可能的大,但是这只是一个感官的度量。现在我们首先从比较简单的二类LDA入手,严谨的分析LDA的原理。

假设我们的数据集D=(x1,y1),(x2,y2),...,((xm,ym))D={(x_1,y_1),(x_2,y_2),...,((x_m,y_m))},其中任意样本xi为n维向量,yiy_i∈{0,1}。我们定义NjN_j(j=0,1)为第j类样本的个数,XjX_j(j=0,1)为第j类样本的集合,而μjμ_j(j=0,1)为第j类样本的均值向量,定义ΣjΣ_j(j=0,1)为第j类样本的协方差矩阵(严格说是缺少分母部分的协方差矩阵)。

μjμ_j的表达式为:

uj=1NjxXjx  (j=0,1)u_j = \frac{1}{N_j}\sum_{x∈X_j}x \ \ (j = 0,1)

j\sum_j的表达式为:

j=xXj(xuj)(xuj)T(j=0,1)\sum_j = \sum_{x∈X_j}(x - u_j)(x - u_j)^T (j = 0,1)

由于是两类数据,因此我们只需要将数据投影到一条直线上即可。假设我们的投影直线是向量w,则对任意一个样本本xix_i,它在直线w的投影为wTxiw^Tx_i,对于我们的两个类别的中心点μ0μ_0,μ1μ_1,在在直线w的投影为wTμ0w^Tμ_0wTμ1w^Tμ_1。由于LDA需要让不同类别的数据的类别中心之间的距离尽可能的大,也就是我们要最大化wTμ0wTμ122||w^Tμ_0−w^Tμ_1||^2_2,同时我们希望同一种类别数据的投影点尽可能的接近,也就是要同类样本投影点的协方差wT0wwT1ww^T\sum_0w和w^T\sum_1w尽可能的小,即最小化wT0wwT1ww^T\sum_0w和w^T\sum_1w。综上所述,我们的优化目标为:

argmaxJ(w)=wTu0wTu122wT0w+wT1w=wT(u0u1)(u0u1)TwwT(0+1)warg max J(w) = \frac{||w^Tu_0 - w^Tu_1||_2^2}{w^T\sum_0w + w^T\sum_1w} = \frac{w^T(u_0 - u_1)(u_0-u_1)^Tw}{w^T(\sum_0 + \sum_1)w}

我们一般定义类内散度矩阵SwS_w为:

Sw=0+1=xX0(xu0)(xu0)T+xX1(xu1)(xu1)TS_w = \sum_0 + \sum_1 = \sum_{x∈X_0}(x-u_0)(x-u_0)^T + \sum_{x∈X_1}(x-u_1)(x-u_1)^T

同时定义类间散度矩阵SbS_b为:

Sb=(u0u1)(u0u1)TS_b = (u_0 - u_1)(u_0 -u_1)^T

这样我们的优化目标重写为:

arg maxJ(w)=wTSbwwTSwwarg \ max J(w) = \frac{w^TS_bw}{w^TS_ww}

仔细一看上式,这不就是我们的广义瑞利商嘛!这就简单了,利用我们第二节讲到的广义瑞利商的性质,我们知道我们的J(w′)最大值为矩阵Sw12SbSw12S_w^{-\frac{1}{2}}S_bS_w^{-\frac{1}{2}}的最大特征值,而对应的w′为Sw12SbSw12S_w^{-\frac{1}{2}}S_bS_w^{-\frac{1}{2}}的最大特征值对应的特征向量! 而Sw1SbS^{−1}_wS_b的特征值和Sw12SbSw12S_w^{-\frac{1}{2}}S_bS_w^{-\frac{1}{2}}的特征值相同,Sw1SbS^{−1}_wS_b的特征向量w和Sw12SbSw12S_w^{-\frac{1}{2}}S_bS_w^{-\frac{1}{2}}的特征向量w′满足w=Sw12ww=S^{−\frac{1}{2}}_ww′的关系!

注意到对于二类的时候,SbwS_bw的方向恒平行于μ0μ1μ_0−μ_1,不妨令Sbw=λ(μ0μ1)S_bw=λ(μ0−μ1),将其带入 (Sw1Sb)w=λw(S^{−1}_wS_b)w=λw,可以得到w=Sw1(μ0μ1)w=S^{−1}_w(μ_0−μ_1), 也就是说我们只要求出原始二类样本的均值和方差就可以确定最佳的投影方向w了。

一定要自己推一下,才好理解!!! 都是一些线代基础知识。

算法的流程

输入:数据集D={(x1,y1),(x2,y2),...,((xm,ym))}D=\{(x_1,y_1),(x_2,y_2),...,((x_m,y_m))\},其中任意样本xi为n维向量,yi{C1,C2,...,Ck}y_i∈\{C_1,C_2,...,C_k\},降维到的维度d。

输出:降维后的样本集DD′

  1. 计算类内散度矩阵SwS_w

  2. 计算类间散度矩阵SbS_b

  3. 计算矩阵Sw1SbS^{−1}_wS_b

  4. 计算Sw1SbS^{−1}_wS_b的最大的d个特征值和对应的d个特征向量(w1,w2,…wd),得到投影矩阵WW

  5. 对样本集中的每一个样本特征xix_i,转化为新的样本zi=WTxiz_i=W^Tx_i

  6. 得到输出样本集D=(z1,y1),(z2,y2),...,((zm,ym))D′={(z_1,y_1),(z_2,y_2),...,((z_m,y_m))}

以上就是使用LDA进行降维的算法流程。实际上LDA除了可以用于降维以外,还可以用于分类。一个常见的LDA分类基本思想是假设各个类别的样本数据符合高斯分布,这样利用LDA进行投影后,可以利用极大似然估计计算各个类别投影数据的均值和方差,进而得到该类别高斯分布的概率密度函数。当一个新的样本到来后,我们可以将它投影,然后将投影后的样本特征分别带入各个类别的高斯分布概率密度函数,计算它属于这个类别的概率,最大的概率对应的类别即为预测类别。


二类LDA例子(与PCA进行对比)

现在有两个数据集:

C1=[(1,2);(2,3);(3,3);(4,5);(5,5)]          C2=[(1,0);(2,1);(3,1);(3,2);(5,3);(6,5)]C_1 = [(1,2);(2,3);(3,3);(4,5);{(5,5)}] \\ \ \ \ \ \ \ \ \ \ \ C_2 = [(1,0);(2,1);(3,1);(3,2);(5,3);(6,5)]

PCA降维

[C1;C2][C_1;C_2]的协方差:

Z={2.7636   2.25452.2545   3.0182}Z = \left\{ \begin{matrix} 2.7636 \ \ \ 2.2545 \\ 2.2545 \ \ \ 3.0182 \end{matrix} \right\}

Z的特征向量和特征值

V={0.72680.68690.68690.7268}V = \left\{ \begin{matrix} -0.7268 & 0.6869 \\ 0.6869 & 0.7268 \end{matrix}\right\}

Z={0.6328005.1490}Z = \left\{ \begin{matrix} 0.6328 & 0\\ 0 & 5.1490\end{matrix} \right\}

看出我们有两个特征值,0.632< 5.1490,选取5.1490所对应的特征向量为投影方向,即:

[0.6869,0.7268]T[0.6869,0.7268]^T


LDA降维

每一类数据的均值:

u1=mean(C1)=[3.0,3.6]Tu_1 = mean(C_1) = [3.0,3.6]^T

u2=mean(C2)=[3.3,2.0]Tu_2 = mean(C_2) = [3.3,2.0]^T

进而算出SBS_B

Sb=(u1u2)(u1u2)T={0.110.530.532.56}S_b = (u_1-u_2)(u_1-u_2)^T = \left\{ \begin{matrix} 0.11 & -0.53 \\ -0.53 & 2.56\end{matrix}\right\}

每一类数据的离散度:

S1=4×cov(C1)={10887.2}S_1 = 4 × cov(C_1) = \left\{ \begin{matrix} 10 & 8 \\8 & 7.2 \end{matrix}\right\}

S2=5×cov(C2)={17.3161616}S_2 = 5 × cov(C_2) = \left\{ \begin{matrix} 17.3 & 16 \\16 & 16 \end{matrix}\right\}

进而算出SwS_w

Sw=S1+S2={27.3242423.2}S_w = S_1 + S_2 = \left\{ \begin{matrix} 27.3 & 24 \\ 24 & 23.2 \end{matrix}\right\}


方法一:

计算Sw1SbS_w^{-1}S_b

Sw1Sb={0.261.270.301.42}S_w^{-1}S_b = \left\{ \begin{matrix} 0.26 & -1.27 \\-0.30 & 1.42\end{matrix}\right\}

Sw1SbS_w^{-1}S_b的特征值和特征向量

V={0.980.670.200.75}V = \left\{ \begin{matrix} -0.98 & 0.67 \\ -0.20 & -0.75 \end{matrix}\right\}

D={0001.69}D = \left\{ \begin{matrix} 0 & 0 \\ 0 & 1.69 \end{matrix}\right\}

看出我们有两个特征值,0<1.69,选取1.69所对应的特征向量为投影方向,即:

[0.6656,0.7463]T[0.6656,-0.7463]^T


方法二:

计算Sw1(u1u2)TS_w^{-1}(u_1-u_2)^T

Sw1(u1u2)T=[0.7936,0.8899]TS_w^{-1}(u_1-u_2)^T = [-0.7936, 0.8899]^T

进行标准化:

[0.6656,0.7463]T[0.6656,-0.7463]^T


PCA和LDA图像比较

不难看出PCA依然是那个PCA,在处理多类数据的降维问题上,效果比LDA差不少!



多类LDA原理

有了二类LDA的基础,我们再来看看多类别LDA的原理。

假设我们的数据集D=(x1,y1),(x2,y2),...,((xm,ym))D={(x_1,y_1),(x_2,y_2),...,((x_m,y_m))},其中任意样本xi为n维向量,yiC1,C2,...,Cky_i∈{C_1,C_2,...,C_k}。我们定义Nj(j=1,2...k)N_j(j=1,2...k)为第j类样本的个数,Xj(j=1,2...k)X_j(j=1,2...k)为第j类样本的集合,而μj(j=1,2...k)μ_j(j=1,2...k)为第j类样本的均值向量,定义Σj(j=1,2...k)Σ_j(j=1,2...k)为第j类样本的协方差矩阵。在二类LDA里面定义的公式可以很容易的类推到多类LDA。

由于我们是多类向低维投影,则此时投影到的低维空间就不是一条直线,而是一个超平面了。假设我们投影到的低维空间的维度为d,对应的基向量为(w1,w2,...wd)(w_1,w_2,...w_d),基向量组成的矩阵为W, 它是一个n×d的矩阵。

此时我们的优化目标应该可以变成为:

WTSbWWTSwW\frac{W^TS_bW}{W^TS_wW}

其中Sb=j=1kNj(uju)(uju)TS_b = \sum^k_{j=1}N_j(u_j - u)(u_j-u)^T,u为所有样本均值向量;Sw=j=1kSwj=j=1kxXj(xuj)(xuj)TS_w = \sum^k_{j=1} S_{wj} = \sum^k_{j=1}\sum_{x∈X_j}(x-u_j)(x-u_j)^T

在这里插入图片描述

但是有一个问题,就是WTSbWW^TS_bWWTSwWW^TS_wW都是矩阵,不是标量,无法作为一个标量函数来优化!也就是说,我们无法直接用二类LDA的优化方法,怎么办呢?一般来说,我们可以用其他的一些替代优化目标来实现。

常见的一个LDA多类优化目标函数定义为:

arg max  J(w)=diagWTSbWdiagWTSwWarg \ max \ \ J(w) = \frac{\prod_{diag}W^TS_bW}{\prod_{diag}W^TS_wW}

其中diagA\prod_{diag}A为A的主对角线元素的乘积,W为n×d的矩阵。

J(W)的优化过程可以转化为:

  J(w)=i=1dWiTSbwii=1dWiTSwwi=i=1dWiTSbwiWiTSwwi\ \ J(w) = \frac{\prod_{i=1}^{d}W_i^TS_bw_i}{\prod_{i=1}^{d}W_i^TS_ww_i} = \prod_{i=1}^{d} \frac{W_i^TS_bw_i}{W_i^TS_ww_i}

仔细观察上式最右边,这不就是广义瑞利商嘛!最大值是矩阵S1wSbS^{−1}wS_b的最大特征值,最大的d个值的乘积就是矩阵S1wSbS^{−1}wS_b的最大的d个特征值的乘积,此时对应的矩阵W为这最大的d个特征值对应的特征向量张成的矩阵。

由于W是一个利用了样本的类别得到的投影矩阵,因此它的降维到的维度d最大值为k-1。为什么最大维度不是类别数k呢?因为SbS_b中每个μjμμ_j−μ的秩为1,因此协方差矩阵相加后最大的秩为k(矩阵的秩小于等于各个相加矩阵的秩的和),但是由于如果我们知道前k-1个μjμ_j后,最后一个μkμ_k可以由前k-1个μjμ_j线性表示,因此SbS_b的秩最大为k-1,即特征向量最多有k-1个。

二类和多类公式的关系

看了上面的二类和多类,是不是会有一个疑问:为什么同一个东西会有两种表达式?两个公式其实是一样的。


LDA vs PCA

LDA用于降维,和PCA有很多相同,也有很多不同的地方,因此值得好好的比较一下两者的降维异同点。

首先我们看看相同点:

  1. 两者均可以对数据进行降维。

  2. 两者在降维时均使用了矩阵特征分解的思想。

  3. 两者都假设数据符合高斯分布。

我们接着看看不同点:

  1. LDA是有监督的降维方法,而PCA是无监督的降维方法

  2. LDA降维最多降到类别数k-1的维数,而PCA没有这个限制。

  3. LDA除了可以用于降维,还可以用于分类。

  4. LDA选择分类性能最好的投影方向,而PCA选择样本点投影具有最大方差的方向。

这点可以从下图形象的看出,在某些数据分布下LDA比PCA降维较优。

当然,某些某些数据分布下PCA比LDA降维较优,如下图所示:


LDA的弊端

看了上面那么多,是不是认为LDA很强大!!!但我们想一种情况,如果多类数据的中心点重合该怎么办,我们的SbS_b直接为0了!这种情况LDA就无法解决!



LDA算法小结

LDA算法既可以用来降维,又可以用来分类,但是目前来说,主要还是用于降维。在我们进行图像识别图像识别相关的数据分析时,LDA是一个有力的工具。下面总结下LDA算法的优缺点。

LDA算法的主要优点有:

  1. 在降维过程中可以使用类别的先验知识经验,而像PCA这样的无监督学习则无法使用类别先验知识。

  2. LDA在样本分类信息依赖均值而不是方差的时候,比PCA之类的算法较优。

LDA算法的主要缺点有:

  1. LDA不适合对非高斯分布样本进行降维,PCA也有这个问题。

  2. LDA降维最多降到类别数k-1的维数,如果我们降维的维度大于k-1,则不能使用LDA。当然目前有一些LDA的进化版算法可以绕过这个问题。

  3. LDA在样本分类信息依赖方差而不是均值的时候,降维效果不好。

  4. LDA可能过度拟合数据。


参考:

-------------本文结束感谢您的阅读-------------
六经蕴籍胸中久,一剑十年磨在手

欢迎关注我的其它发布渠道