看到这份标题第一反应是很典型的“课程作业/入门项目”配置——k-medoids 聚类算法的 MATLAB 源代码数据导入和画图部分带中文注释。别看它结构简单这套东西背后把聚类算法的几个关键环节全串起来了数据读入、距离计算、中心点更新、结果可视化。我已经帮多位同学把这个模型从“能跑通”调到“能讲清”今天就把这份代码的完整拆解、踩坑记录和可复现步骤一起整理出来。先说清楚适用范围这份代码适合入门机器学习、模式识别课程设计或者需要在 MATLAB 里快速做无监督聚类分析的开发者。如果你的数据量在几千条以内、维度不超过 5 维这份代码的训练速度和内存占用都能接受。要处理百万级样本那就得换 Mini-Batch 方案或上 Python 版本了。1. 项目背景与 k-medoids 算法原理拆解1.1 为什么选 k-medoids 而不是 k-meansk-means 是大多数人第一个接触的聚类算法它的核心计算是“簇内均值”。均值这个东西有两个天然缺陷第一它会被离群点带偏。比如某个簇里混入一个极端值簇中心会被硬生生拉过去导致边界划分错位。第二均值点是一个虚拟点不一定对应真实存在的样本后续分析时想“看看这个簇的代表长什么样”得额外找最近的真实点很麻烦。k-medoids 的做法完全不同它强行要求每个簇的中心必须是簇内真实存在的样本点这个点就叫 medoid中心样本。每次迭代只在“现有样本”里选代表天然免疫离群值而且聚类结果可以直接解释成“这个点就是本簇最典型的用户/商品/案例”。如果你处理的是金融交易记录、传感器异常数据、图像特征向量这类离群值占比高的数据k-medoids 的稳定性优势非常明显。当然代价也有——每次更新中心点要做多次相异度比较计算量比 k-means 大一截。1.2 PAM 算法的两步迭代逻辑MATLAB 源码里最常见的 k-medoids 实现是 PAMPartitioning Around Medoids算法它的思路拆开就两步构建阶段的目的是选初始中心。做法是遍历所有样本选出一个对象作为第一个中心使得它到其余所有样本的相异度总和最小然后反复挑选能使“整体相异度增量”最小的新样本加入中心集直到攒够 k 个中心。交换阶段则是循环优化对每一对“当前中心点 m_i”和“非中心点 h”计算把 m_i 换成 h 之后代价的变化量。如果代价减小就执行替换重复这个操作直到任何替换都无法再降低代价。这里的“代价”通常定义为每个普通样本到其所属最近中心点的距离之和。代码里核心循环就是围绕这个代价函数转的。理解了这两个阶段你就知道了源码里 200 行左右的逻辑主线是什么——它不像神经网络那样有复杂的反向传播所有核心都集中在“替换中心点是否划算”这一个判断上。提示PAM 算法对初始中心比较敏感不同初始点导致最终聚类结果差异明显。代码里如果没做多次随机初始化取最优的操作建议自己加上后面我会给具体做法。2. 代码框架与数据导入部分的实现解读2.1 完整代码模块划分整份源码按功能可以切成四个模块从上到下依次是数据准备模块、算法初始化模块、PAM 迭代主循环模块、结果可视化模块。数据准备和可视化是注释最密集的区段这是有意的安排——这两块高度依赖外部环境文件路径、坐标系、数据范围不写清楚使用者根本没法改。数据准备模块的核心任务就一句话把存储在外部文件里的矩阵读进 MATLAB 工作区并做初步清洗。这里涉及三个问题文件格式、路径写法、空值处理。源码里默认支持 .xlsx、.txt 和 .mat 三种格式用一行 switch 分支控制。2.2 数据加载区注释翻译与设计要点先看导数据这一段典型的源码注释和实现逻辑大致是这样% 导入数据部分 % 请根据你的实际数据文件格式修改文件名 file_name 和路径 file_path % 数据格式要求每一行是一个样本前 d 列是特征最后一列可以是标签可省略 file_path ./data/; file_name iris_data.xlsx; data readmatrix([file_path, file_name]); % 如果文件最后一列是标签先剥离标签列仅保留特征矩阵 if size(data, 2) size(data, 2) - 1 % 占位写法实际以代码为准 X data(:, 1:end-1); true_labels data(:, end); else X data; end设计要点有三个。第一是readmatrix这个函数它在较新版本的 MATLABR2019b 之后里能自动识别数值和文本混合的 Excel 表格老代码里常用的xlsread已经标了 deprecated 状态新写代码建议直接上readmatrix。第二是“剥离标签列”这个操作很多公开数据集最后带真实标签聚类是无监督算法不需要标签参与训练但后面对比聚类准确率时要用所以注释里专门提示用户保留。第三是./data/这种相对路径写法它假设你当前工作目录和控制脚本在同一级避免了绝对路径在换电脑后失效的坑。这里面最容易被忽略的是数据标准化。如果不同特征列的数值范围差异很大比如一列是 0~1另一列是 0~10000距离计算会被大数值列主导聚类效果几乎必然不理想。源码注释里如果没有提醒这一步我强烈建议你导入数据后加一行% 标准化处理每个特征列减去均值并除以标准差消除量纲影响 X (X - mean(X)) ./ std(X);把标准化注释也写进中文注释里之后调参时不容易忘。3. 核心聚类过程的 MATLAB 实现3.1 初始化策略PAM 的初始化有两种常见写法。第一种是随机抽 k 个样本做 medoid实现最简单但结果方差大。第二种是“最大最小法”第一个中心选离所有样本最远的点之后每个新中心都选“离已有中心最近距离最大”的样本这样中心散布较均匀比纯随机稳定得多。源码里通常提供的是随机初始化注释里会提示运行多次取最优。我的建议是直接改成最大最小法代码量增加不超过十行% 初始化 medoid采用最大最小距离策略使初始中心尽可能分散 rng(42); % 固定随机种子保证结果可复现 medoid_idx zeros(k, 1); % 先选距离矩阵总和最大的样本作为第一个中心简化为随机选第一个 medoid_idx(1) randi(n); for i 2:k % 计算每个非中心样本到已有最近中心的距离取最小值中的最大值 dist_to_center min(dist_mat(:, medoid_idx(1:i-1)), [], 2); [~, new_idx] max(dist_to_center); medoid_idx(i) new_idx; end固定随机种子比如rng(42)是个好习惯。它能保证你每次运行结果一致方便调试和比较实验——尤其在做参数调优时不确定的初始点会让前后结果没法公平对比。3.2 分配与更新循环主循环的代码逻辑是容易绕晕的地方我把伪代码结构抽出来看主循环分为两个阶段分配阶段根据当前 medoid 集把每个样本归到距离最近的 medoid 所属簇。这一步用距离矩阵的行比较实现。源码里最常见的做法是先预计算一个n x k的矩阵第 i 行存的是样本 i 到每个 medoid 的距离然后取最小值所在列索引。更新阶段对每个簇遍历簇内所有样本计算把它设为新 medoid 后簇内总距离的变化量。若总距离更小则替换。这个阶段最耗时因为涉及“簇内每个点都试一遍”的穷举式搜索。更新阶段的计算量能占到整个程序的 70% 以上。如果样本量大这一段建议改成向量化实现——用矩阵运算替掉内层循环而不是一个一个算距离% 更新 medoid 的向量化近似写法 for c 1:k members (assignments c); cluster_data X(members, :); % 计算簇内两两距离矩阵选行和最小的样本作为新 medoid intra_dist pdist2(cluster_data, cluster_data, squaredeuclidean); [~, new_medoid_in_cluster] min(sum(intra_dist, 2)); new_idx find(members); medoid_idx(c) new_idx(new_medoid_in_cluster); end这样写的好处有二一是pdist2内部用 BLAS 优化过的矩阵乘法算距离比自己写双重 for 循环快不少二是行和最小的样本理论上就是到簇内所有样本距离和最小的点正符合 medoid 的定义。代价是额外开了一整块n^2的内存数据量上万时留意内存占用。3.3 距离矩阵计算的向量化技巧距离度量直接决定聚类形态。源码里常见写法是支持三种距离欧氏距离、曼哈顿距离、余弦距离。欧氏适合连续实数特征曼哈顿L1在高维稀疏数据上更稳它不像欧氏那样对大数值差敏感余弦适合文本向量或用户行为向量这类重视方向、不重视绝对大小的数据。三种距离的向量化 MATLAB 写法供参考% 欧氏距离||a-b||^2 ||a||^2 ||b||^2 - 2*a*b dist_mat sqrt(sum(X.^2, 2) sum(X.^2, 2) - 2 * (X * X)); % 曼哈顿距离按列分别计算绝对差再求和用循环实现 % 余弦距离1 - 向量夹角的余弦值 dist_cos 1 - (X * X) ./ (sqrt(sum(X.^2, 2)) * sqrt(sum(X.^2, 2)));欧氏距离的展开式是一个很常用的技巧原理是二次展开避免双循环网上搜“距离矩阵向量化”几乎都会看到这一行。注意X * X得到的矩阵要小心对角线数值误差——因为浮点计算可能导致对角线出现 -1e-15 这种极小负数开方前先max(0, dist_sq)钳一下。4. 可视化模块画图部分的注释与配色设计4.1 二维散点图绘制聚类结果的可视化通常针对二维数据源码里的画图部分是把每个簇画成不同颜色散点中心点用大号星标标出。核心用scatter函数分散参数设置如下S点大小设为 30太小会被整体淹没太大在样本密集处会粘连。filled填充散点颜色比空心圆清晰很多。MarkerFaceAlpha设置透明度为 0.6样本重叠多时能看出密度层次。LineWidth和MarkerEdgeColor中心点要突出边框加粗、颜色用黑边。画图部分的主要代码结构类似这样% 画图部分 % 画出聚类结果不同簇用不同颜色区分medoid 用星号突出标记 figure; hold on; for c 1:k scatter(X(assignments c, 1), X(assignments c, 2), 30, colors(c, :), filled, MarkerFaceAlpha, 0.6); end scatter(X(medoid_idx, 1), X(medoid_idx, 2), 200, k, p, LineWidth, 2); title([k-medoids 聚类结果k , num2str(k)]); xlabel(特征 1); ylabel(特征 2); grid on; hold off;这里有个注释要特别提醒colors这个变量是代码内部定义的一个k x 3矩阵每行是一个 RGB 三元组。如果你的 k 大于 7直接用 MATLAB 默认的lines或parula色带抽取颜色即可否则手工写 10 种以上的 RGB 值很容易撞色。4.2 中心点突出显示与图例画中心点有个细节星标p只有 200 point size 时才足够醒目但也会遮挡周围样本。实测最好的做法是分两层画——先画一个带透明度的浅色大圆作为光环再在上面画实心星标。这样既突出中心又不挡信息% 画中心点光环半透明的大空心圆让星标更突出且不遮挡附近样本 scatter(X(medoid_idx, 1), X(medoid_idx, 2), 400, colors(c, :), o, LineWidth, 1.5, MarkerFaceAlpha, 0.15); scatter(X(medoid_idx, 1), X(medoid_idx, 2), 200, k, p, LineWidth, 2);图例方面如果是人工标注簇名称的场景可以直接用legend加字符串元胞数组。但聚类结果通常没有先验名称建议在标题里写清楚是第几次运行、k 值和距离度量图例反而可以省略。可视化最重要的目标是让读者一眼判断“簇分得开吗”“边界在哪”而不是精确读出每个点的编号。高维数据没法直接画二维散点这时候可以在画图模块前加一步 PCA 降维保留前两主成分再画注释里记得说明这属于“降维可视化不等同于聚类输入”。5. 常见问题排查与调试笔记5.1 问题速查表我在实际使用中把这些年跑 k-medoids 遇到过的问题做成速查表基本覆盖了新手会踩的坑现象原因解决方案中文注释显示为乱码MATLAB 编辑器的编码和源文件不一致用 MATLAB 的“打开”重新选择 UTF-8 编码或命令行执行edit重开文件readmatrix报错无法读取文件路径不对或文件在 MATLAB 的搜索路径之外先用which file_name检查是否在路径内或把cd切到数据目录聚类结果每次都不同初始化随机没有固定随机种子顶部加rng(42)并注释说明该行作用某个簇为空初始中心选得太近或数据分布极不均匀换最大最小距离初始化或检查标准化是否遗漏画图时中文标签乱码MATLAB 字体不支持中文字符在代码里加set(0,DefaultAxesFontName,Microsoft YaHei)距离矩阵维度溢出数据量过大导致矩阵内存不够用pdist2分块计算或把样本量降到一万以内测试结果始终是一坨没有区分度特征量纲不一致检查是否做了标准化避免某列数值主导距离迭代次数过多收敛慢距离计算用了循环而非向量化把内层循环改pdist2sum的向量化写法5.2 防止局部最优的三个操作习惯PAM 是局部搜索算法理论上不保证全局最优实操层面有三个有效习惯第一个是多次随机重启。同一个 k 值多跑 20 次每次用不同随机种子记录目标函数总距离和最小的那一次作为最终结果。代价是耗时变成 20 倍但聚类质量在大部分数据集上明显稳定。第二个是“先用 k-means 找中心位置再替换成最近的真实点”作为初始化。具体做法是把 k-means 跑一遍得到 k 个均值中心然后每个均值找离它最近的原始样本作为 k-medoids 的初始中心。这个方案在 sklearn 的KMedoids里叫k-medoids实测收敛一致性大幅提高。第三个是兜底策略朴素地检查每次迭代是否有中心点重复。如果第 i 个和第 j 个簇的中心相同那就是初始化失败或者 k 选大了果断换初始化重新跑别让代码带着明显错误继续跑完所有迭代。5.3 中文注释和显示乱码的处理这份源码最大的特色就是中文注释所以乱码问题必须单独讲。MATLAB 的编辑器默认按本地编码通常是 GBK打开文件但源码文件保存成 UTF-8 时中文字符在旧版本里就可能显示成乱码。好用的处理步骤是先用记事本打开源码文件另存为“带有 BOM 的 UTF-8”编码再用 MATLAB 打开重试。更省事的办法是在 MATLAB 命令行执行feature(DefaultCharacterSet, UTF-8)切换默认编码但这只对当前会话有效下次启动要重新执行。另外建议在源码这样写注释% 功能对输入样本矩阵 X 执行 k-medoids 聚类 % 输入X —— n x d 样本矩阵n 为样本数d 为特征维度 % k —— 聚类簇数须为小于 n 的正整数 % 输出medoid_idx —— k x 1 向量中心样本的行索引 % assignments —— n x 1 向量每个样本所属的簇编号 % total_cost —— 数值最终簇内总相异度越小说明聚类越紧密把人话写进注释几个月后回看既有“现在在干什么”也有“为什么这么干”。这也是我改别人代码时最看重的一点——不是注释越多越好而是关键决策处必须写清动机。总结一句这份代码值得留作基础模板把距离度量参数扩展到三种初始化换成最大最小距离迭代过程改成向量化再做中文注释就是一个能放进项目文档和课程报告的完整方案。