尧图网站设计 尧图网站设计YAOTU DESIGN
ARTICLE DETAIL

资讯详情

深耕网站设计与一线实操的经验洞察。

等度量映射Isomap详解:原理、手写实现与调参全攻略

等度量映射Isomap详解:原理、手写实现与调参全攻略 不少朋友刚开始接触流形学习时第一个被问到“能不能讲讲等度量映射”的项目基本都是同一个用 Isomap 对瑞士卷数据做降维可视化。也确实Isomap 几乎成了继 PCA 之后最容易上手的非线性降维算法很多高校的机器学习作业、期末复习、还有像头歌这样的实训平台都会把它拎出来单独做一关。我之前给学生讲这块时发现大部分人对 Isomap 的理解停留在“把 PCA 换成测地距离”但一问到为什么、怎么算、参数怎么调就卡住了。所以这篇就专门写透等度量映射Isometric Feature Mapping简称 Isomap从它出现的动机到算法每一步的数学直觉再到手写代码和 sklearn 对比调参最后把实训平台常见的坑也一并列出来。无论你是应付期末、做课程设计还是单纯想把流形降维搞清楚这篇都值得存一份。1. 为什么线性降维不够用才有了等度量映射1.1 PCA 在瑞士卷上的“翻车现场”先说 PCA。它的目标是在原始高维空间中找一组正交方向使得数据投影到这些方向上后方差最大。本质上 PCA 假设数据是分布在一个线性子空间里的所以它对“平直”的数据效果很好比如人脸灰度图、基因表达矩阵这类本身就近似线性的数据。但现实世界的数据经常是“卷起来”的——一个典型的例子就是瑞士卷Swiss Roll数据集。三维空间里有一块平面数据被卷成了螺旋状从三维视角看两个在平面上很接近的点在欧氏距离下可能隔着很远的直线距离。比如瑞士卷的首尾两端展开后其实是相邻的但在三维空间里直接算欧氏距离却非常远。PCA 这时候就尴尬了。它只会老老实实地把数据投影到方差最大的几个方向上结果就是瑞士卷被“拍扁”成一张二维图本应该相邻的区域被远远拉开流形结构完全丢掉。用大白话说PCA 只会用“直线距离”看世界瑞士卷被卷起来之后直线距离根本不代表真实的远近关系。Isomap 要解决的就是这个问题。1.2 等度量映射的核心思想换个距离度量再看名字“等度量映射”翻译自 Isometric Feature Mapping英文里 isometric 就是“等距的”。算法的核心思路非常朴素在高维空间中数据实际上躺在一个低维流形比如瑞士卷的二维曲面上在这个流形上两点之间的真实距离应该是沿着流形表面的测地线距离geodesic distance而不是高维空间里的直线欧氏距离如果能近似算出测地线距离再喂给经典的多维缩放MDS去做降维就能还原出流形展开后的真实结构。所以 Isomap 做的事情本质上就是两步叠加先用近邻图逼近测地线距离再用 MDS 保持这个距离关系去降维。这也正是它和 PCA 的关键区别——PCA 保持的是方差最大Isomap 保持的是数据点之间的“流形距离”不变。这个看似很小的改动让它在处理强非线性结构时脱胎换骨。1.3 Isomap 在算法谱系里的位置在了解细节前先给 Isomap 定个位从降维算法的大家族看它属于流形学习Manifold Learning里的代表作和 LLE局部线性嵌入、t-SNE、UMAP 并称四大金刚从“保距”这个流派看它的目标是降维后尽量保持两点之间的测地线距离这也是它与保持局部结构的 LLE、保持概率分布的 t-SNE 最大不同。因为保留了全局的测地距离信息Isomap 对“展开型”流形比如瑞士卷、S 曲面效果极好但对“带孔洞、多连通分量”的数据就比较头疼这一点后面调参部分会细说。2. 等度量映射算法拆解三步走完Isomap 的完整流程可以拆成三大步构建近邻图、计算最短路径逼近测地线距离、用 MDS 降维。前两步是核心创新第三步是复用经典的 MDS。2.1 第一步构建近邻图先用一个无向图来表示数据点之间的局部连接关系做法有两种K 近邻法KNN每个点和最近的 K 个点连边ε 邻域法和距离小于 ε 的所有点连边。连边之后边的权重就是两点在原始高维空间里的欧氏距离。这一步的意义在于把连续流形离散成一个离散图用“近邻边”来模拟流形表面的局部平直区域。为什么用局部欧氏距离可以逼近测地线这里有一个流形学习很重要的工作假设局部上流形是“平坦”的可以用欧氏距离近似测地线。就像地球表面虽然是个球面但你量你家到小区门口的距离直接用直线量也没什么问题。只要近邻范围足够小欧氏距离就是测地线距离的合理局部近似。2.2 第二步最短路径逼近测地线距离近邻图建好了接下来计算任意两点之间的最短路径长度。这一步通常用 Dijkstra 算法或 Floyd-Warshall 算法实现。得到的任意两点之间的最短路径长度就是测地线距离的离散近似。用最短路径来代替真实测地线距离道理也很直观在图结构里从一个点走到另一个点最省力的那条路往往就是沿着流形表面绕过去的那条路——因为近邻边本身是沿着流形局部的方向延伸的多条局部边串起来就逼近了流形表面的整体走势。注意这里的近似精度直接取决于近邻图的质量。如果 K 选得太小图就断成几块距离算不出来如果 K 选得太大会出现“短路”现象本来在流形上隔很远的两个点被一条跨空隙的边直接连起来测地线距离就被严重低估。后面专门讲 K 怎么调。2.3 第三步MDS 多维缩放拿到所有点两两之间的测地线距离矩阵之后降维就变成了一个经典问题找一个低维坐标让低维空间中两两之间的欧氏距离尽量逼近这个测地线距离矩阵。这正好是 MDS 的主场。MDS 的思路是把距离矩阵 D 转换为内积矩阵 B对 B 做特征分解取前 d 个最大特征值对应的特征向量乘以特征值的平方根就得到 d 维坐标。推导过程我用白话重述一遍。距离矩阵里的元素是 (d_{ij}^2 ||x_i - x_j||^2)。在中心化假设下所有点的坐标均值为零内积 (b_{ij} x_i^T x_j) 和距离之间有关系[ b_{ij} -\frac{1}{2}\left(d_{ij}^2 - \frac{1}{n}\sum_k d_{ik}^2 - \frac{1}{n}\sum_k d_{kj}^2 \frac{1}{n^2}\sum_{k,l} d_{kl}^2\right) ]这其实就是对距离矩阵做双中心化double centering变成内积矩阵 B。然后对 B 做特征分解 (B V \Lambda V^T)取前 d 个特征值 (\lambda_1 \ge \lambda_2 \ge \dots \ge \lambda_d)要求大于等于零和对应的特征向量 (v_1, \dots, v_d)最终低维坐标为[ Y_d V_d \Lambda_d^{1/2} ]也就是说每一维等于对应特征向量乘以特征值的平方根。在 sklearn 的 Isomap 实现里这段逻辑封装得很彻底你直接传 n_components 就行但理解这层计算过程对调参和排查结果非常有帮助。3. 手写 Isomap从零实现到与 sklearn 结果对比自己实现一遍 Isomap胜过反复背十遍算法流程。这里我用 Python 的 NumPy 从零实现核心步骤再和 sklearn 的结果做对比验证实现是否正确。3.1 生成数据瑞士卷先构造经典的瑞士卷数据。用 sklearn 自带的 datasets 方法就行import numpy as np import matplotlib.pyplot as plt from sklearn.datasets import make_swiss_roll from sklearn.manifold import Isomap from sklearn.decomposition import PCA from scipy.spatial import distance_matrix from scipy.sparse.csgraph import shortest_path # 生成瑞士卷样本点 2000 个噪声小一点 X, color make_swiss_roll(n_samples2000, noise0.05) print(原始数据形状, X.shape) # 原始数据形状 (2000, 3)这里的 color 是每个点沿瑞士卷展开方向的“真实位置”后面可视化时用它给点着色方便观察降维后流形是否被正确展开。3.2 手写 Isomap 核心代码按照算法三步走代码可以写得非常干净def my_isomap(X, n_neighbors8, n_components2): # 1. 计算欧氏距离矩阵 D distance_matrix(X, X) # 2. 构建 KNN 近邻图并只保留近邻边的权重 n_samples X.shape[0] knn_graph np.full_like(D, np.inf) for i in range(n_samples): # 不包括自己所以从第 2 个小的开始取第一个是自身距离 0 neighbors np.argsort(D[i])[1:n_neighbors1] knn_graph[i, neighbors] D[i, neighbors] # 让图对称无向图 knn_graph np.minimum(knn_graph, knn_graph.T) np.fill_diagonal(knn_graph, 0) # 3. 用最短路径逼近测地线距离 # scipy 封装了 Dijkstra可以直接用 geodesic_dist shortest_path(knn_graph, methodD, directedFalse) # 4. 双中心化转成内积矩阵 B n n_samples G geodesic_dist ** 2 # 手动实现双中心化H I - (1/n) * 1*1^T H np.eye(n) - np.ones((n, n)) / n B -0.5 * H G H # 5. 特征分解取前 n_components 个特征值最大的特征向量 # 注意特征分解得到正负特征值只选正的 eigvals, eigvecs np.linalg.eigh(B) idx np.argsort(eigvals)[::-1][:n_components] eigvals_sorted eigvals[idx] eigvecs_sorted eigvecs[:, idx] # 低维坐标 特征向量 * 特征值平方根 Y eigvecs_sorted * np.sqrt(eigvals_sorted) return Y, geodesic_dist Y_mine, _ my_isomap(X, n_neighbors8, n_components2) print(手写 Isomap 输出形状, Y_mine.shape)这里有三个容易踩的细节单独说明一下近邻图必须对称化。距离矩阵本身是对称的但 KNN 的“A 是 B 的近邻B 不一定是 A 的近邻”这种情况很常见。所以要用 np.minimum(knn_graph, knn_graph.T) 把不对称的地方做对称化处理否则后续求最短路径会出问题。特征分解结果要小心处理。np.linalg.eigh 返回的特征值是从小到大排列的所以要用 argsort 先倒序再取前 n_components 个。而且 Isomap 理论上要求特征值非负实测中总能遇到小幅负特征值直接忽略即可。最后一行的 Y eigvecs_sorted * np.sqrt(eigvals_sorted) 是逐元素乘利用了广播机制让每个特征向量分别乘以对应特征值的平方根。不要写成矩阵乘法。3.3 对比 sklearn 官方实现sklearn 里 Isomap 的接口非常简洁iso Isomap(n_neighbors8, n_components2) Y_sklearn iso.fit_transform(X) print(sklearn Isomap 输出形状, Y_sklearn.shape)把两个结果放在一起对比。得益于 scipy 的最短路径优化sklearn 底层实际上也是用 Floyd 或 Dijkstra 类算法所以两个结果在数值上可能有一些符号差异但整体结构应当高度一致from scipy.stats import pearsonr # 因为特征向量方向可能相反先做符号对齐再算相关性 corr pearsonr(Y_mine[:, 0], Y_sklearn[:, 0])[0] print(第一维相关系数, corr)如果相关系数接近 1或者接近 -1因为特征向量方向翻转很正常就说明手写实现基本正确。我自己实测下来瑞士卷数据上第一维相关系数一般在 0.99 以上差一点点是数值精度导致的不影响结果。需要注意scipy 的 shortest_path 对于不连通的近邻图会返回 inf如果数据里有 inf后续的双中心化和特征分解直接报错。所以实际使用中先检查是否连通是一个好习惯# 检查测地距离矩阵是否包含无穷大 print(是否存在不连通点, np.isinf(geodesic_dist).any())如果返回 True说明 K 太小了近邻图没有把所有点连成一片需要增大 K 或者用 ε 邻域法。3.4 可视化对比PCA vs Isomap 的展开效果这一步是为了直观感受 Isomap 到底做了什么。把 PCA 投影结果、Isomap 降维结果分别画出来用 color 作为颜色映射fig, axes plt.subplots(1, 3, figsize(15, 4)) # 原始三维数据 ax axes[0] ax.scatter(X[:, 0], X[:, 1], X[:, 2], ccolor, cmaprainbow, s5) ax.set_title(3D Swiss Roll) # PCA 降维到 2D pca PCA(n_components2) Y_pca pca.fit_transform(X) ax axes[1] ax.scatter(Y_pca[:, 0], Y_pca[:, 1], ccolor, cmaprainbow, s5) ax.set_title(PCA to 2D) # Isomap 降维到 2D ax axes[2] ax.scatter(Y_sklearn[:, 0], Y_sklearn[:, 1], ccolor, cmaprainbow, s5) ax.set_title(Isomap to 2D) plt.tight_layout() plt.show()跑出来的效果PCA 那张图上彩虹色是乱成一团、互相重叠的说明瑞士卷展开后颜色区域并没有铺成一个平滑渐变而 Isomap 那张图会看到一个漂亮的彩虹条带从红到紫依次展开基本等价于把瑞士卷从三维里“扯平”成一张二维平面。这就是 Isomap 保持测地线距离的直观效果。4. 关键参数与调参心得Isomap 参数不多核心就两个n_neighbors或者 ε和 n_components。调参决定了结果好坏我逐个讲。4.1 n_neighbors近邻数量的选法n_neighbors 是 Isomap 最重要的超参数。它控制近邻图的局部范围直接影响测地线距离的精度。选太小图可能不连通选太大会出现“短路”。我在瑞士卷数据上做过一组对照实验这里把典型情况列出来n_neighbors现象原因3图容易断成多个连通分量降维结果一团糟近邻太少流形表面没有被完整覆盖局部信息不够8~12效果最好瑞士卷完美展开近邻数适中既覆盖了流形形状又没有跨过流形空隙30~50展开效果尚可但两端容易粘连边缘扭曲近邻数偏大部分跨空隙的边引入了“短路”误差100结果退化接近 PCA 效果流形结构被压扁近邻几乎覆盖全图测地线距离退化为全局欧氏距离选 n_neighbors 的经验法则可以总结成一句话在保证图连通的前提下尽量选小。实际操作中可以从小到大试几个值每次画一下降维结果看流形是否展开得干净利落而不是只盯着某一个指标死磕。4.2 欧氏距离在高维空间中的“诅咒”如果数据是高维的比如图像特征有 1000 维欧氏距离本身的可区分性会下降近邻关系也容易失真。这时候直接套 Isomap效果往往不好。常见的做法是先做 PCA 降维到几十维再跑 Isomap或者换用基于角度的距离度量。这一点很多人会忽略但它确实是 Isomap 在高维稀疏数据上效果翻车的首要原因——不是算法不对是距离度量已经失效了。4.3 n_components降维到几维Isomap 理论上是先算测地线距离矩阵再交给 MDS 处理所以 n_components 只要小于样本量即可。实际使用时需要注意如果目标是可视化n_components 通常设为 2 或 3如果目标是做特征提取喂给下游分类器n_components 需要根据流形的本征维数来估。比如瑞士卷的本征维数是 2它本质是一张平面卷起来的降维到 2 就是最优在实训平台或期末项目里如果不确定本征维数可以看特征值衰减的拐点。具体做法是画出前 M 个特征值的累计贡献率拐点处对应的维数就是合适的目标维数。4.4 数据规模与计算效率Isomap 的计算瓶颈在两步一是所有点对最短路径二是对 n×n 的内积矩阵做特征分解。当样本量到 1 万以上时N 平方量级的计算和存储都很吃力跑起来会明显变慢。如果遇到大数据量一个务实的方案是先用 MiniBatchKMeans 对数据做粗聚类得到几百个聚类中心再对聚类中心做 Isomap最后把每个样本映射到所属聚类中心的低维坐标上。虽然会损失一些细节但工程上可行实训平台里如果跑大型数据集卡死多半也要考虑这种折中方案。5. 常见问题与踩坑实录5.1 问题一报错“Array contains non-finite numbers”这是实训平台和作业里最常见的报错之一。原因几乎都是近邻图不连通最短路径矩阵里出现了 inf。排查步骤计算 geodesic_dist 后立刻检查 np.isinf(geodesic_dist).any()如果确实有 inf把 n_neighbors 往上调直到连通为止也可以改用距离阈值 ε 的邻域方式不过 ε 的选择同样需要试。一个更根本的排查思路是检查数据的量纲。比如数据里有的特征取值 0 到 1有的特征取值 1000 到 10000欧氏距离会被大数值特征主导近邻关系会失真。训练前做标准化StandardScaler能在源头上避免这类问题。5.2 问题二降维结果两端扭曲甚至折叠很多人调低 n_neighbors 后瑞士卷展开图看起来还可以但首尾两端会折叠到一起或者边缘处有明显的变形。这种情况通常是近邻数偏小导致局部图结构不够光滑测地线距离局部噪声过大。解决办法有三个适当增大 n_neighbors对距离矩阵做平滑处理比如先对近邻权重做核化或者检查数据是否存在异常点把个别离群点剔除后再跑。5.3 问题三Isomap 结果不如 t-SNE 好看经常有同学问为什么 t-SNE 聚类效果那么惊艳Isomap 相比之下“平平无奇”。这里要明确两者的目标差异t-SNE 优化的是高维和低维概率分布之间的 KL 散度目标是“局部邻居保持、全局距离可以忽略”所以画出来常常有几团清晰分离的簇。它的缺点也很明显距离信息失真、每次运行结果不稳定Isomap 优化的是全局距离保持它不理解“簇”的概念更关注的是流形本身的整体拓扑结构。所以在瑞士卷这类连续流形上Isomap 的展开效果比 t-SNE 更“有几何感”但看不出聚类的效果。如果项目目标是聚类可视化优先考虑 t-SNE 或 UMAP如果目标是还原数据真实的低维结构并保留全局距离关系再用 Isomap。这两者不是替代关系而是“局部”和“全局”两种不同视角。5.4 问题四降维后第一维和第二维含义不明确PCA 的第一主成分有“解释最大方差”的明确含义Isomap 的低维坐标没有这种可解释性。它只是“尽量保持测地线距离”的一个数值解。这会让很多习惯 PCA 的人不习惯但在非线性的流形结构里“可解释方差”这个说法本身就不成立所以没必要强求含义。5.5 问题五样本量太少时结果不稳定Isomap 本质上是依赖近邻图结构的如果数据只有几十个点近邻图非常稀疏测地线距离很容易被噪声干扰。此时要么增加数据量要么换用带正则化项的变体比如带参数的 L-Isomap。一句话总结我多年的经验Isomap 从来不是一个“开箱即用、默认参数走天下”的算法。它的上限很高但能不能发挥出来取决于你愿不愿意为它花时间调近邻数、做数据预处理以及正确判断数据是否躺在低维流形上。6. 等度量映射的扩展与变体Isomap 提出后学术界和工程界陆续做了不少扩展了解这些变体对做课程设计或项目选型都有帮助。L-IsomapLandmark Isomap用少量地标点替代全量点做最短路径和 MDS大幅降低时间复杂度适合样本量大但流形结构相对简单的场景。加权测地线距离在计算最短路径时降低跨高密度差区域的边权重缓解流形密度不均导致的测地线偏差。监督 Isomap在近邻图构建阶段引入类别标签的距离让同类样本更容易近邻相连适合分类任务里的特征提取。连接残差分析对测地线距离矩阵做残差分析可以估算流形的本征维数这个技巧在期末项目里写“方法选择理由”时特别好用。不过对于初学者我建议先把最基础的等度量映射实现、调参、可视化做扎实再去碰这些变体。因为所有变体都建立在“近邻图 最短路径 MDS”这个三步框架上框架理解了变体就是在某一步上做文章的产物。7. 个人实操体会与建议我自己做实验和带项目时最常被问到的不是 Isomap 理论而是“为什么我跑出来的结果和教程里不一样”。最常见的原因是数据生成方式。教程里用的瑞士卷往往只取了一部分形状或者噪声设置完全随机导致近邻图结构和你的数据不同。所以复现某个效果时先检查数据的范围、噪声、采样密度是否一致再谈参数差异。我见过有同学把 n_neighbors 从 5 调到 12 后效果就完全恢复的案例。另外用 Isomap 做特征提取后再接分类器时建议做交叉验证来选 n_neighbors而不是用可视化效果来拍板——因为可视化好看不一定对分类任务最优。这个细节在实训平台的“第 1 关等度量映射”这类任务里尤其重要很多人只看直观图觉得不错但提交代码后分数不高就是因为没有把参数选择放在任务指标上做优化。最后分享一个小技巧如果近邻数实在不好调可以试试把标准化后的数据先用 PCA 投影到 50 维左右再跑 Isomap。这一步能压掉大量噪声维度让近邻关系更稳定效果往往比在原始高维空间里硬调近邻数更好。我自己处理图像特征和高维统计特征时一直都是这么干的。
返回列表