简介这是一份Matlab实现kMedoids聚类的示例代码资源适合刚接触聚类分析、想要快速理解算法原理的数据挖掘学习者。压缩包内共1个m文件总体约2KB内容精炼代码脚本以少量随机生成的多维数据为例完整演示了从构造特征矩阵、调用内置kmedoids函数到指定聚类数量、距离度量方式如欧氏距离或曼哈顿距离的流程。读者还可以修改InitMethod、MaxIter、Tol等参数观察不同初始化和收敛条件对聚类结果的影响直观体会medoids以真实样本为中心、相比K-means质心更抗离群点的特点。目前已有237人学习下载对希望快速上手kMedoids聚类应用的Matlab用户来说是一份轻量而实用的参考将其替换为自己的数据矩阵即可进一步验证算法效果也便于在此基础上扩展市场细分、社区检测或图像分割等应用。1. kMedoids 在 matlab 里的真实用法从 Kmedios.m 看聚类能不能直接落到你的数据上谈到 kmedoids 聚类matlab 用户手里最常出现的资源就是类似 Kmedios.rar 里这样一个 Kmedios.m 脚本。脚本本身不难跑生成随机数据、调用内置 kmedoids、画散点图三步就出结果。但我在实际项目里帮人排查聚类问题时发现多数人卡住的从来不是函数调用而是「为什么每次结果都不一样」「K 到底取几」「有离群点时 medoid 和质心差在哪」。这份资源给了一个能用真实样本当簇中心的完整示例适合刚接触 kMedoids、想把聚类落到数据分群、图像分割或异常检测场景的工程师。我把脚本拆开逐段讲再把常见翻车点罗列出来。2. 为什么选 medoids 而不是 centroids原理差异与三种最适合换算法的场景2.1 medoids 是真实样本centroids 是虚拟质心kMedoids 和 K-means 的核心差别就一句话K-means 用算出来的均值当中心kMedoids 用数据集里真实存在的一个样本当中心。这个真实样本叫做 medoid它是在当前簇内找一个点使它到簇内其他所有点的距离之和最小。K-means 的中心是虚拟的它把簇内所有样本在每个维度上求平均。问题是均值天生怕离群点一个极端值就能把质心拉偏。kMedoids 的思路更接近中位数它不关心均值只关心「谁最能代表这群样本」所以对离群点的容忍度明显更高。另一个实际好处是medoid 就是原始数据里的一行你可以直接回溯到它对应的用户 ID、订单号或图像块这在业务上非常好解释。从目标函数看K-means 最小化的是簇内平方误差和SSE而 kMedoids 最小化的是样本到 medoid 的绝对距离总代价。这意味着它不要求距离是欧氏距离的平方你甚至可以换成曼哈顿距离、余弦距离这在 K-means 里就很难办因为 K-means 的均值计算隐含了欧氏空间假设。Matlab 内置的kmedoids函数正是围绕这个思想实现的而 Kmedios.m 脚本里的调用方式就是最标准的那一套。对比项K-meanskMedoids中心定义簇内均值虚拟点簇内真实样本medoid对离群点敏感均值被拉偏鲁棒用样本代表示范可解释性中心不是任何实例中心可回溯原始 ID支持距离主要是欧氏距离欧氏、曼哈顿、余弦等计算复杂度低适合大规模高交换阶段代价大典型场景图像压缩、海量日志分桶用户画像、社区检测、含噪数据2.2 三种最适合换算法的场景第一种是数据里明显有离群点。比如传感器数据、交易金额、用户时长这类带长尾的字段K-means 会把质心拖向离群点导致簇边界失真。kMedoids 因为中心是真实样本离群点只会成为它自己的一个孤簇不会污染别的簇。第二种是业务上需要「代表性个体」。做过用户分群的都知道平均用户往往不存在但「最典型的用户」真实存在。kMedoids 找出来的 medoid 可以直接做用户画像背书比如拿一个真实 VIP 用户代表整个高价值簇运营那边一眼就能看懂。K-means 给出一组浮点数的坐标业务侧完全无法落地。第三种是特征不是连续数值或者空间结构不适合欧氏距离。比如文本向量、类别型编码、社区网络里的节点K-means 在这些场景下算均值没有意义。kMedoids 配合曼哈顿或余弦距离仍然能给出「真实存在的代表节点」。图像分割里也常有人用 kMedoids 选代表色块因为 medoid 色块是真实像素组合出来的不会出现 K-means 那种不存在的混合色。2.3 计算复杂度PAM 的代价与应对手段kMedoids 的经典实现是 PAMPartitioning Around Medoids它在每个迭代里尝试用非 medoid 样本替换当前 medoid然后重新计算总代价。这个交换过程的理论复杂度是 O(k(n-k)²)比 K-means 的 O(nkd) 高出不少。也就是说几万行数据以内跑起来还舒服百万行级别就别指望直接呼救了。Matlab 内置的kmedoids在内部也是 PAM 思想但它支持Replicates做多次随机初始化并且InitMethod可以选PlusPlus用 K-means 的思路挑彼此距离远的初始 medoid能显著减少到达局部最优的概率。数据量大时我的常用做法是先随机抽样几千行跑一遍 kmedoids拿到 medoid 后再把它作为全量数据的中心做一次性分配。这样精度损失不大但速度能快一个量级。关于这一点Kmedios.m 里没写但实际项目里基本离不开它。3. Kmedios.m 逐段拆解随机数据生成、kmedoids 调用与 7 个参数的实际含义3.1 先生成一份能复现的随机数据集% 固定随机种子保证每次运行结果一致 rng(42); % 三个高斯簇的中心 mu [0 0; 6 6; -4 5]; % 簇内协方差矩阵两个维度带一点相关性 sigma [1 0.8; 0.8 1]; % 每个簇生成 200 个样本拼成 600x2 的数据矩阵 data [mvnrnd(mu(1,:), sigma, 200); mvnrnd(mu(2,:), sigma, 200); mvnrnd(mu(3,:), sigma, 200)];rng(42)的作用是把随机数生成器的状态固定下来这样任何人运行这段代码拿到的data都完全一样。做聚类调试时这一点特别重要否则你没法判断结果变化到底是算法问题还是数据变了。mvnrnd的第一个参数是均值向量第二个是协方差矩阵第三个是样本数。协方差矩阵里对角线是方差非对角线 0.8 表示两个维度正相关让簇的形状是倾斜的椭圆而不是正圆。这样更接近真实数据也能在可视化时看出 medoid 和几何中心的区别。这里三个簇中心距离分得比较开K3 是非常明显的答案。脚本里故意先给一个「一眼能看出答案」的数据目的是让你先验证调用流程没写错再替换成自己的数据。3.2 调用内置 kmedoids核心参数逐个说K 3; [idx, medoids_idx] kmedoids(data, K, ... Distance, euclidean, ... InitMethod, PlusPlus, ... MaxIter, 100, ... Tol, 1e-4, ... Replicates, 5);kmedoids最少只需要传数据矩阵和聚类数但实际使用我会把五个参数全写出来。第一个返回值idx是每个样本所属簇的编号长度和样本数一致第二个返回值medoids_idx是 medoid 在数据矩阵里的行号。注意它返回的是索引不是坐标本身要拿原始样本得再用data(medoids_idx, :)去取。下面这张参数表是我日常用的推荐值参数作用建议值Distance样本间距离度量方式优先euclidean有离群点换cityblockInitMethod初始 medoid 选择策略数据量不大用PlusPlus量大用random配 ReplicatesMaxIter最大迭代次数100~300不用贪大Tol收敛容差默认 1e-4 即可精度要求高可设 1e-6Replicates重复运行的次数返回最优结果5~10局部最优风险大时加InitMethod选PlusPlus是 K-means 的思路第一个 medoid 随机选后面的 medoid 倾向于选离已有 medoid 更远的样本。这样初始化质量高收敛快局部最优概率小。random就是完全随机挑 K 个样本当初始 medoid如果数据分布不够均匀很容易卡在局部最优。Replicates是后悔药机制跑 5 次取总代价最小的一次成本翻倍但结果稳定很多。MaxIter不是越大越好我见过有人设 5000结果运行时间长了一倍聚类效果和 200 次差不多。迭代到后期每次交换带来的代价改善极小Tol就是在检测这个「改善幅度」。当连续两次迭代的目标值变化小于Tol算法就认为收敛提前退出。3.3 从 medoids_idx 取回真实样本medoid_samples data(medoids_idx, :); for j 1:K fprintf(簇 %d 的 medoid 是第 %d 行样本: [%.2f, %.2f]\n, ... j, medoids_idx(j), medoid_samples(j,:)); end这一步是 kMedoids 最有价值的地方。medoids_idx给出的是数据矩阵的行号data(medoids_idx, :)能把这一行的完整特征取出来。如果你的原始数据附带了用户 ID、时间戳、业务标签直接data_table(medoids_idx(j), user_id)就能拿到「这个簇的代表性用户是谁」。K-means 输出的是一个浮点坐标你永远不知道它对应哪个人kMedoids 输出的是一个真实样本业务同事可以直接拿着这个样本去访谈、去打标签。做用户画像和社区检测时这一步几乎等于把聚类结果从技术语言翻译成了业务语言。3.4 用轮廓系数和 CH 指数评估聚类质量sil silhouette(data, idx); mean_sil mean(sil); fprintf(平均轮廓系数: %.4f\n, mean_sil); ch calinskiHarabasz(data, idx); fprintf(CH 指数: %.2f\n, ch);轮廓系数对每个样本计算一个值范围在 [-1, 1] 之间。接近 1 说明这个样本离自己簇很近、离相邻簇很远聚类结构清晰接近 0 说明在边界上负值说明可能被分错了簇。看mean(sil)的同时我强烈建议看一眼sil的分布如果大部分样本在 0.5 以上但有一小撮是负值问题往往出在离群点或 K 选太大。Calinski-Harabasz 指数是簇间离散度和簇内离散度的比值越大越好。它没有上界但适合在同一种距离度量下比较不同 K 值。CH 指数计算很快我在做 K 值扫描时一般先用它粗筛再用轮廓系数确认。4. 聚类结果可视化与质量评估散点图、silhouette 与 K 值扫描4.1 用 gscatter 染出簇归属用黑叉标出 medoidfigure; gscatter(data(:,1), data(:,2), idx, rgb, o, 8); hold on; plot(medoid_samples(:,1), medoid_samples(:,2), kx, ... MarkerSize, 12, LineWidth, 2); legend(簇 1, 簇 2, 簇 3, Medoid); title(kMedoids 聚类结果× 为 medoid); grid on;gscatter按idx的取值自动给每个簇分配颜色第三个参数是分组变量后面的rgb指定颜色顺序o指定点的形状8 是点的大小。plot把三个 medoid 用黑色叉号画在最上层。这样一张图能同时回答两个问题数据明显分成几堆、每个堆的代表样本在哪个位置。一个容易被忽视的细节是hold on。不写这一行第二个plot会开新图覆盖掉gscatter的画布你只能看到一个干巴巴的黑叉。Matlab 绘图里几乎所有「两个图形叠在一起」的用法都离不开hold on这是新手最常见的一处漏写。4.2 统计每个簇的样本数和簇内平均距离for j 1:K members data(idx j, :); avg_dist mean(pdist2(members, medoid_samples(j,:))); fprintf(簇 %d: 样本数 %d, 到 medoid 平均距离 %.3f\n, ... j, size(members, 1), avg_dist); end簇样本数能暴露聚类是否均衡。如果某个簇只有两三个样本大概率是离群点被单独成簇了这种情况就要回去处理数据而不是调参数。pdist2(members, medoid_samples(j,:))计算每个成员到该簇 medoid 的欧氏距离mean取平均。这个平均距离是簇内紧凑度的简单度量数值越小说明簇越紧。4.3 K 值扫描从 2 到 6 循环对比评估指标Ks 2:6; score_matrix zeros(length(Ks), 2); for i 1:length(Ks) [idx_i, ~] kmedoids(data, Ks(i), ... Distance, euclidean, ... InitMethod, PlusPlus, ... Replicates, 3); score_matrix(i, 1) mean(silhouette(data, idx_i)); score_matrix(i, 2) calinskiHarabasz(data, idx_i); fprintf(K%d Silhouette%.4f CH%.2f\n, ... Ks(i), score_matrix(i, 1), score_matrix(i, 2)); end聚类数和聚类效果不是单调关系K 越大簇内越紧凑但模型越复杂、越容易过拟合。用循环把 K 从 2 到 6 全部跑一遍对比 silhouette 和 CH 的走势是最务实的选 K 办法。理想情况下 silhouette 在某个 K 值出现峰值或者 CH 在某个 K 值出现明显拐点再往后增长放缓就取那个 K。有一点要提醒silhouette函数对大数据集很慢如果样本量超过几万可以改成silhouette(data, idx, Distance, euclidean)配合抽样或者直接用evalclusters里的并行选项。5. kMedoids 高频避坑与排查5 个实测翻车点每条都有解决办法5.1 报错 Undefined function or variable kmedoids现象是代码一行没跑就开始报错提示找不到kmedoids函数。我第一次遇到还以为是 Matlab 没装好折腾了半天重装环境。原因是kmedoids是 R2018a 才加入 MATLAB 的函数而且属于 Statistics and Machine Learning Toolbox。旧版本或者精简安装没带这个工具箱都会报同样的错。先运行ver(stats)看看工具箱版本如果显示Statistics and Machine Learning Toolbox检查一下版本号如果根本没显示说明工具箱没装上。解决办法分两种。一是装工具箱或者升版本这取决于你的 license。二是换一个不依赖工具箱的纯脚本实现思路是把 PAM 算法用朴素循环写出来。function [idx, medoid_rows] simple_pam(data, K, max_iter) n size(data, 1); rng(0); medoid_rows randsample(n, K); idx zeros(n, 1); for iter 1:max_iter % 每个样本归到最近的 medoid d pdist2(data, data(medoid_rows, :)); [~, idx] min(d, [], 2); % 在每个簇内找新的 medoid到簇内其他样本距离和最小 new_med medoid_rows; for j 1:K members find(idx j); if isempty(members), continue; end sub_d pdist2(data(members, :), data(members, :)); [~, pos] min(sum(sub_d, 2)); new_med(j) members(pos); end if isequal(new_med, medoid_rows), break; end medoid_rows new_med; end d pdist2(data, data(medoid_rows, :)); [~, idx] min(d, [], 2); end这段实现抓的是 PAM 的核心先随机挑 K 个样本当 medoid然后循环做两件事——分配样本到最近 medoid再在每个簇内挑一个「到簇内其他点距离总和最小」的样本当新 medoid。pdist2计算两两距离矩阵sum(sub_d, 2)对每一行求和最小的那行对应的就是簇内的新 medoid。性能远不如内置函数但几千行数据足够用。5.2 每次运行结果都不一样K 值扫描曲线剧烈抖动现象是连续跑两次同一个脚本聚类结果完全对不上medoids_idx变化很大连散点图的簇颜色都变了。原因是InitMethod默认或显式设成了random初始 medoid 每次都是随机选的而 PAM 的交换过程容易落到局部最优。不同的初始点会走向不同的局部最优表现就是每次运行结果不同。这也是 kMedoids 最容易被吐槽「玄学」的地方但问题不在算法在初始化。解决方法是固定随机种子并加大重复次数。rng(42)放在调用前Replicates设为 5 或 10让 Matlab 从多个初始点里挑代价最小的结果。如果数据量允许InitMethod换成PlusPlus也能明显降低抖动。从那以后我养成的习惯是任何聚类脚本第一行写rng(0)跑批结果必须能复现否则后面所有分析都建立在不稳定的地基上。5.3 数据量一大就跑不动甚至直接卡死现象是数据从几千行换成几万行脚本执行时间从秒级变成分钟级再往上直接内存吃满卡死。原因是 PAM 的复杂度是 O(k(n-k)²)n 是样本数。这个平方项非常致命n 翻一倍时间接近翻四倍。另外pdist2会生成 n×n 的距离矩阵几万行数据就是几十 GB 内存直接爆掉。这是 kMedoids 的结构性代价不是优化两行代码能解决的。解决思路是降维打击先对数据做随机抽样比如抽 5000 行跑kmedoids拿到 medoid然后用knnsearch或pdist2把剩余样本分配到最近的 medoid 上。抽样只影响 medoid 的选择精度不影响整体结构速度却能快一个量级。这一步我习惯写成脚本的一部分而不是单独的手工操作因为只要数据更新就会再跑一遍。另外还可以换个距离度量。cityblock曼哈顿距离的计算比欧氏距离少了平方和开方对pdist2的性能有一定改善。不过这是次要手段采样才是主力。5.4 轮廓系数很高但业务上分出来的簇完全不合理现象是 silhouette 到 0.7 以上感觉聚类质量非常好一细看却发现高频用户和低活跃用户被分到了同一个簇或者某个簇里全是异常样本业务根本无法解释。原因是对聚类质量的评估太依赖数值指标而 silhouette 和 CH 都是纯几何指标完全不考虑业务语义。数据里的离群点经常会被单独分出来当一簇这在几何上「很清晰」在业务上却毫无意义。另外如果特征里有一个数值波动特别大它会在距离计算里占主导把真正重要的业务字段淹没了。解决办法是先做特征标准化再聚类。zscore把每个特征缩放到均值 0、方差 1避免量纲大的字段主导距离。然后在评估前先看每个簇的样本数和业务特征分布再决定要不要剔除离群点。我常用的一招是聚类完了先跑一段描述性统计看看每个簇的均值、中位数、业务占比几何指标只做参考业务解释不通就回退调参。5.5 数据里有 NaN聚类结果莫名出现空洞或报错现象是调用不报错但idx里有些样本被分到了明显错误的簇或者silhouette直接报NaN相关的错误还有的版本会直接终止运行。原因是kmedoids对NaN的处理不透明。部分旧版本会把NaN样本直接忽略导致输出idx长度和输入对不上有些版本的NaN会被当成一个特殊值参与距离计算结果整个簇都被污染。聚类函数不是fillmissing它不会好心替你处理缺失值。解决方法是进聚类之前主动清洗数据。对数值型特征用fillmissing(data, linear)做插值对离散型特征可以用rmmissing删掉有缺失的行。注意不要只清洗一个变量rmmissing默认删除的是整行有缺失的样本如果缺失率超过 5%删行会损失太多信息改成fillmissing更稳妥。清洗完画一下缺失值占比确认没有列还带着大量空洞再往下走。6. 进阶技巧把 silhouette 扫描封装成函数告别手动试 K6.1 一个通用的自动选 K 函数基于前面第 4 章的循环思路我习惯把它封装成一个可复用函数放进项目公共工具目录里。function [bestK, bestIdx, bestMedoids, stats] auto_kmedoids(data, Krange) stats zeros(length(Krange), 2); bestK Krange(1); bestScore -inf; for i 1:length(Krange) K Krange(i); [idx, medRow] kmedoids(data, K, ... InitMethod, PlusPlus, Replicates, 3); s mean(silhouette(data, idx)); c calinskiHarabasz(data, idx); stats(i, :) [s, c]; % 优先看 silhouetteCH 只做辅助参考 if s bestScore bestScore s; bestK K; bestIdx idx; bestMedoids data(medRow, :); end end end调用时传数据矩阵和一个 K 的候选范围比如[2 3 4 5 6]返回分数最高的 K 以及对应的聚类结果。这个函数的价值在于把「扫 K、看指标、选最优」三件事合成一步数据更新了直接重新调用不用每次都复制粘贴循环代码。判断逻辑只取 silhouette 峰值因为 CH 指数没有上界不同 K 之间的差值不好设定阈值实际对比里 silhouette 更直观。6.2 替换距离度量的两个边界条件kmedoids的Distance参数支持euclidean、cityblock、cosine等选项。用欧氏距离时 medoid 是「几何上最居中」的真实样本换成曼哈顿距离后medoid 会更接近簇内的中位数样本对重尾分布更友好。做高维稀疏数据比如文本向量或用户行为计数矩阵时cosine是更合理的选择因为它只看方向不看模长。但要提醒一点换距离度量后silhouette 的Distance参数也要同步换否则评估用的距离和聚类用的距离不一致指标会失真。6.3 结尾你的验证习惯是什么这套流程我前前后后用在不同项目里踩过的坑比上面写的还多。现在每次拿到一份新数据我已经形成了固定动作第一件事rng(0)第二件事清洗 NaN第三件事跑一遍auto_kmedoids扫 K最后才看业务解释。这三步做完聚类结果基本不会翻车。kMedoids 不是银弹但它给的是「真实存在的代表样本」光这一点在业务沟通里就值回票价。希望帮到你。本文还有配套的精品资源点击获取