ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

Matlab实现SOM聚类:从Excel读取到可视化完整指南

Matlab实现SOM聚类:从Excel读取到可视化完整指南 很多人一提到聚类脑子里第一个蹦出来的就是K-means。但真拿一批业务数据上手做几轮之后你会很快发现K-means的几个限制K值你得事先拍脑袋初始中心一换结果就跟着漂而且它对数据形态的假设也比较理想化。这时候自组织特征映射SOM反倒是一个更顺手的选择。这篇文章想聊的就是我用Matlab实现的一个SOM数据聚类程序。数据从Excel读入跑完直接输出每个样本的聚类标签整个过程不依赖任何商业统计软件也不需要额外的数据库只要电脑上装了Matlab就能复现。SOM最大的特点是无监督、不需要预先指定类别个数而且训练完还能把高维数据铺到一个二维网格上看分布。很多做探索性数据分析、客户分群、故障诊断、图像分割的人都用它来先“摸一摸”数据的大致结构再决定下一步怎么做。这篇文章不只给代码还会把参数怎么调、Excel怎么整理、U-Matrix怎么看这些实操里绕不开的细节一并讲清楚。适合正在做课程作业、毕业论文或者工作中刚接到一批数据还不知道从哪儿下手的读者。1. 为什么聚类要选SOM原理和选型思路1.1 K-means的三个痛点先说K-means。它确实简单收敛快代码到处都能抄但用它的前提是你要知道自己想聚成几类。现实中大部分时候数据里到底有几类你根本不知道。比如电商用户分群你可能预期分4类跑出来之后发现第2类里还能拆出两类可K-means一旦定了K再怎么调初始中心也很难跳出那个固定框架。第二个痛点是初始中心敏感。K-means用随机初始化不同的随机种子得到的局部最优解可能差别很大。你上午跑出一版结果下午换台电脑又跑出一版老板问你怎么两边对不上你只能解释“随机种子不同”。这在实际交付里很尴尬。第三个痛点K-means默认类别是凸的、球形簇。数据分布如果是长条形、月牙形、嵌套形K-means基本无能为力。相比之下SOM没有这些先验假设它靠神经元之间的拓扑关系去“贴”数据的形态边界不规则的簇也能表达出来。1.2 SOM的本质把高维数据铺到二维网格上SOM的思想很有意思。它把若干个神经元排成一个网格常见的是二维矩形网格每个神经元都有一组权重向量维度跟样本特征数一致。训练开始后每条样本逐个送进来跟所有神经元算距离距离最小的那个神经元获胜叫做BMUBest Matching Unit最佳匹配单元。但SOM不只是让赢家学习。它以BMU为中心把周围邻域内的神经元也拉过来一起更新离得越近的神经元权重被拽得越厉害离得远的只是轻微动一下。这个“协作”过程是关键经过反复训练网格上相邻神经元的权重向量会变得越来越像。最终结果就是原本高维空间里相近的样本会落到网格上相近的神经元附近高维空间的簇在网格上也会自然形成一块连通的区域。这就是SOM的拓扑保持特性。它相当于把高维分布“摊平”到二维平面上聚类结构肉眼可见。这也是为什么很多人喜欢拿SOM做探索性分析你不用提前设定类别数训练完看一眼U-Matrix就能判断数据大致有几群。1.3 为什么用Matlab、用Excel选Matlab不是因为它比Python高级而是因为矩阵运算确实省事。SOM的核心操作是反复计算样本与权重向量的距离这在Matlab里就是一条矩阵减法加sum。而且Matlab对Excel的读写支持非常成熟readtable一行就能把表读进来writetable一行就能把结果写回去对大多数习惯用Excel管理数据的人来说交互成本很低。再者很多高校实验室和传统工科背景下Matlab就是默认工具。神经网络工具箱里虽然也提供了selforgmap函数但装上工具箱的人未必多。我用的是一个几十行的纯手写版本不依赖任何工具箱换个版本也能跑这在实际交付时能省掉很多环境兼容的麻烦。2. 数据准备与参数设计先想清楚再动手2.1 Excel数据的格式约定与读取我把话放前面SOM程序本身不难80%的问题都出在数据格式上。程序里的数据读取逻辑很简单但Excel表格如果不按约定整理后面跑出来一堆NaN或报错调试起来相当费劲。先约定格式。第一行必须放表头每个特征占一列每行一个样本。例如一列用户ID后面几列是消费金额、活跃天数、平均客单价等数值特征。表头建议用英文或拼音用中文也能读但会出现变量名带中文的问题处理起来反而不顺手。Excel文件本身建议保存为xlsx格式别用xls。聚类时只让数值特征参与训练ID、姓名、日期这些列不能进训练矩阵。你在读数据时可以把这些列单独留着最后输出结果时再拼回去这样聚类标签能和原始记录对上号。读取时用readtable就够了% 读取Excel data readtable(data.xlsx); % 查看表结构确认每列类型 summary(data)如果表里有文本列data{:, :}会把整表变成cell之后再算距离就会报错。一种稳妥的办法是按列号取数值列比如% 假设第1列是ID第2~6列是特征 X data{:, 2:6};2.2 特征归一化min-max还是z-scoreSOM用的是欧氏距离而欧氏距离对量纲极其敏感。拿用户数据举例消费金额可能是几千活跃天数可能是几十如果不归一化“消费金额”这一列会在距离计算里一票否决其他特征模型基本就退化成了只看单维度的排序问题。我通常先用min-max归一化把每个特征压缩到[0,1]区间% 逐列min-max归一化 X_min min(X); X_max max(X); X_norm (X - X_min) ./ (X_max - X_min);min-max的缺点是受离群点影响大。如果某列有一个极端大值其他正常值会被压缩到很窄的区间细节全丢。这时候可以考虑z-score标准化% z-score标准化 X_norm (X - mean(X)) ./ std(X);z-score的好处是不怕离群点但标准化之后数据不再落在[0,1]对SOM影响不大因为距离计算是无偏的。我自己的习惯是先看每列的最大最小值明显有极端值的用z-score分布相对均匀的用min-max。另外要注意如果某一列是常数min-max归一化会出现除以0的问题算出来全是NaN。我实际踩过这个坑后来在代码里加了判断% 处理常数列统一置0 rangeX X_max - X_min; rangeX(rangeX 0) 1; X_norm (X - X_min) ./ rangeX;2.3 神经元网格多大才合适网格大小直接决定聚类分辨率。网格太小Topology结构容易被压扁原本不同的簇可能被挤到同一块区域网格太大会出现大量没有样本落到上面的“空神经元”聚类结果太碎反而不好解读。我一般按样本量来估网格规模。经验范围是让神经元数量落在5sqrt(n)到10sqrt(n)之间n是样本条数。举个例子300条样本sqrt(300)约等于175倍是8510倍是170所以5x5到12x12的网格都在合理区间。实际跑的时候我会从偏小的网格开始比如5x5先看U-Matrix的整体轮廓再逐步放大。下面这张表是我在几个数据集上试出来的参考范围不是严格公式但可以當做起点样本量建议网格说明小于2005x5或6x6网格大了容易空转200~10008x8或10x10最常用区间1000~500010x10~15x15视特征复杂度调整大于500015x15或20x20配合批处理训练网格形状也可以选六边形拓扑比矩形网格更自然但在Matlab里手写复杂度会高一些。初期用矩形网格就够用了。2.4 学习率、邻域半径和迭代次数三个参数怎么调才不翻车SOM训练里有三个关键参数学习率、邻域半径、迭代次数。这三个参数不需要特别精确但方向不能错。学习率控制每次权重更新的步长。初始值我习惯设在0.5到0.8之间随着训练轮数逐渐衰减到接近0.01。衰减太快神经元还没来得及展开就定型了衰减太慢后期权重还在大幅摆动训练曲线会来回震荡。线性衰减是够用的方案learnRate initLearnRate * (1 - epoch / totalEpochs);邻域半径控制每次更新时“拉拢”多大的范围。初始半径我一般设为网格最大边长的一半比如10x10的网格初始半径是5。训练过程中半径也从大往小缩最后稳定在1左右。半径大时神经元之间互相牵连整体结构快速展开半径小时每个神经元才开始精细地贴向各自对应的样本。迭代次数就比较直观了。每次把所有样本完整过一遍训练叫做一个epoch。200个epoch是一个比较稳的起点样本量少、特征维度低时可以降到100。判断是否收敛有个土办法每个epoch结束后计算所有样本到各自BMU的平均距离如果这个值后面基本不再下降就说明模型已经稳定了。3. 核心实现从Excel到聚类标签的完整代码3.1 数据读取与归一化的Matlab写法先给出一段能直接跑的读取和预处理代码。这里我假设Excel第一列是样本ID之类的文本信息第2到第6列是数值特征% SOM聚类 - 数据读取与预处理 clc; clear; close all; % 1. 读取Excel data readtable(data.xlsx); % 2. 取出数值特征假设第2~6列为特征 X data{:, 2:6}; % 如果表格里没有非数值列可直接用 X data{:, :}; % 3. 检查是否有缺失值 fprintf(缺失值数量%d\n, sum(isnan(X(:)))); if any(isnan(X(:))) error(数据存在NaN请先处理缺失值); end % 4. min-max归一化带常数列保护 X_min min(X, [], 1); X_max max(X, [], 1); rangeX X_max - X_min; rangeX(rangeX 0) 1; X_norm (X - X_min) ./ rangeX; fprintf(数据读取完毕%d条样本%d个特征\n, size(X, 1), size(X, 2));运行到这一步X_norm就是可以用来训练的矩阵。如果你的原始数据列数不同把列号换成实际的位置就行。3.2 SOM训练核心思路竞争、协作、更新训练循环只有三个动作竞争、协作、更新。竞争就是找BMU协作是计算邻域内每个神经元的影响权重更新是把这些神经元的权重向量朝当前样本拉近。为了方便计算模型权重W用一个二维矩阵存放每一行是一个神经元的权重向量。同时维护一个网格坐标矩阵gridPos记录每个神经元在网格上的位置。这样可以避免用循环嵌套去处理二维索引代码更简洁也更容易套用到不同尺寸的网格。训练循环的伪逻辑是这样的外层循环控制epoch内层循环遍历每条样本。每个epoch先把样本顺序随机打乱再逐条送入网络寻找BMU并更新权重。打乱顺序能避免样本本身排列顺序影响训练结果这个小细节我后面还会再提。对应代码如下% SOM训练核心参数 mapSize [10 10]; % 网格尺寸 10x10 totalEpochs 200; % 迭代轮数 initLearnRate 0.5; % 初始学习率 initRadius max(mapSize) / 2; % 初始邻域半径 n size(X_norm, 1); % 样本数 d size(X_norm, 2); % 特征维度 % 初始化网格坐标和权重 [I, J] ind2sub(mapSize, (1:prod(mapSize))); gridPos [I, J]; W rand(prod(mapSize), d); % 权重矩阵每个神经元一行 % 训练 for epoch 1:totalEpochs % 学习率与半径线性衰减 learnRate initLearnRate * (1 - epoch / totalEpochs); radius initRadius * (1 - epoch / totalEpochs) 0.5; % 打乱样本顺序 order randperm(n); for k order x X_norm(k, :); % 竞争找BMU diff W - x; dist2 sum(diff .^ 2, 2); [~, bmu] min(dist2); % 协作计算邻域高斯权重 diffPos gridPos - gridPos(bmu, :); distGrid sqrt(sum(diffPos .^ 2, 2)); influence exp(-distGrid .^ 2 / (2 * radius ^ 2)); % 更新所有神经元按影响力向样本移动 W W learnRate .* influence .* (x - W); end % 每50轮输出一次平均距离观察收敛 if mod(epoch, 50) 0 dist2 sum((W - X_norm).^2, 1); avgDist mean(sqrt(dist2)); fprintf(Epoch %d, 平均距离 %.4f\n, epoch, avgDist); end end这段代码里的邻域计算是整个SOM的精华。influence是一个向量值落在0到1之间BMU所在位置影响最大离得远的位置接近0。权重更新公式写成W W learnRate * influence .* (x - W)相当于每个神经元都往样本x方向走一小步但走的幅度由学习率和邻域影响共同决定。循环跑完之后W就是训练好的SOM模型。3.3 聚类标签提取从BMU到最终类别训练结束后每条样本都要确定自己的归属。对每一条样本在全部分布的权重矩阵里找距离最小的神经元这个神经元的编号就是样本的BMU索引% 为每条样本分配BMU sampleBmu zeros(n, 1); for k 1:n diff W - X_norm(k, :); dist2 sum(diff .^ 2, 2); [~, sampleBmu(k)] min(dist2); end到这里已经有了每个样本的神经元标签。但直接用神经元编号当聚类标签有个问题10x10网格有100个神经元但实际数据可能只有5~6类如果直接把100个神经元编号输出成标签类别太碎不实用。所以常规做法是对神经元的权重向量做二次聚类。用K-means对W做一次聚类指定你最终想要的类别数K然后把每个神经元归到一个大类里最后把样本的BMU映射到对应大类上。这一步相当于先用SOM把数据降维到网格再用K-means在网格上划大类能获得比直接K-means更稳定的结果% 二次聚类将神经元分为K类 K 3; neuronLabel kmeans(W, K, Replicates, 10); % 映射到样本 label neuronLabel(sampleBmu);kmeans里的Replicates设为10表示用10组不同的初始中心跑10次选最优可以减轻K-means随机性的影响。这里的K需要你根据业务或U-Matrix的观察结果来定。3.4 U-Matrix可视化让聚类边界自己显形U-Matrix是SOM最直观的可视化方法。它计算每个神经元和相邻神经元之间的平均距离把结果映射成灰度图像。距离大的区域在图上呈现亮色说明相邻神经元之间差异大这些都是聚类边界的候选位置距离小的区域呈现暗色说明神经元之间相似度高对应数据稠密的核心区。U-Matrix的实现不复杂一个双循环就能搞定。为了可读性我先把权重矩阵重新排列成网格形状再对每个神经元取四邻域求平均距离% 计算U-Matrix Wmap reshape(W, [mapSize, d]); uMat zeros(mapSize); for i 1:mapSize(1) for j 1:mapSize(2) neighborDist []; % 四邻域上、下、左、右 for di [-1 0 1 0] dj 0; % 这里用一个更直观的邻域遍历 continue; end % 上面是占位实际用下面的写法更清晰 neighborDist []; for di -1:1 for dj -1:1 if abs(di) abs(dj) ~ 1 continue; end ni i di; nj j dj; if ni 1 ni mapSize(1) nj 1 nj mapSize(2) dW squeeze(Wmap(i, j, :)) - squeeze(Wmap(ni, nj, :)); neighborDist(end1) sqrt(sum(dW .^ 2)); %#okSAGROW end end end uMat(i, j) mean(neighborDist); end end % 显示U-Matrix figure; imagesc(uMat); colormap(bone); colorbar; title(U-Matrix);看U-Matrix有个技巧亮色的带状区域就是潜在的分类边界。如果图像上能看到两三块清晰的暗色区域被亮色带隔开那数据大概就是两三类的结构。这时候再去设二次聚类的K值心里就有底了。3.5 完整程序整合把结果写回Excel把前面所有代码串起来再加一段结果输出逻辑这就是一个能直接交付的完整程序。输出部分我用writetable写回Excel把原始数据和聚类标签放在同一张表里% 结果输出 data.Label label; data.BMUNeuron sampleBmu; writetable(data, clustered_result.xlsx); % 统计每个类别的样本量 tabulate(label)这里直接修改了data表把Label列和BMUNeuron列加进去。如果原始数据第一列是ID导出后也能对上号。我自己在项目里执行完这一句之后还会顺手用一句histogram画一下各类别的数量分布确认没有哪个类别样本量少到离谱。这里要提醒一个版本问题writetable是较新版本Matlab才有的函数如果你的环境比较老可以把这行换成xlswrite(clustered_result.xlsx, data)但旧函数不支持表格类型需要先把数据转成cell数组。如果装的是近几年版本直接用writetable就行。4. 常见问题与排查技巧实录4.1 Excel读不进来问题出在哪我遇到最多的报错是readtable读取失败错误信息五花八门但原因通常就三类文件路径不对、Excel文件被占用、xls格式兼容问题。文件路径这个最容易踩。Matlab在Windows下的路径分隔符是反斜杠复制文件夹地址时经常带着中文或空格readtable解析起来容易出问题。我的习惯是把Excel文件和脚本放在同一个目录下直接用文件名读取不写绝对路径。如果一定要读其他路径先把当前目录切过去cd(D:\project\data); data readtable(data.xlsx);还有一次同事发来的文件实际上是CSV文件改了扩展名readtable虽然能读但中文表头编码全乱。遇到这种情况我会先用Excel打开文件“另存为”重新保存成标准xlsx再做后续处理。4.2 训练结果每次跑都不一样SOM虽然比K-means稳定但权重矩阵是随机初始化的每次运行结果会有细微差别。如果只是细微扰动不影响聚类结论这可以接受。但如果你希望程序结果可以复现就在脚本最前面固定随机种子rng(42);固定种子之后每次运行的结果是完全一致的。输出报告、做论文复现的时候这一步别漏。另外注意训练循环里我用了randperm打乱顺序这个操作本身也受随机种子影响。固定了rng之后整体训练过程就是确定性的。4.3 训练不收敛平均距离反复震荡训练过程中如果fprintf打印出来的平均距离没有下降趋势而是忽高忽低多半是学习率衰减得太慢了。到训练后期学习率应该已经降得很低如果初值设得太大且衰减公式不对后期权重还在大幅摆动距离自然稳不下来。排查方法是把平均距离的曲线画出来% 记录每轮平均距离画曲线 avgDistHistory(epoch) mean(sqrt(dist2));如果曲线总体下降但尾部有小幅波动这是正常的。如果曲线整个像锯齿就把initLearnRate调小到0.3或者把totalEpochs加长到300。实测下来200轮对大多数标准数据集已经够用但高维稀疏数据可能需要更多轮次。4.4 出现大量空神经元网格上有些神经元从头到尾没当过任何样本的BMU这就是空神经元。造成空神经元的原因通常是网格设得太大或者数据本身类别数远小于神经元数。空神经元本身不影响聚类结果但会让U-Matrix上出现大片的暗区干扰判断。处理办法很简单一是缩小网格二是做二次聚类时把空神经元也归入最近的簇。我一般用后者因为网格大小还要兼顾U-Matrix的分辨率不想为了消除空神经元而牺牲可视化细节。4.5 数据量大时训练太慢SOM逐样本更新复杂度天然是O(epochs * n * m)n是样本数m是神经元数。几千条样本跑起来还好到了几万条就明显迟钝了。我的优化策略有两种一是随机采样一部分代表性样本参与训练得到权重后再把所有样本映射到BMU二是把外层epoch循环减少到50通过增大初始学习率快速展开拓扑结构。对于Excel里几十万行的大表我其实不推荐直接上SOM。先把数据抽样到5000条以内做探索性分析确定聚类模式和K值再对大样本跑一次K-means这种组合拳在实际项目里的效果很好。4.6 问题速查表问题现象可能原因解决方案readtable报错路径有中文/文件被占用把文件放在脚本目录直接使用文件名读取训练结果NaN原始数据有缺失值提前检查isnan并处理缺失值权重全是NaN常数列归一化除以0给rangeX加保护除0位置置1结果每次不一致随机初始化未固定脚本开头加rng(42)平均距离震荡学习率衰减过慢降低initLearnRate或延长epochsU-Matrix边界模糊网格太大且训练不足缩小网格、增加迭代次数空神经元过多网格尺寸不合理缩小网格或二次聚类时归并大数据训练慢样本量过大抽样训练或先减少epochs这些坑我基本都踩过一遍。SOM本身并不复杂复杂的是数据准备和参数之间的配合。真要说有什么经验值得强调那就是别一上来就追求完美的大网格先用5x5跑一版看看U-Matrix的轮廓再逐步放大网格、微调参数。我实际做过的几个项目里SOM最有价值的一部分不是输出标签本身而是那张U-Matrix图——对着它讲数据结构和分类结果比如你在向业务方或导师解释“这批数据为什么分成这几类”时一张拓扑保留的灰度图比一堆聚类中心坐标有说服力得多。最后再分享一个小技巧给SOM训练前的归一化参数min、max、mean、std单独存一份到Excel或mat文件里。之后来新样本时用同一套参数做归一化再映射到训练好的SOM网格上就能直接预测新样本属于哪个簇。这个流程把SOM从“离线的聚类分析工具”升级成了“可复用的分类器”在很多实际业务里都是加分项。
返回列表