1. 概要
线性判别分析是有监督的降维算法,也是一个分类算法,核心思想是:寻找一个投影矩阵使得类内样本的的投影尽可能近,类间样本投影尽可能远。
2. 理论基础
2.1 二分类
首先,先看该算法流程,再进一步理解,以下介绍均在二分类中考虑
结合线性判别分析的重要思想,我们的目标是
J = ∥ w T μ 0 − w T μ 1 ∥ 2 2 w T Σ 0 w + w T Σ 1 w = w T ( μ 0 − μ 1 ) ( μ 0 − μ 1 ) T w w T ( Σ 0 + Σ 1 ) w (1) J=\frac{\Vert\bm{w}^{\rm{T}}\bm{\mu}_0-\bm{w}^{\rm{T}}\bm{\mu}_1\Vert_2^2}{\bm{w}^{\rm{T}}\Sigma_0\bm{w}+\bm{w}^{\rm{T}}\Sigma{_1}\bm{w}} =\frac{\bm{w}^{{\rm{T}}}(\bm{\mu}_0-\bm{\mu}_1)(\bm{\mu}_0-\bm{\mu_1})^{\rm{T}}\bm{w}}{\bm{w}^{\rm{T}}(\Sigma_0+\Sigma_1)\bm{w}} \tag{1} J=wTΣ0w+wTΣ1w∥wTμ0−wTμ1∥22=wT(Σ0+Σ1)wwT(μ0−μ1)(μ0−μ1)Tw(1)
类内散度矩阵:
S w = Σ 0 + Σ 1 = ∑ x ∈ X 0 ( x − μ 0 ) ( x − μ 0 ) T + ∑ x ∈ X 1 ( x − μ 1 ) ( x − μ 1 ) T (2) S_w=\Sigma_0+\Sigma_1=\sum_{\bm{x}\in{X_0}}(\bm{x}-\bm{\mu}_0)(\bm{x}-\bm{\mu}_0)^{\rm{T}}+\sum_{\bm{x}\in{X_1}}(\bm{x}-\bm{\mu}_1)(\bm{x}-\bm{\mu}_1)^{\rm{T}} \tag{2} Sw=Σ0+Σ1=x∈X0∑(x−μ0)(x−μ0)T+x∈X1∑(x−μ1)(x−μ1)T(2)
类间散度矩阵:
S b = ( μ 0 − μ 1 ) ( μ 0 − μ 1 ) T (3) S_b= (\bm{\mu}_0-\bm{\mu}_1)(\bm{\mu}_0-\bm{\mu}_1)^{\rm{T}}\tag{3} Sb=(μ0−μ1)(μ0−μ1)T(3)
因此公式(1)可以改写为:
J = w T S b w w T S w w (4) J=\frac{\bm{w^{\rm{T}}S_bw}}{\bm{w^{\rm{T}}S_w}\bm{w}}\tag{4} J=wTSwwwTSbw(4)
求解该问题具体细节参考西瓜书,直接给出结果: w = S w − 1 ( μ 0 − μ 1 ) \bm{w}=S^{-1}_w(\bm{\mu}_0-\bm{\mu_1}) w=Sw−1(μ0−μ1).
西瓜书中介绍这里的矩阵求逆,考虑到数值解的稳定性,通常采用奇异值分解进行求解:
S w = U Σ V T (5) S_w=\rm{U}\Sigma V^{\rm{T}}\tag{5} Sw=UΣVT(5)
因此, S w − 1 = V Σ − 1 U T S_w^{-1}=\rm{V\Sigma^{-1}}U^{\rm{T}} Sw−1=VΣ−1UT,带入公式即可求解 w \bm{w} w.
2.1 多分类
当应用到多分类时,定义全局散度矩阵:
S t = S b + S w = ∑ i = 1 m ( x i − μ ) ( x i − μ ) T (6) S_t=S_b+S_w=\sum_{i=1}^m(\bm{x_i-\mu})(\bm{x_i}-\bm{\mu})^{\rm{T}}\tag{6} St=Sb+Sw=i=1∑m(xi−μ)(xi−μ)T(6)
μ \bm{\mu} μ 指的全部示例的均值向量
定义类内散度矩阵:
S w = ∑ i = 1 N S w i (7) S_w=\sum_{i=1}^NS_{w_i}\tag{7} Sw=i=1∑NSwi(7)
直接理解为全局的散度矩阵就是每一类的散度矩阵求和即可
定义全局类间散度矩阵:
S b = S t − S w (8) S_b=S_t-S_w \tag{8} Sb=St−Sw(8)
总体优化目标:
max W t r ( W T S w W ) t r ( W T S b W ) (9) \max_{\mathrm{W}}\frac{\mathrm{tr}(\mathrm{W}^{\rm{T}}S_w\mathrm{W})}{\mathrm{tr}(\mathrm{W}^{\rm{T}}S_b\mathrm{W})}\tag{9} Wmaxtr(WTSbW)tr(WTSwW)(9)
tr表示矩阵的迹
广义特征值问题求解:
S b W = λ S w W (10) S_b\mathrm{W} = \lambda S_w\mathrm{W}\tag{10} SbW=λSwW(10)
W \mathrm{W} W 的闭式解为 S w − 1 S b S_w^{-1}S_b Sw−1Sb 的 d ′ d^\prime d′ 个最大非零广义特征值对应的特征向量组成的矩阵, d ′ ≤ N − 1 d^{\prime}\le N-1 d′≤N−1代码实现就是使用直接的矩阵运算求解,构建投影矩阵,最后得到投影后的样本矩阵。注意:进行线性判别分析时,降低的维度最大是(类别数-1),如果类别数是10,样本数是100,则降低的最高维数是9.
**理解:**起初,可以理解二分类的线性判别分析,多分类难以理解,我选择这样理解,多分类就是多个w拼成一个W,求解的问题也变成了求W矩阵而不是求解w向量问题,从这里可以直接印证二分类问题将样本维度降低到1,就是样本矩阵乘以一个投影向量即得到投影后的维度为1,也就是将样本全部投影到一条直线上使得类间的距离尽可能大,使得类内距离尽可能小
3. 代码详解
算法流程:
- 计算全局散度矩阵 S t S_t St.
- 计算类内散度矩阵 S w S_w Sw.
- 计算类间散度矩阵 S b S_b Sb.
- 计算矩阵 S w − 1 S b S_w^{-1}S_b Sw−1Sb.
- 计算 S w − 1 S b S_w^{-1}S_b Sw−1Sb 的最大的d个特征值和对应的d个特征向量 ( w 1 , w 2 , ⋯ , w d ) (\bm{w}_1,\bm{w}_2,\cdots, \bm{w}_d) (w1,w2,⋯,wd).
- 对样本集中每个样本进行变换,构建新的样本集 Z = W T X Z=W^{\rm{T}}X Z=WTX.
朴素版本:
import numpy as np
from sklearn import datasets
# 准备数据
iris = datasets.load_iris()
X = iris.data
Y = iris.target
target_names = iris.target_names
# 1. 全局散度矩阵
mean = X.mean(axis=0)
St = np.dot((X - mean).T, X-mean)
# 2. 计算类内散度矩阵
Sw = np.zeros((X.shape[1], X.shape[1]))
for i in set(Y):
Sw += np.dot((X[Y==i, :] - X[Y==i, :].mean(axis=0)).T, X[Y==i, :] - X[Y==i, :].mean(axis=0))
# 3. 类间散度矩阵
Sb = St - Sw
# 4. Sw^-1Sb矩阵
temp = np.dot(mat(Sw).I, Sb)
# 5. 求特征值与特征向量
eigenvalue, featurevector = np.linalg.eig(temp)
# 6. 直接得到特征向量组成投影矩阵W进行样本矩阵投影变换
temp_index = np.argsort(eigenvalue) # 从小到大排序
W = featurevector[:, temp_index[-3:-1]] # 取最大的两个特征值对应的特征向量构成投影矩阵
trans_X = np.array(np.dot(X, W)) # 得到降维结果
高级版本:
from sklearn import datasets
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
# 准备数据
iris = datasets.load_iris()
X = iris.data
y = iris.target
target_names = iris.target_names
# 直接使用api
lda = LinearDiscriminantAnalysis(solver="eigen")
X_r2 = lda.fit(X, y).transform(X)
如果感兴趣可以查看sklearn源码,说明的非常详细,其中关于线性判别分析部分使用了三种方法,分别是特征值分解,svd求解,最小二乘求解。简单版本是参考西瓜书进行梳理产生,仅供参考。
参考链接:《机器学习(西瓜书)》