转载：https://www.cnblogs.com/pinard/p/6266408.html?utm_source=itdadao&utm_medium=referral#!comments
局部线性嵌入(Locally Linear Embedding) 是非常重要的姜维方法。与传统的PCA，LDA等关注样本方差的降维方法相比，LLE关注降维时保持样本局部的线性特征，由于这个特性，它广泛运用于图像识别，高维数据可视化等领域。
LLE属于流形学习(Manifold Learning)的一种。流形学习是一大类基础流形的框架。数学意义上的流形比较抽象，不过我们可以认为LLE中的流形是一个不闭合的曲面。这个流形曲面有数据分布比较均匀，且比较稠密的特征。基于流行的降维算法就是将流形从高维到低维的降维过程，在降维的过程中我们希望流形在高维的一些特征可以得到保留。
一个形象的流形降维过程如下图。我们有一块卷起来的布，我们希望将其展开到一个二维平面，我们希望展开后的布能够在局部保持布结构的特征，其实也就是将其展开的过程，就想两个人将其拉开一样。
在局部保持布结构的特征，或者说数据特征的方法有很多种，不同的保持方法对应不同的流形算法。比如等距映射(ISOMAP)算法在降维后希望保持样本之间的测地距离而不是欧式距离，因为测地距离更能反映样本之间在流形中的真实距离。 但是等距映射算法的一个问题是要找所有样本全局的最优解，当数据量很大，样本维度很高时，计算非常耗时，鉴于这个问题，LLE通过放弃所有样本全局最优的降维，只是通过保证局部最优来降维。同时假设样本集在局部是满足线性关系的，进一步减少的降维的计算量。LLE首先假设数据在较小的局部是线性的，也就是说，某一个数据可以由它领域中的几个样本点来线性表示。比如：一个样本 x 1 x_1 x1，我们在它原始高维空间中用K-近邻思想找到和它最近的三个点 x 2 , x 3 , x 4 . x_2,x_3,x_4. x2,x3,x4. 然后假设 x 1 x_1 x1 可以由 x 2 , x 3 , x 4 x_2,x_3,x_4 x2,x3,x4 线性表示， 即：
x 1 = w 12 x 2 + w 13 x 3 + w 14 x 4 x_1 = w_{12}x_2 + w_{13}x_3 + w_{14}x_4 x1=w12x2+w13x3+w14x4 其中， w 12 , w 13 , w 14 w_{12},w_{13},w_{14} w12,w13,w14 为权重系数。在通过LLE降维后，我们希望 x 1 x_1 x1 在低维空间对应的投影 x 1 ′ x'_1 x1′ 和 x 2 , x 3 , x 4 x_2,x_3,x_4 x2,x3,x4 对应的投影 x 2 ′ , x 3 ′ , x 4 ′ x'_2,x'_3,x'_4 x2′,x3′,x4′ 线性表示，即： x 1 ′ ≈ w 12 x 2 ′ + w 13 ′ x 3 ′ + w 14 x 4 ′ x'_1 \approx w_{12}x'_2 + w'_{13}x'_3 + w_{14}x'_4 x1′≈w12x2′+w13′x3′+w14x4′ 也就是说，投影前后线性关系的权重系数 w 12 , w 13 , w 14 w_{12},w_{13},w_{14} w12,w13,w14 是尽量不变或者改变最小的。从上可知，线性关系只在样本的附近起作用，离样本远的样本对局部的线性关系没有影响，因此降维的复杂度降低了很多。
对于LLE算法，我们首先要确定邻域大小的选择，即我们需要多少个邻域样本来线性表示某个样本。假设这个值为k。我们可以通过和KNN一样的思想通过距离度量比如欧式距离来选择某样本的k个最近邻。
在寻找到某个样本 x i x_i xi 的k个最近邻之后我们就需要找到 x i x_i xi 和这k个邻居之间的线性关系，也就是要找到线性关系的权重系数。找线性关系，这显然是一个回归问题。假设我们有m个n维样本 { x 1 , x 2 , ⋅ ⋅ ⋅ , x m } \left \{ x_1,x_2,···,x_m \right \} {x1,x2,⋅⋅⋅,xm},我们可以用均方差作为回归问题的损失函数：即： J ( w ) = ∑ i = 1 m ∥ x i − ∑ j ∈ Q ( i ) w i j x j ∥ 2 2 J(w)=\sum_{i=1}^{m}\left \|x_i- \sum_{j\in Q(i)}w_{ij}x_j \right \|^{2}_2 J(w)=i=1∑m∥∥∥∥∥∥xi−j∈Q(i)∑wijxj∥∥∥∥∥∥22 其中 Q ( i ) Q(i) Q(i) 表示 i i i 的K个近邻样本集合。一般我们也会对权重系数 w i j w_{ij} wij做归一化的限制，即权重系数需要满足
∑ j ∈ Q ( i ) w i j = 1 \sum_{j\in Q(i)} w_{ij} = 1 j∈Q(i)∑wij=1 对于不在样本 x i x_i xi 邻域内的样本 x j x_j xj，我们令对应的 w i j = 0 w_{ij} =0 wij=0 ，这样可以把 w w w扩展到整个数据集的维度。 也就是我们需要通过上面两个式子求出我们的权重系数。一般我们可以通过矩阵和拉格朗日乘子法来求解这个最优化问题。对于第一个式子，我们先将其矩阵化： J ( w ) = ∑ i = 1 m ∥ x i − ∑ j ∈ Q ( i ) w i j x j ∥ 2 2 = ∑ i = 1 m ∥ ∑ j ∈ Q ( i ) w i j x i − ∑ j ∈ Q ( i ) w i j x j ∥ 2 2 = ∑ i = 1 m ∥ ∑ j ∈ Q ( i ) w i j ( x i − x j ) ∥ 2 2 = ∑ i = 1 m W i T ( x i − x j ) ( x i − x j ) T W i J(w)=\sum_{i=1}^{m}\left \|x_i- \sum_{j\in Q(i)}w_{ij}x_j \right \|^{2}_2\\ = \sum_{i=1}^{m}\left \| \sum_{j\in Q(i)}w_{ij}x_i-\sum_{j\in Q(i)}w_{ij}x_j \right \|^2_2\\ = \sum_{i=1}^{m}\left \| \sum_{j\in Q(i)}w_{ij}(x_i-x_j) \right \|^2_2\\ = \sum_{i=1}^m W_i^T(x_i-x_j)(x_i-x_j)^TW_i J(w)=i=1∑m∥∥∥∥∥∥xi−j∈Q(i)∑wijxj∥∥∥∥∥∥22=i=1∑m∥∥∥∥∥∥j∈Q(i)∑wijxi−j∈Q(i)∑wijxj∥∥∥∥∥∥22=i=1∑m∥∥∥∥∥∥j∈Q(i)∑wij(xi−xj)∥∥∥∥∥∥22=i=1∑mWiT(xi−xj)(xi−xj)TWi 的编辑器似乎不支持LaTeX的换行对齐，勉强看一下吧
其中 W i = ( w i 1 , w i 2 , w i 3 , ⋅ ⋅ ⋅ w i k ) T W_i=(w_{i1},w_{i2},w_{i3},···w_{ik})^T Wi=(wi1,wi2,wi3,⋅⋅⋅wik)T。
我们令矩阵 Z i = ( x i − x j ) ( x i − x j ) T , j ∈ Q ( i ) Z_i = (x_i-x_j)(x_i-x_j)^T,j\in Q(i) Zi=(xi−xj)(xi−xj)T,j∈Q(i)则上式进一步简化为 J ( W ) = ∑ i = 1 k W i T Z i W i J(W) = \sum_{i=1}^k W_i^TZ_iW_i J(W)=i=1∑kWiTZiWi 对于权重系数的归一化限制，可以写成： ∑ j ∈ Q ( i ) w i j = W I T 1 k = 1 \sum_{j\in Q(i)}w_{ij} = W_I^T 1_k=1 j∈Q(i)∑wij=WIT1k=1其中 1 k 1_k 1k为k维全1向量。 现在将以上两个矩阵化之后的式子用拉格朗日乘子法合成为一个优化目标： L ( W ) = ∑ i = 1 K W i T Z i W i + λ ( W i T 1 k − 1 ) L(W) = \sum_{i=1}^KW_i^TZ_iW_i + \lambda (W_i^T1_k-1) L(W)=i=1∑KWiTZiWi+λ(WiT1k−1)对 W W W 进行求导并令其值为0，得到 2 Z i W i + λ 1 k = 0 2Z_iW_i+\lambda 1_k =0 2ZiWi+λ1k=0即有: W i = λ ′ Z i − 1 1 k W_i=\lambda'Z_i^{-1}1_k Wi=λ′Zi−11k其中 λ ′ = − 1 2 λ \lambda' =-\frac{1}{2}\lambda λ′=−21λ 为一个常数。利用 W i T 1 k = 1 W_i^T1_k=1 WiT1k=1，对 W i W_i Wi归一化，那么最终我们的权重系数 W i W_i Wi 为： W i = Z i − 1 1 k 1 k T Z i − 1 1 k W_i = \frac{Z_i^{-1}1_k}{1_k^TZ_i^{-1}1_k} Wi=1kTZi−11kZi−11k 现在我们得到了高维的权重系数，那么我们希望这些权重系数对应的线性关系在降维后的低维一样得到保持。假设我们的n维样本集 { x 1 , x 2 , ⋅ ⋅ ⋅ , x m } \left \{ x_1,x_2,···,x_m \right \} {x1,x2,⋅⋅⋅,xm}在低维的d维度对应投影为 { y 1 , y 2 , ⋅ ⋅ ⋅ , y m } \left \{ y_1,y_2,···,y_m \right \} {y1,y2,⋅⋅⋅,ym}, 则我们希望保持线性关系，也就是希望对应的均方差损失函数最小，即最小化损失函数 J ( Y ) J(Y) J(Y)如下： J ( y ) = ∑ i = 1 m ∥ y i − ∑ j = 1 m w i j y j ∥ 2 2 J(y) = \sum_{i=1}^m \left \| y_i- \sum_{j=1}^mw_{ij}y_j \right \|^{2}_2 J(y)=i=1∑m∥∥∥∥∥yi−j=1∑mwijyj∥∥∥∥∥22可以看到这个式子和高维损失函数几乎相同。在高维时，高维数据已知，目标是最小值对应的权重系数 W W W ，而在低维是权重系数 W W W 已知，求对应的低维数据。注意，这里的 W W W已经是 m × m m×m m×m 维度，之前的 W W W是 m × k m×k m×k 维度，我们将那些不在邻域位置的W的位置取值为0，将 W W W 扩充到 m × m m×m m×m 维度。
为了得到标准化的低维数据，一般我们也加入约束条件如下: ∑ i = 1 m y i = 0 ; 1 m ∑ i = 1 m y i y i T = I \sum_{i=1}^my_i=0; \frac{1}{m}\sum_{i=1}^my_iy_i^T = I i=1∑myi=0;m1i=1∑myiyiT=I首先我们将目标损失函数矩阵化： J ( y ) = ∑ i = 1 m ∥ y i − ∑ j = 1 m w i j y j ∥ 2 2 = ∑ i = 1 m ∥ Y I i − Y W i ∥ 2 2 = t r ( Y ( I − W ) ( I − W ) T Y T ) J(y) = \sum_{i=1}^m \left \| y_i- \sum_{j=1}^mw_{ij}y_j \right \|^{2}_2\\=\sum_{i=1}^m\left \| YI_i-YW_i \right \|^{2}_2\\=tr(Y(I-W)(I-W)^TY^T) J(y)=i=1∑m∥∥∥∥∥yi−j=1∑mwijyj∥∥∥∥∥22=i=1∑m∥YIi−YWi∥22=tr(Y(I−W)(I−W)TYT)如果我们令 M = ( I − W ) ( I − W ) T M=(I-W)(I-W)^T M=(I−W)(I−W)T，则优化函数可写成 J ( Y ) = t r ( Y M Y T ) J(Y)=tr(YMY^T) J(Y)=tr(YMYT)，tr为矩阵求迹。约束函数矩阵化为： Y Y T = m I YY^T=mI YYT=mI
如果熟悉谱聚类和PCA，会发现这里的优化过程几乎一样。其实最小化 J ( Y ) J(Y) J(Y) 对应的 Y Y Y 就是 M M M 的最小的d个特征值所对应的d个特征向量组成的矩阵。当然我们同样可以通过拉格朗日乘子法来求解函数： L ( Y ) = t r ( Y M Y T + λ ( Y Y T − m I ) ) L(Y) = tr(YMY^T+\lambda(YY^T-mI)) L(Y)=tr(YMYT+λ(YYT−mI))上式对 Y Y Y求导另其为0，得到 2 M Y T + 2 λ Y T = 0 2MY^T + 2\lambda Y^T=0 2MYT+2λYT=0，即 M Y T = λ Y T MY^T=\lambda Y^T MYT=λYT，这样我们就清楚要得到最小的d维数据集，需要求出矩阵 M M M 最小的d个特征值所对应的d个特征向量组成的矩阵 Y = ( y 1 , y 2 , ⋅ ⋅ ⋅ , y d ) Y = (y_1,y_2,···,y_{d}) Y=(y1,y2,⋅⋅⋅,yd) 即可。
一般的，由于M的最小特征值为0不能反应数据特征，此时对应的特征向量全为1.我们通常选择M的第2到d+1个最小特征值对应的特征向量 Y = ( y 2 , y 3 , ⋅ ⋅ ⋅ , y d + 1 ) Y = (y_2,y_3,···,y_{d+1}) Y=(y2,y3,⋅⋅⋅,yd+1) 来得到最终的Y。为什么M的最小特征值为0呢？这是因为 W T e = e W^Te = e WTe=e, 得到 ∣ W T − I ∣ e = 0 |W^T−I|e=0 ∣WT−I∣e=0,由于 e ≠ 0 e \neq 0 e̸=0,所以只有 W T − I = 0 W^T−I=0 WT−I=0,即 ( I − W ) T = 0 (I−W)^T=0 (I−W)T=0,两边同时左乘 I − W I−W I−W ，即可得到 ( I − W ) ( I − W ) T e = 0 e (I−W)(I−W)^Te=0e (I−W)(I−W)Te=0e,即M的最小特征值为0.
在上一节我们已经基本推导了LLE降维的整个流程，现在我们对算法过程做一个总结。整个LLE算法用一张图可以表示如下：
从图中可以看出，LLE算法主要分为三步，第一步是求K近邻的过程，这个过程使用了和KNN算法一样的求最近邻的方法。第二步，就是对每个样本求它在邻域里的K个近邻的线性关系，得到线性关系权重系数W，具体过程在第三节第一部分。第三步就是利用权重系数来在低维里重构样本数据，具体过程在第三节第二部分。具体过程如下：
输入：样本集 D = { x 1 , x 2 , ⋅ ⋅ ⋅ , x m } D=\left \{ x_1,x_2,···,x_m \right \} D={x1,x2,⋅⋅⋅,xm} 最近邻数k，降维到的维数d输出： 低维样本集矩阵 D ′ D' D′ for i 1 to m， 按欧式距离作为度量，计算和xi最近的的k个最近邻 ( x i 1 , x i 2 , ⋅ ⋅ ⋅ , x i k ) ( x_{i1},x_{i2},···,x_{ik} ) (xi1,xi2,⋅⋅⋅,xik)for i 1 to m，求出 局部协方差矩阵 Z i = ( x i − x j ) T ( x i − x j ) Z_i=(x_i-x_j)^T(x_i-x_j) Zi=(xi−xj)T(xi−xj)，并求出对应的权重系数向量: W i = Z i − 1 1 k 1 k T Z i − 1 1 k W_i = \frac{Z_i^{-1}1_k}{1_k^TZ_i^{-1}1_k} Wi=1kTZi−11kZi−11k有权重系数向量 W i W_i Wi 组成权重系数矩阵 W W W，计算矩阵 M = ( I − W ) ( I − W ) T M=(I-W)(I-W)^T M=(I−W)(I−W)T计算矩阵M的前d+1个特征值，并计算这d+1个特征值对应的特征向量 { y 1 , y 2 , ⋅ ⋅ ⋅ y d + 1 } \left \{y_1,y_2, ···y_d+1 \right \} {y1,y2,⋅⋅⋅yd+1}由第二个特征向量到第d+1个特征向量所张成的矩阵即为输出低维样本集矩阵 D ′ = ( y 2 , y 3 , ⋅ ⋅ ⋅ y d + 1 ) D′=(y_2,y_3, ···y_d+1) D′=(y2,y3,⋅⋅⋅yd+1)LLE算法很简单高效，但是却有一些问题，比如如果近邻数k大于输入数据的维度时，我们的权重系数矩阵不是满秩的。为了解决这样类似的问题，有一些LLE的变种产生出来。比如：Modified Locally Linear Embedding(MLLE)和Hessian Based LLE(HLLE)。对于HLLE，它不是考虑保持局部的线性关系，而是保持局部的Hessian矩阵的二次型的关系。而对于MLLE，它对搜索到的最近邻的权重进行了度量，我们一般都是找距离最近的k个最近邻就可以了，而MLLE在找距离最近的k个最近邻的同时要考虑近邻的分布权重，它希望找到的近邻的分布权重尽量在样本的各个方向，而不是集中在一侧。
另一个比较好的LLE的变种是Local tangent space alignment(LTSA)，它希望保持数据集局部的几何关系，在降维后希望局部的几何关系得以保持，同时利用了局部几何到整体性质过渡的技巧。
这些算法原理都是基于LLE，基本都是在LLE这三步过程中寻求优化的方法。具体这里就不多讲了。
LLE是广泛使用的图形图像降维方法，它实现简单，但是对数据的流形分布特征有严格的要求。比如不能是闭合流形，不能是稀疏的数据集，不能是分布不均匀的数据集等等，这限制了它的应用。下面总结下LLE算法的优缺点。 LLE算法的主要优点有：
可以学习任意维的局部线性的低维流形算法归结为稀疏矩阵特征分解，计算复杂度相对较小，实现容易。LLE算法的主要缺点有：
算法所学习的流形只能是不闭合的，且样本集是稠密均匀的。算法对最近邻样本数的选择敏感，不同的最近邻数对最后的降维结果有很大影响。