• 拉普拉斯特征映射(Laplacian Eigenmaps)


    1 介绍

      拉普拉斯特征映射(Laplacian Eigenmaps)是一种不太常见的降维算法,它看问题的角度和常见的降维算法不太相同,是从局部的角度去构建数据之间的关系。具体来讲,拉普拉斯特征映射是一种基于图的降维算法,它希望相互间有关系的点(在图中相连的点)在降维后的空间中尽可能的靠近,从而在降维后仍能保持原有的数据结构。

    2 推导

      拉普拉斯特征映射通过构建邻接矩阵为 W 的图来重构数据流形的局部结构特征。

      其主要思想是:如果两个数据实例 ij 很相似,那么 ij 在降维后目标子空间中应该尽量接近。设数据实例的数目为 n ,目标子空间即最终的降维目标的维度为 m 。 定义 n×m 大小的矩阵 Y ,其中每一个行向量 yTi 是数据实例 i 在目标 m 维子空间中的向量表示(即降维后的数据实例 i )。我们的目的是让相似的数据样例 ij 在降维后的目标子空间里仍旧尽量接近,故拉普拉斯特征映射优化的目标函数如下:

        mini,jyiyj2Wij

      下面开始推导:

        ni=1nj=1yiyj2Wij=ni=1nj=1(yTiyi2yTiyj+yTjyj)Wij=ni=1(nj=1Wij)yTiyi+nj=1(ni=1Wij)yTjyj2ni=1nj=1yTiyjWij=2ni=1DiiyTiyi2ni=1nj=1yTiyjWij=2ni=1(Diiyi)T(Diiyi)2ni=1yTi(nj=1yjWij)=2trace(YTDY)2ni=1yTi(YW)i=2trace(YTDY)2trace(YTWY)=2trace[YT(DW)Y]=2trace(YTLY)

      其中 W 是图的邻接矩阵,对角矩阵 D 是图的度矩阵 (Dii=nj=1Wij)L=DW 为图的拉普拉斯矩阵。

      变换后的拉普拉斯特征映射优化的目标函数如下:

        mintrace(YTLY) s.t. YTDY=I

      限制条件  s.t.YTDY=I 为了消除低维空间中的缩放因子,也为了保证 Dii 值较大的样本点在低维空间中更为重要,下面用拉格朗日乘子法对目标函数求解:

        f(Y)=tr(YTLY)+tr[Λ(YTDYI)]

        f(Y)Y=LY+LTY+DTYΛT+DYΛ=2LY+2DYΛ=0

        LY=DYΛ

      其中,Λ 为一个对角矩阵,另外 LD 均为实对称矩阵,其转置与自身相等。对于单独的 y 向量,上式可写为: Ly=λDy,这是一个广义特征值问题。通过求得 m 个最小非零特征值所对应的特征向量,即可达到降维的目 的。

      关于这里为什么要选择 m 个最小非零特征值所对应的特征向量。将 LY=DYΛ 带回到 mintrace(YTLY) 中,由于有着约束条件 YTDY=I 的限制,可以得到 mintrace(YTLY)=mintrace(Λ) 。即为特 征值之和。我们为了目标函数最小化,要选择最小的 m 个特征值所对应的特征向量。

    3 步骤

      使用时算法具体步骤为:

      步骤1:构建图

        使用某一种方法来将所有的点构建成一个图,例如使用KNN算法,将每个点最近的K个点连上边。K是一个预先设定的值。

      步骤2:确定权重

        确定点与点之间的权重大小,例如选用热核函数来确定,如果点 i 和点 j 相连,那么它们关系的权重设定为:

        Wij=exixj2t

        另外一种可选的简化设定是 Wij=1 如果点 ij 相连,否则 Wij=0

      步骤3:特征映射

        计算拉普拉斯矩阵 L 的特征向量与特征值: Ly=λDy

        使用最小的 m 个非零特征值对应的特征向量作为降维后的结果输出。

     

    Code:

    # coding:utf-8
    
    import numpy as np
    import matplotlib.pyplot as plt
    from sklearn.datasets import load_digits
    from mpl_toolkits.mplot3d import Axes3D
    
    def make_swiss_roll(n_samples=100, noise=0.0, random_state=None):
        #Generate a swiss roll dataset.
        t = 1.5 * np.pi * (1 + 2 * np.random.rand(1, n_samples))
        x = t * np.cos(t)
        y = 83 * np.random.rand(1, n_samples)
        z = t * np.sin(t)
        X = np.concatenate((x, y, z))
        X += noise * np.random.randn(3, n_samples)
        X = X.T
        t = np.squeeze(t)
        return X, t
    
    def rbf(dist, t = 1.0):
        '''
        rbf kernel function
        '''
        return np.exp(-(dist/t))
    
    def cal_pairwise_dist(x):
    
        '''计算pairwise 距离, x是matrix
        (a-b)^2 = a^2 + b^2 - 2*a*b
        '''
        sum_x = np.sum(np.square(x), 1)
        dist = np.add(np.add(-2 * np.dot(x, x.T), sum_x).T, sum_x)
        #返回任意两个点之间距离的平方
        return dist
    
    def cal_rbf_dist(data, n_neighbors = 10, t = 1):
    
        dist = cal_pairwise_dist(data)
        dist[dist < 0] = 0
        n = dist.shape[0]
        rbf_dist = rbf(dist, t)
    
        W = np.zeros((n, n))
        for i in range(n):
            index_ = np.argsort(dist[i])[1:1+n_neighbors]
            W[i, index_] = rbf_dist[i, index_]
            W[index_, i] = rbf_dist[index_, i]
    
        return W
    
    def le(data,
              n_dims = 2,
              n_neighbors = 5, t = 1.0):
        '''
        :param data: (n_samples, n_features)
        :param n_dims: target dim
        :param n_neighbors: k nearest neighbors
        :param t: a param for rbf
        :return:
        '''
        N = data.shape[0]
        W = cal_rbf_dist(data, n_neighbors, t)
        D = np.zeros_like(W)
        for i in range(N):
            D[i,i] = np.sum(W[i])
    
        D_inv = np.linalg.inv(D)
        L = D - W
        eig_val, eig_vec = np.linalg.eig(np.dot(D_inv, L))
    
        sort_index_ = np.argsort(eig_val)
    
        eig_val = eig_val[sort_index_]
        print("eig_val[:10]: ", eig_val[:10])
    
        j = 0
        while eig_val[j] < 1e-6:
            j+=1
    
        print("j: ", j)
    
        sort_index_ = sort_index_[j:j+n_dims]
        eig_val_picked = eig_val[j:j+n_dims]
        print(eig_val_picked)
        eig_vec_picked = eig_vec[:, sort_index_]
    
        # print("L: ")
        # print(np.dot(np.dot(eig_vec_picked.T, L), eig_vec_picked))
        # print("D: ")
        # D not equal I ???
        print(np.dot(np.dot(eig_vec_picked.T, D), eig_vec_picked))
    
        X_ndim = eig_vec_picked
        return X_ndim
    
    if __name__ == '__main__':
        X, Y = make_swiss_roll(n_samples = 2000)  #生成瑞士卷数据集
        X_ndim = le(X, n_neighbors = 5, t = 20)
        
        fig = plt.figure(figsize=(12,6))
        ax1 = fig.add_subplot(121, projection='3d')
        ax1.scatter(X[:, 0], X[:, 1], X[:, 2], c = Y)
        
        ax2 = fig.add_subplot(122)
        ax2.scatter(X_ndim[:, 0], X_ndim[:, 1], c = Y)
        plt.show()
    
        X = load_digits().data
        y = load_digits().target
    
        dist = cal_pairwise_dist(X)
        max_dist = np.max(dist)
        print("max_dist", max_dist)
        X_ndim = le(X, n_neighbors = 20, t = max_dist*0.1)
        plt.scatter(X_ndim[:, 0], X_ndim[:, 1], c = y)
        plt.savefig("LE2.png")
        plt.show()

     


    __EOF__

  • 本文作者: Blair
  • 本文链接: https://www.cnblogs.com/BlairGrowing/p/16029398.html
  • 关于博主: 评论和私信会在第一时间回复。或者直接私信我。
  • 版权声明: 本博客所有文章除特别声明外,均采用 BY-NC-SA 许可协议。转载请注明出处!
  • 声援博主: 如果您觉得文章对您有帮助,可以点击文章右下角推荐一下。
  • 相关阅读:
    table的展开折叠按钮操作
    Android CAN 简记
    JVM 别和我说你还不知道这几种垃圾回收器?Serial |Parallel|ParNew|CMS|G1|ZGC
    计算机毕业设计Java婚纱摄影网站(源码+系统+mysql数据库+lw文档)
    《自然语言处理实战:利用Python理解、分析和生成文本》读书笔记:第3章 词中的数学
    spring boot项目一次性能测试的总结
    Microsemi Libero SoC 教程2 (点亮LED闪烁)
    图纸管理制度《一》
    kubernetes code-generator使用
    1、web自动化环境安装(selenium安装、卸载和Chromedriver下载、安装)
  • 原文地址:https://www.cnblogs.com/BlairGrowing/p/16029398.html