特征抽取篇之LLE(Locally Linear Embedding)

前言

与等度量映射（Isomap）相同，局部线性嵌入（Locally Linear Embedding, 简称LLE）也是流形学习算法中的一员，是一种非线性的降维方法。LLE在2000年被发明出来，并发表在Sience期刊上，截止2019年5月，被引量高达1.3万。它与Isomap试图保持样本间的距离不同，LLE算法试图保持的是邻域内样本之间的线性关系。如下图所示：

如图，高维空间中样本 $x_i$ 可由它的近邻样本 $xj\boldsymbol x_j$ ， $xk\boldsymbol x_k$ ， $xl\boldsymbol x_l$ 通过线性组合而重构出来，即
$xi=wijxj+wikxk+wilxl(1)\boldsymbol x_i=w_{ij}\boldsymbol x_j+w_{ik}\boldsymbol x_k+w_{il}\boldsymbol x_l(1)$
LLE算法希望式(1)中的关系在低维空间中得以保持。

算法推导

LLE先为每个样本 $x_i$ 找近邻下标集合 $Q_i$ ，然后计算基于 $Q_i$ 中的样本点对 $x_i$ 进行线性重构的系数 $w_i$ ，优化目标为重构误差最小：
$min⁡W∑i=1m∥xi−∑j∈Qiwijxj∥22(2)s.t.∑j∈Qiwij=1\begin{aligned} \min\limits_{\boldsymbol W}&\sum^m_{i=1}\| \boldsymbol x_i-\sum_{j \in Q_i}w_{ij}\boldsymbol x_j \|^2_2(2)\\ s.t.&\sum_{j \in Q_i}w_{ij}=1 \end{aligned}$
转换为
$∑i=1m∥xi−∑j∈Qiwijxj∥22=∑i=1m∥∑j∈Qiwijxi−∑j∈Qiwijxj∥22=∑i=1m∥∑j∈Qiwij(xi−xj)∥22=∑i=1mWiT(xi−xj)(xi−xj)TWi=∑i=1mWiTCiWi\begin{aligned} \sum^m_{i=1}\| \boldsymbol x_i-\sum_{j \in Q_i}w_{ij}\boldsymbol x_j \|^2_2 &=\sum^m_{i=1}\| \sum_{j \in Q_i}w_{ij}\boldsymbol x_i- \sum_{j \in Q_i}w_{ij}\boldsymbol x_j \|^2_2 \\ &=\sum^m_{i=1}\| \sum_{j \in Q_i}w_{ij}(\boldsymbol x_i- \boldsymbol x_j) \|^2_2\\ &=\sum^m_{i=1}\boldsymbol W^T_i(\boldsymbol x_i-\boldsymbol x_j)(\boldsymbol x_i-\boldsymbol x_j)^T\boldsymbol W_i\\ &=\sum^m_{i=1}\boldsymbol W^T_i\boldsymbol C_i\boldsymbol W_i \end{aligned}$
其中， $Wi=(wi1,wi2,⋅⋅⋅,wik)T\boldsymbol W_i=(w_{i1},w_{i2},\cdot\cdot\cdot,w_{ik})^T$ ， $k$ 是 $Q_i$ 集合的长度， $Ci=(xi−xj)(xi−xj)T\boldsymbol C_i=(\boldsymbol x_i-\boldsymbol x_j)(\boldsymbol x_i-\boldsymbol x_j)^T$ ， $\in Q_i$ 。
$∑j∈Qiwij=WiT1k=1\sum_{j\in Q_i}w_{ij}=\boldsymbol W_i^T\boldsymbol 1_k=1$
其中， $1k\boldsymbol 1_k$ 为k维全1向量。
运用拉格朗日乘子法可得，
$J(W)=∑i=1mWiTCiWi+λ(WiT1k−1)∂J(W)∂Wi=2CiWi+λ1k\begin{aligned} J(\boldsymbol W)&=\sum^m_{i=1}\boldsymbol W^T_i\boldsymbol C_i\boldsymbol W_i+\lambda(\boldsymbol W_i^T\boldsymbol 1_k-1)\\ \cfrac{\partial J(\boldsymbol W)}{\partial \boldsymbol W_i} &=2\boldsymbol C_i\boldsymbol W_i+\lambda\boldsymbol 1_k \end{aligned}$
令 $∂J(W)∂Wi=0\cfrac{\partial J(\boldsymbol W)}{\partial \boldsymbol W_i}=0$ ，故
$Wi=−12λCi−11kWi=λCi−11k\begin{aligned} \boldsymbol W_i&=-\cfrac{1}{2}\lambda\boldsymbol C_i^{-1}\boldsymbol 1_k\\ \boldsymbol W_i&=\lambda\boldsymbol C_i^{-1}\boldsymbol 1_k\\ \end{aligned}$
其中，由于 $λ\lambda$ 原本为一个常数，加入 $−12这个放缩因子-\cfrac{1}{2}这个放缩因子$ 并无影响。利用 $WiT1k=1\boldsymbol W^T_i\boldsymbol 1_k=1$ ，对 $Wi\boldsymbol W_i$ 归一化，可得
$Wi=Ci−11k1kCi−11k\boldsymbol W_i=\cfrac{\boldsymbol C^{-1}_i\boldsymbol 1_k}{\boldsymbol 1_k\boldsymbol C^{-1}_i\boldsymbol 1_k}$
在得到高维空间样本点之间的关系之后，按照LLE的算法思想（在低维空间中保持 $w_i$ 不变），于是 $x_i$ 降维后的样本 $z_i$ 应达到的优化目标如下式(2)：
$min⁡Z∑i=1m∥zi−∑j∈Qiwijzj∥22(3)\min\limits_{\boldsymbol Z}\sum^m_{i=1}\| \boldsymbol z_i-\sum_{j \in Q_i}w_{ij}\boldsymbol z_j \|^2_2(3)$
与求解 $W$ 一样的过程：
$min⁡Z∑i=1m∥zi−∑j∈Qiwijzj∥22=∑i=1m∥ZIi−ZWi∥22=∑i=1m∥Z(Ii−Wi)∥22=∑i=1m(Z(Ii−Wi))TZ(Ii−Wi)=∑i=1m(Ii−Wi)TZTZ(Ii−Wi)=tr((I−W)TZTZ(I−W))=tr(Z(I−W)(I−W)TZT)=tr(ZMZT)\begin{aligned} \min\limits_{\boldsymbol Z}\sum^m_{i=1}\| \boldsymbol z_i-\sum_{j \in Q_i}w_{ij}\boldsymbol z_j \|^2_2&=\sum^m_{i=1}\|\boldsymbol Z\boldsymbol I_i-\boldsymbol Z\boldsymbol W_i\|^2_2\\ &=\sum^m_{i=1}\|\boldsymbol Z(\boldsymbol I_i-\boldsymbol W_i)\|^2_2\\ &=\sum^m_{i=1}(\boldsymbol Z(\boldsymbol I_i-\boldsymbol W_i))^T\boldsymbol Z(\boldsymbol I_i-\boldsymbol W_i)\\ &=\sum^m_{i=1}(\boldsymbol I_i-\boldsymbol W_i)^T\boldsymbol Z^T\boldsymbol Z(\boldsymbol I_i-\boldsymbol W_i)\\ &=tr((\boldsymbol I-\boldsymbol W)^T\boldsymbol Z^T\boldsymbol Z(\boldsymbol I-\boldsymbol W))\\ &=tr(\boldsymbol Z(\boldsymbol I-\boldsymbol W)(\boldsymbol I-\boldsymbol W)^T\boldsymbol Z^T)\\ &=tr(\boldsymbol Z\boldsymbol M\boldsymbol Z^T) \end{aligned}$
其中， $M=(I−W)(I−W)T\boldsymbol M=(\boldsymbol I-\boldsymbol W)(\boldsymbol I-\boldsymbol W)^T$ 。
再加上约束条件 $ZTZ=I\boldsymbol Z^T\boldsymbol Z=\boldsymbol I$ ，保证得到标准化（标准正交空间）的低维数据。
式(2)重写为：
$min⁡Ztr(ZMZT)(4)s.t.ZTZ=I\begin{aligned} &\min\limits_{\boldsymbol Z}tr(\boldsymbol Z\boldsymbol M\boldsymbol Z^T)(4) \\ &s.t.\boldsymbol Z^T\boldsymbol Z=\boldsymbol I \end{aligned}$
式(3)可以通过特征值分解求解： $M$ 最小的 $d^{'}$ 个特征值对应的特征向量组成的矩阵即为 $ZT\boldsymbol Z^T$ 。

算法步骤

输入：样本集 $D={x1,x2,⋅⋅⋅,xm}D=\{\boldsymbol x_1,\boldsymbol x_2,\cdot\cdot\cdot,\boldsymbol x_m\}$ ；近邻参数 $k$ ；低维空间维数 $d^{'}$ 。
步骤：
1：for i=1,2,…,m do
2：确定 $xi\boldsymbol x_i$ 的 $k$ 近邻；
3：从式(2)求得 $w_{ij}$ ， $\in Q_i$ ；
4：对于 $\notin Q_i$ ，令 $w_{ij}=0$ ；
5：end for
6：从式(3)得到 $M\boldsymbol M$ ；
7：对 $M\boldsymbol M$ 进行特征值分解；
8：return $M\boldsymbol M$ 的最小 $d^{'}$ 个特征值对应的特征向量
输出：样本集 $D$ 在低维空间的投影 $Z={z1,z2,⋅⋅⋅,zm}Z=\{ \boldsymbol z_1,\boldsymbol z_2,\cdot\cdot\cdot,\boldsymbol z_m \}$ 。

实验效果

调用sklearn库

"""
LLE : Locally Linear Embedding
Refercences :
[1]周志华.机器学习[M].清华大学出版社,2016:425.
[2]http://scikit-learn.sourceforge.net/dev/modules/generated/sklearn.manifold.LocallyLinearEmbedding.htmlAuthor : Ggmatch
Date : 2019/5/14
"""
from time import timeimport matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.ticker import NullFormatterfrom sklearn import manifold, datasets# 制造样本
n_points = 1000
X, color = datasets.samples_generator.make_s_curve(n_points, random_state=0)
n_neighbors = 10fig = plt.figure(figsize=(6, 4))  #画板
gs = fig.add_gridspec(1,2)  #共2副子图
ax1 = fig.add_subplot(gs[0,0], projection='3d')  #第一幅子图表示原始样本分布
ax1.scatter(X[:, 0], X[:, 1], X[:, 2], c=color, cmap=plt.cm.Spectral)# LLE降维
n_components = 2t0 = time()  #计时开始
Y = manifold.LocallyLinearEmbedding(n_neighbors, n_components).fit_transform(X)
t1 = time()  #计时结束
ax2 = fig.add_subplot(gs[0,1])
ax2.scatter(Y[:, 0], Y[:, 1], c=color, cmap=plt.cm.Spectral)  #第2副子图表示降维后样本分布
ax2.set_title("LLE (%.2g sec)" % (t1 - t0))
ax2.xaxis.set_major_formatter(NullFormatter())
ax2.yaxis.set_major_formatter(NullFormatter())plt.show()

实验效果

参考文献

[1]Roweis S T, Saul L K. Nonlinear dimensionality reduction by locally linear embedding[J]. science, 2000, 290(5500): 2323-2326.
[2]周志华.机器学习[M].清华大学出版社,2016:425.
[3]http://scikit-learn.sourceforge.net/dev/modules/generated/sklearn.manifold.LocallyLinearEmbedding.html