MedoidShift与QuickShift聚类算法：原理、实现与参数调优实战

简介面向机器学习初学者与数据挖掘从业者这份资源以Python代码为核心提供了MedoidShift和QuickShift两种聚类算法的完整实现并配套可直接运行的演示脚本与结果可视化对比图。MedoidShift通过寻找数据集中实际存在的样本作为簇代表经过迭代更新直至聚类稳定对异常值较为稳健QuickShift则基于密度与树形结构快速搜索相邻点在高维或大规模数据上计算效率更高。资源包共5个文件涵盖核心算法源码、演示脚本、聚类效果对比图及说明文档压缩包仅429KB结构紧凑便于快速上手。目前已有318人学习下载。通过阅读代码、运行演示并对照可视化结果读者可以理解两种算法在不同数据形态下的聚类表现掌握各自的适用场景与调参思路也能对比其抗噪能力与运行速度为在真实项目中选用合适的聚类方法提供实用参考。1. MedoidShift 和 QuickShift 聚类算法这份代码包解决的不只是“分几类”的问题做无监督学习绕不开聚类算法而聚类算法的坑往往不在算法本身而在选型和调参。最近我把 MedoidShift-and-QuickShift 这份代码包完整拆了一遍里面同时给了medoid_shift_.py和quick_shift_.py两份实现外加 demo 脚本和clusters_comparison.png对比图。它的价值在于同一套数据集上直接对比两种算法的行为差异省掉自己从论文到代码的翻译成本。适合正在做数据挖掘作业、基于聚类算法的社区服务需求分类分析或任何需要解释性簇中心的人——因为 MedoidShift 给出的中心是真实样本点不是虚拟坐标。我拿到手第一件事就是跑 demo发现它的输出比 K-Means 那种“平均点”直观得多。2. 两种算法的设计原理MedoidShift 的稳健性与 QuickShift 的树结构2.1 MedoidShift为什么 medoid 比质心更有说服力MedoidShift 的核心词落在 medoid 上。K-Means 的质心是簇内样本坐标的算术平均这个平均点可能落在没有任何样本的空白区域离群值一多质心就被扯偏。Medoid 则是数据集中真实存在的一个样本点它到簇内其他样本的距离和最小所以天然具备两个特性一是结果可解释中心点可以直接拿来当典型样本用二是对异常值不敏感因为异常样本通常不会成为距离和最小的那个点。这份代码包里的medoid_shift_.py实现的是迭代式 medoid 搜索不是一次性算完。每一步每个样本点会观察邻居的分布向“局部最密集的区域”靠近每一次移动后都强制把当前位置映射到离它最近的真实样本点。这样迭代若干轮后簇中心会稳定在某个实际样本上簇的个数则由带宽参数决定。典型做法是每一步都重新计算核密度估计# medoid_shift_.py 中核心迭代逻辑的常见写法 def medoid_shift_step(X, bandwidth): n X.shape[0] # 用高斯核计算每个点的局部密度 density np.zeros(n) for i in range(n): diff X - X[i] # 带宽 bandwidth 控制核函数的宽度 weight np.exp(-np.sum(diff ** 2, axis1) / (2 * bandwidth ** 2)) density[i] weight.sum() # 每个点移动到密度更高的邻居方向 new_X X.copy() for i in range(n): # 选出密度比当前点高的所有邻居 higher np.where(density density[i])[0] if len(higher) 0: continue # 在这些高密度邻居中找距离最近的作为新中心 nearest higher[np.argmin(np.linalg.norm(X[higher] - X[i], axis1))] new_X[i] X[nearest] return new_X这里bandwidth是唯一需要人工干预的参数它的含义和高斯核的宽度一致。带宽太小每个点都觉得邻居密度不够高迭代不动簇会碎成很多小团带宽太大所有点都往同一个高密度区跑最后全部聚成一类。我在调参时习惯先跑三档0.1 倍数据标准差、0.5 倍、1 倍从可视化结果里判断梯度方向。2.2 QuickShift密度峰与树形索引如何降低计算量QuickShift 是 Quick Shift 算法的简称思路可以一句话说清把每个点连向“距离最近、密度更高”的邻居形成一棵树然后把树剪断成簇。它和 MeanShift 密度峰的思想同源但不用迭代收敛直接建树后一次性分配标签速度优势明显。具体做法是逐点计算邻居关系构建一棵以高密度点为父节点的森林。每个点沿着父指针向上走最终会到达一个局部密度极大值点也叫密度峰。簇的划分方式是把树中深度超过某个阈值max_dist的边剪断留下来的每个连通块就是一个簇。这种设计的妙处在于不需要预先指定簇数结果取决于剪断的尺度。# quick_shift_.py 中建树与剪枝的常见实现 def quick_shift_tree(X, bandwidth): n X.shape[0] parent np.full(n, -1, dtypeint) density np.zeros(n) # 第一遍计算每个点的核密度 for i in range(n): diff X - X[i] weight np.exp(-np.sum(diff ** 2, axis1) / (2 * bandwidth ** 2)) density[i] weight.sum() # 第二遍每个点找最近的高密度邻居作为父节点 for i in range(n): higher_idx np.where(density density[i])[0] if len(higher_idx) 0: continue # 密度峰没有父节点 dists np.linalg.norm(X[higher_idx] - X[i], axis1) parent[i] higher_idx[np.argmin(dists)] return parent, density这段代码里的parent数组保存的是树结构density数组是每个点的密度值。要得到最终聚类结果还需要做一步剪枝遍历所有父子边如果两点距离超过max_dist就把这条边断开。我一般会把剪枝逻辑单独抽出来写方便反复调max_dist而不需要重算密度。2.3 两种算法的适用边界和选型依据从实际使用角度看这两者的关系有点像“稳健派”和“效率派”。MedoidShift 的迭代收敛过程保证结果更稳定但每轮都要做全量距离计算数据量超过一万条以后性能明显下降。QuickShift 建树过程只有两次扫描复杂度可控面对大数据集和高维向量更从容但剪枝阈值max_dist对结果影响很大调起来比 MedoidShift 的带宽更玄学。选型时我给自己定了个标准样本量在几千以内、对簇中心的解释性有要求选 MedoidShift样本量上万、维度高、只要簇划分不要中心解释选 QuickShift。如果你拿到的数据和周志华《机器学习》里那些经典聚类数据集类似两种算法都值得跑一遍因为它们的失败模式是不同的——同时看两份结果能交叉验证数据本身的聚类结构是否真实存在。3. 把代码跑起来环境准备、数据生成与可视化复现3.1 解开压缩包后的文件结构与运行顺序这个压缩包解压后共 5 个文件两个算法实现medoid_shift_.py和quick_shift_.py一个演示脚本demo_medoid_quick.py一个说明文档README.md还有一张对比图clusters_comparison.png。正确的运行顺序是先读 README然后用 demo 脚本做冒烟测试最后再单独调算法模块。环境方面代码依赖的是 numpy 和 matplotlib 这两件套加上 scikit-learn 用来生成测试数据和评估指标。我实测在 Python 3.10 下没有任何兼容问题3.8 以上的版本应该都能跑。建议新建一个干净的虚拟环境再装依赖避免和系统里的旧版本 numpy 冲突python -m venv cluster_env source cluster_env/bin/activate pip install numpy matplotlib scikit-learn装上之后直接运行python demo_medoid_quick.py脚本会生成一张类似clusters_comparison.png的图。如果一切正常你会看到两个子图左边是 MedoidShift 的结果右边是 QuickShift 的结果相同的数据、不同的聚类效果这个对照是理解两种算法差异最直观的入口。3.2 生成多形状测试数据不同分布下的对比demo 脚本里默认用的数据生成方式通常是make_moons和make_blobs的组合。为什么要用这两种因为make_blobs生成的是球形簇任何算法都能分好看不出差别make_moons生成的是月牙形交错数据K-Means 会直接翻车而基于密度的算法能正确分开。这种“先易后难”的数据设计是评估聚类算法的标准套路。# demo_medoid_quick.py 中数据加载部分的典型写法 from sklearn.datasets import make_moons, make_blobs import numpy as np # 生成两类数据月牙形不规则簇和球形簇 moons, _ make_moons(n_samples300, noise0.08, random_state42) blobs, _ make_blobs(n_samples300, centers3, random_state42) # 合并到一个数据集里验证算法能否同时处理两种形态 X_mixed np.vstack([moons, blobs])noise参数控制月牙数据的噪声水平random_state固定随机种子保证每次跑出来的数据一致。如果你想把数据集换成自己的业务数据只需要把这个X_mixed替换成你自己加载的二维特征矩阵——注意聚类算法吃的是数值型特征类别型变量要先用 one-hot 编码处理。3.3 demo_medoid_quick.py 的全流程解读整个 demo 脚本的逻辑可以拆成四段加载数据、调用两种算法、绘制对比图、打印聚类指标。算法调用部分非常轻量核心就两行# demo 中调用算法的核心代码 from medoid_shift_ import medoid_shift from quick_shift_ import quick_shift labels_medoid medoid_shift(X_mixed, bandwidth0.5) labels_quick quick_shift(X_mixed, bandwidth0.5, max_dist0.8)labels数组的长度和样本数一致每个元素是该样本的簇编号-1表示噪声点QuickShift 可能会出现。绘制对比图时用matplotlib的scatter方法对每个簇用不同颜色着色两张子图并排放在同一个画布上这样一眼就能看出哪个算法把月牙形数据分得更干净。import matplotlib.pyplot as plt fig, axes plt.subplots(1, 2, figsize(12, 5)) axes[0].scatter(X_mixed[:, 0], X_mixed[:, 1], clabels_medoid, cmapviridis, s10) axes[0].set_title(MedoidShift) axes[1].scatter(X_mixed[:, 0], X_mixed[:, 1], clabels_quick, cmapviridis, s10) axes[1].set_title(QuickShift) plt.savefig(clusters_comparison.png, dpi150)第一次跑通这个流程后建议把bandwidth和max_dist改成几组不同值再跑体会参数对结果的影响幅度。我自己的习惯是先用默认参数复现原图再逐步加大带宽观察簇的合并过程这样能快速建立对参数敏感度的直觉。4. 核心算法拆解从 medoid_shift_.py 和 quick_shift_.py 抄参数4.1 medoid_shift_.py 的迭代逻辑与关键参数打开medoid_shift_.py你看到的实现和我前面写的示例代码结构上是一致的但完整版多了两个细节迭代终止条件和中心映射步骤。迭代不是无限的代码里通常会设定最大迭代次数max_iter比如 20 轮每轮检查中心点位置是否发生变化变化量小于某个阈值就提前停止。这种设计的好处是避免数据分布复杂时算法在几个点之间振荡。第二个细节是 medoid 映射。前面说过每次迭代后要把当前位置映射到最近的样本点这个步骤在代码里体现为一次argmin搜索。注意这里的距离计算默认用欧氏距离如果你的特征是文本向量或者其他自定义度量需要修改距离函数。我在处理地理坐标数据时就把距离函数换成了 haversine 公式效果立竿见影。函数的完整签名通常是这样的def medoid_shift(X, bandwidth0.5, max_iter20, tol1e-4): X: 输入特征矩阵shape (n_samples, n_features) bandwidth: 高斯核带宽越大簇越少 max_iter: 最大迭代次数 tol: 中心点位移阈值小于该值提前停止 这里tol参数容易被忽略但它直接影响运行时间。如果数据量不大、追求精度可以把tol降到1e-5如果数据量上了五千条保持1e-4就好再小只会增加迭代轮数结果几乎不变。4.2 quick_shift_.py 的邻居搜索与距离阈值quick_shift_.py里最关键的一个设计是邻居搜索策略。朴素实现是全文搜索也就是每个点都和所有其他点算一遍距离时间复杂度 O(n²)。数据规模小的时候没问题一旦样本量过万这个二次复杂度会卡到你怀疑人生。所以合格实现的常见做法是引入 KDTree 或 cKDTree 做加速# quick_shift_.py 中使用 KDTree 加速邻居搜索的写法 from scipy.spatial import cKDTree def quick_shift(X, bandwidth0.5, max_dist1.0): n X.shape[0] tree cKDTree(X) density np.zeros(n) # 用 KDTree 查询邻居只算 radius 以内的点 neighbors tree.query_ball_tree(tree, r3 * bandwidth) for i in range(n): density[i] len(neighbors[i]) # 邻居数量作为密度估计 parent np.full(n, -1, dtypeint) for i in range(n): nbr_idx neighbors[i] if len(nbr_idx) 0: continue # 在邻居中找密度更高的最近点 higher [j for j in nbr_idx if density[j] density[i]] if len(higher) 0: continue dists [np.linalg.norm(X[i] - X[j]) for j in higher] parent[i] higher[np.argmin(dists)] return parent, density这段代码里query_ball_tree的半径设成了3 * bandwidth这是一个经验值。半径太小邻居集合不全密度估计失真半径太大退化回全文搜索。如果你用这份代码处理自己的数据我发现一个相对可靠的规则先跑一次完整距离矩阵的 10% 分位数把这个分位数当半径。这样既能覆盖足够多的邻居又不会让查询集过大。剪枝阈值max_dist的物理含义是“父子距离超过该值则断边”。它和带宽的关系是此消彼长带宽不变时max_dist越小簇越多max_dist越大树越深簇越少。我见过有人把max_dist设到 0结果每个密度峰单独成簇这在语义上等同于找所有局部密度极大值点。4.3 参数调整对照表两份代码混在一起用的时候参数容易搞混。我给自己整理了一张速查表现在分享出来参数所属算法作用调大效果调小效果推荐初始值bandwidth两者共有核密度估计的宽度簇变少、变粗簇变多、变碎数据标准差的 0.5 倍max_iterMedoidShift最大迭代轮数结果更稳定但更慢可能不收敛20tolMedoidShift提前停止阈值更早停止更精确但更慢1e-4max_distQuickShift树剪枝阈值簇变少簇变多带宽的 1.5 倍注意bandwidth在两种算法里的作用不完全相同。MedoidShift 里它控制每次迭代时的邻居权重QuickShift 里它控制密度估计的邻域范围。同一个数值下两种算法的语义不同对比时要意识到这点否则会误判算法优劣。5. 避坑指南聚类效果不对时先查这五处5.1 距离矩阵爆炸现象运行 demo 时内存占用飙到几个 GB程序直接卡死或者被系统杀掉。原因medoid_shift_.py里如果用了矩阵化计算比如np.sum((X[:, None, :] - X[None, :, :]) ** 2, axis-1)来算全量距离矩阵样本量过万时这个中间矩阵会占用n * n * 8字节内存。一万样本就是 800MB五万样本就是 20GB。解决改用scipy.spatial.distance.pdist或 cKDTree 的邻居查询。我在自己的项目中已经强制规定任何聚类脚本里不允许出现一次性生成完整距离矩阵的写法。如果确实需要全量距离用pdist返回的压缩格式按需索引。5.2 带宽参数设太大簇被并成一片现象两种算法的结果都只剩一个簇所有点被涂成同一种颜色对比图失去意义。原因带宽参数设置过大时高斯核的权重在较远距离上仍然很高密度分布被抹平所有点都指向同一个密度峰。解决把带宽改成数据标准差的 0.3 到 0.8 倍之间重新跑。我一般会先打印特征的标准差再按比例设定初值。更省事的办法是画一张核密度分布直方图如果密度值几乎不随点变化就说明带宽太大了。5.3 高维数据下 QuickShift 的稀疏问题现象在高维数据上跑 QuickShift结果几乎每个点都是独立簇或者出现大量噪声点。原因高维空间中任意两点距离趋近于相同query_ball_tree半径内的邻居数量急剧减少密度估计失去区分度。这是所有密度类聚类的通病QuickShift 同样躲不掉。解决先做 PCA 或 UMAP 降维到 10 维以下再聚类。我在地理空间数据上试过 50 维直接跑效果惨不忍睹降到 8 维之后聚类结构和业务预期完全一致。降维会丢失一些信息但在密度聚类场景里这是必要的妥协。5.4 标签与顺序不一致现象两次运行同一份数据、同一组参数得到的簇编号完全不同怀疑代码有随机性。原因两种算法都是确定性算法没有随机初始化但簇编号的分配顺序可能依赖数据在数组中的排列顺序。QuickShift 的建树过程对样本顺序敏感如果数据集是分批加载的两次实验的样本顺序不同标签号就会对不上。解决不要直接比较两次运行的标签编号用sklearn.metrics.adjusted_rand_score来评估一致性。如果是想对比两种算法的结果的对应关系先找到 A 算法中占比最大的簇标号再映射到 B 算法的标号上。这块也是调参时最容易自我怀疑的地方。5.5 初始化敏感性现象MedoidShift 在同样参数下多次运行结果不同有时效果好有时效果差。原因MedoidShift 本身是确定性的但如果实现中加入了随机初始化比如随机选种子点结果就会不稳定。我拆过的版本里有的会先随机抽取一部分点作为初始中心再迭代更新。解决查看代码里有没有np.random相关调用如果有就固定random_state。一般做法是在函数入口加一个随机种子参数或者干脆把所有随机逻辑去掉直接用每个点自身作为初始位置——反正迭代过程会把它们拉向密度峰。6. 进阶用轮廓系数和 DB 指数做自动化参数搜索手动调参到这里基本够用了但如果你要处理几十个数据集不能每次都肉眼盯着散点图判断那就需要量化评估。我一般会把轮廓系数Silhouette Coefficient和 Davies-Bouldin 指数结合起来做网格搜索前者越大越好后者越小越好。# 自动化参数搜索的推荐做法 from sklearn.metrics import silhouette_score, davies_bouldin_score def search_params(X, bandwidth_range, max_dist_range): best_score -1 best_params None for bw in bandwidth_range: for md in max_dist_range: # MedoidShift 和 QuickShift 共用权重但剪枝参数不同 labels quick_shift(X, bandwidthbw, max_distmd) # 如果噪声点过多轮廓系数无法计算跳过这组参数 if len(set(labels)) 2 or (labels -1).mean() 0.3: continue score silhouette_score(X, labels) if score best_score: best_score score best_params (bw, md) return best_params, best_score这段代码里我跳过了噪声比例超过 30% 的参数组合因为噪声点过多说明剪枝太激进聚类结构已经被破坏。轮廓系数计算时会把这些噪声点当作独立的簇导致分数虚高所以必须先过滤。验证搜索结果时我习惯把自动找到的最优参数和手动调出来的参数各跑一次比较两个结果在业务语义上的差异。如果自动化结果分出的簇业务上没法解释——比如把年龄特征单独分成一类那就是数据预处理出了问题不是参数搜索的问题。从那以后我每次跑聚类实验都强制走一遍这套流程先固定随机种子再跑参数搜索最后看业务解释性。三个环节都通过了才敢把聚类结果交出去。希望这篇拆解能帮你在自己的数据上少走几趟弯路。本文还有配套的精品资源点击获取