简介这份资料面向图像处理入门者与Matlab算法实践者围绕K-means聚类在图像特征分割中的应用展开帮助读者理解无监督聚类如何将像素按RGB或灰度特征划分为若干类别并解决分割结果受初始质心影响、K值选取等实际问题。压缩包共10个文件以7个m脚本为核心另含2个png与1个jpg示例图像整体约487KB脚本覆盖初始化、距离计算、聚类主流程与结果验证等环节图像文件可直接用于运行与对比实验。目前已有1381人学习下载说明该案例在同类资源中具备一定参考价值。读者可借助完整代码复现从数据预处理、质心迭代更新到分割图像输出的全过程并通过调整K值与初始质心观察分割差异进而掌握参数优化与结果评估思路适合作为课程设计、算法练习或图像分割项目的起步参考。1. 从一张过曝的叶片照片说起K-means 图像特征分割到底在切什么去年帮一个做植物表型的朋友处理一批温室叶片图像背景是白色传送带叶片本身有反光。用固定阈值分割反光区域直接被判成背景叶面积算出来少了将近三成。换成 K-means 聚类算法做图像特征分割把 RGB 三个通道当成三维特征空间里的坐标让像素自己抱团反光区域因为颜色接近高亮白被单独聚成一类反而好剥离了。这件事让我意识到K-means 在图像分割里的价值不在于它多先进而在于它不需要你事先知道阈值是多少——你只需要告诉它我想分成几堆。这篇要讲的就是基于 K-means 聚类算法的图像特征分割配套 Matlab 代码怎么落地。适合两类人一类是刚接触图像处理、想找一个能跑通、能改参数、能看懂每一步在干什么的入门方案另一类是做过阈值分割或边缘检测但遇到光照不均、背景复杂就翻车的从业者。核心问题只有一个怎么把一张图里的像素按颜色或纹理特征自动分成有意义的几组并且让这个分组结果对后续的测量、识别、压缩真正有用。Matlab 在这个方向上的优势很直接矩阵运算天然适配 K-means 的向量化实现Image Processing Toolbox 里的imread、rgb2lab、imsegkmeans能省掉大量底层读写和色彩空间转换的代码。但我不打算只调一个imsegkmeans就交差那样你遇到内存溢出或者聚类结果抖动时完全不知道从哪下手。下面从原理到代码到踩坑一层层拆。2. K-means 做图像分割的原理与 Matlab 选型为什么不是 GMM 也不是超像素2.1 把像素当成数据点特征空间里的距离度量一张 M×N 的 RGB 图像展开后是 M×N 个三维向量每个向量就是该像素的 (R, G, B) 值。K-means 要做的事是在这个三维空间里找 K 个中心点使得每个像素到离它最近的中心点的距离平方和最小。这个距离通常用欧氏距离但如果你把图像转到 Lab 色彩空间L 通道代表亮度a、b 代表色度欧氏距离就更接近人眼感知的色差。迭代过程只有两步反复执行到收敛分配步每个像素归到最近的中心更新步每个中心改成它名下所有像素的均值。收敛条件一般是中心点移动量小于某个容差或者达到最大迭代次数。这里有个容易被忽略的点K-means 对初始中心敏感不同的初始化可能收敛到不同的局部最优。Matlab 的imsegkmeans默认用 k-means 初始化比纯随机选点稳得多这也是我一般不自己手写随机初始化的原因。特征不限于 RGB。你可以把像素坐标 (x, y) 拼到特征向量后面这样空间上靠近的像素更容易被分到一类分割结果更连续。代价是坐标的量纲和颜色量纲不一致需要归一化否则坐标项会主导距离计算。我通常把坐标缩放到 [0, 1]颜色也缩放到 [0, 1]再拼起来。2.2 Matlab 里三种实现路径的取舍第一种是直接调imsegkmeans一行代码出标签图。适合快速验证 K 值选得对不对但它内部做了降采样加速对小目标不友好而且你拿不到中间迭代过程。第二种是kmeans函数把像素矩阵传进去指定Distance为sqeuclideanReplicates设为 3 到 5 次重复取最优。这个函数返回中心点和每个像素的归属索引你能完全控制特征矩阵的构造。缺点是它把每个像素当独立样本大图会吃内存。一张 4000×3000 的图有 1200 万像素kmeans直接跑大概率内存爆掉。第三种是自己写迭代循环配合blockproc或分块处理。这是我在处理大幅遥感影像时的做法把图切成 512×512 的块每块单独聚类块之间用重叠区域做中心点对齐。代码量大但内存可控而且能并行。选型建议教学和中小图1000×1000 以内用第二种能看清每一步快速出图用第一种工业级大图用第三种。下面给的代码以第二种为主因为它最能暴露参数和坑。2.3 特征构造RGB、Lab 还是灰度加纹理只做颜色分割RGB 够用。但 RGB 三个通道相关性高欧氏距离在 RGB 空间里和人眼感知偏差大。转到 Lab 后L 通道单独控制亮度a、b 控制色度分割结果对阴影更鲁棒。Matlab 转换用rgb2lab注意输入必须是 double 或 uint8uint8 会自动按 sRGB 规则转。如果图像里不同区域的色差很小但纹理差异大纯颜色特征会失效。这时候要加纹理特征常用的是局部二值模式LBP或者灰度共生矩阵GLCM的对比度、能量。做法是先用extractLBPFeatures或graycomatrix算出每个像素邻域的纹理描述子和颜色特征拼在一起再做 K-means。特征维度上去了K 值也要相应调整否则聚类会过于碎片化。提示特征拼接前一定要归一化。颜色值范围是 0 到 255 或 0 到 1纹理特征范围可能是 0 到几十不归一化的话量纲大的特征会完全主导距离。3. 用 Matlab 跑通 K-means 图像分割的最小代码从读图到标签着色3.1 完整可复现脚本与逐段说明下面这段代码在 Matlab R2020b 及以上版本可以直接运行不需要额外工具箱只依赖 Image Processing Toolbox 的rgb2lab和label2rgb。我特意把特征构造、聚类、后处理拆成独立段落方便你替换其中任何一步。% kmeans_seg_demo.m % 基于 K-means 的 RGB 图像分割最小示例 clear; clc; close all; % 1. 读图并降采样控制内存 img imread(leaf_sample.jpg); % 替换为你的图片路径 img imresize(img, 0.5); % 长宽各减半像素数降为 1/4 [h, w, ~] size(img); % 2. 构造特征矩阵Lab 颜色 归一化坐标 lab rgb2lab(img); % 转换到 Lab 空间 lab double(lab); L lab(:,:,1) / 100; % L 范围约 0-100 A (lab(:,:,2) 128) / 255; % a 范围约 -128-127 B (lab(:,:,3) 128) / 255; % b 范围约 -128-127 [X, Y] meshgrid(1:w, 1:h); X (X - 1) / (w - 1); % 坐标归一化到 0-1 Y (Y - 1) / (h - 1); % 特征向量[L, a, b, X, Y]坐标权重用 wxy 控制 wxy 0.3; % 坐标项权重0 表示纯颜色 features [L(:), A(:), B(:), wxy*X(:), wxy*Y(:)]; % 3. K-means 聚类重复 3 次取最优 K 4; % 聚类数先试 4 rng(42); % 固定随机种子保证可复现 [idx, C] kmeans(features, K, ... Distance, sqeuclidean, ... Replicates, 3, ... MaxIter, 200, ... Display, final); % 4. 标签图还原为二维图像 labelMap reshape(idx, [h, w]); % 5. 用聚类中心颜色着色便于目视检查 centerColors C(:, 1:3); % 取 Lab 中心的前三维 centerColors centerColors .* [100, 255, 255] - [0, 128, 128]; centerColors lab2rgb(centerColors); % 转回 RGB 用于显示 centerColors max(0, min(1, centerColors)); segImg label2rgb(labelMap, centerColors, k, shuffle); % 6. 显示原图、标签图、着色图 figure(Position, [100 100 1200 400]); subplot(1,3,1); imshow(img); title(原图); subplot(1,3,2); imagesc(labelMap); axis image off; title(标签图); subplot(1,3,3); imshow(segImg); title(聚类着色);这段代码的逻辑链条是读图降采样控制规模转 Lab 让颜色距离更符合感知拼坐标让空间连续kmeans做聚类最后把一维索引还原成二维标签图并着色。rng(42)这行别省K-means 有随机初始化不固定种子的话每次跑出来的标签编号和边界都可能不一样调试时你会怀疑人生。3.2 关键参数怎么调K 值、坐标权重、重复次数K 值是第一变量。K 太小不同区域被合并K 太大同一区域碎成多块。我一般从 3 开始试每次加 1看标签图里目标区域是否被完整切出。定量判断可以用轮廓系数Matlab 没有内置但可以自己算对每个像素算它到同簇其他点的平均距离 a 和到最近异簇的平均距离 b轮廓系数为 (b-a)/max(a,b)全图取均值。值越接近 1 越好但计算量大我通常只在 K 候选值之间比较时用而且会再降采样一次。坐标权重wxy控制空间连续性。设成 0纯颜色分割同色但分离的区域会被归为一类设成 0.3 到 0.5空间上靠近的像素更容易同簇边界更平滑但可能把颜色不同但相邻的区域强行合并。我的经验是目标区域颜色区分度高时wxy取 0.1 到 0.2颜色接近、需要靠空间位置区分时取 0.5 以上。Replicates设为 3 到 5表示用不同初始化跑多次取距离和最小的那次。次数越多越稳但耗时线性增长。MaxIter默认 100我一般设 200因为图像像素多100 次有时还没收敛。Display设为final能看到最终迭代次数和总距离和如果迭代次数顶到MaxIter说明还没收敛要加大。3.3 分割结果的后处理去噪、填洞、边界提取K-means 输出的标签图往往有孤立小区域像椒盐噪声。用medfilt2对标签图做中值滤波窗口 3×3 或 5×5能去掉大部分孤立点。注意medfilt2要求输入是二维数值矩阵标签图正好符合。如果目标区域内部有空洞用imfill填洞。先对每个标签做二值化bw labelMap k再bw2 imfill(bw, holes)最后把填洞后的区域合并回标签图。这一步在细胞分割、缺陷检测里特别常用因为聚类可能把目标内部的高光或阴影单独聚成另一类形成空洞。边界提取用bwperim或者对标签图求梯度。bwperim返回二值边界叠加到原图上就能看到分割轮廓。如果要亚像素边界用edge配合 Canny但那是另一条路和 K-means 是互补关系不是替代。4. 避坑与排查K-means 图像分割翻车的五个真实场景4.1 现象跑大图直接内存溢出Matlab 卡死原因kmeans函数把每个像素当独立样本特征矩阵大小是 像素数 × 特征维度。一张 4000×3000 的 RGB 图像素数 1200 万特征维度 5矩阵就是 1200 万 × 5 的 double约 480 MB加上kmeans内部的距离矩阵和临时变量实际内存占用可能是这个数的好几倍。如果图更大或者特征维度更高爆内存是必然的。解决先降采样。imresize到长边 1000 以内像素数降到百万级内存和速度都能接受。如果必须保留原分辨率用分块处理blockproc把图切成 512×512 的块每块单独聚类块之间重叠 32 像素用重叠区域的中心点做对齐。或者改用imsegkmeans它内部做了降采样加速但小目标会丢适合大区域分割。4.2 现象每次运行分割结果都不一样标签编号乱跳原因K-means 的初始中心是随机选的不同的初始化收敛到不同的局部最优。标签编号本身没有语义第 1 类这次是叶片下次可能是背景所以直接比较两次的标签图没有意义。解决固定随机种子rng(42)放在kmeans之前。这样每次初始中心一样结果可复现。但要注意固定种子只是让结果稳定不保证结果最优。如果换一张图种子不变但数据分布变了可能收敛到更差的解。所以固定种子用于调试和对比正式跑的时候还是用Replicates多次取最优。另外比较两次分割结果时不要直接比标签编号要比区域重叠度比如用 IoU 或者混淆矩阵。4.3 现象分割边界锯齿严重目标边缘像狗啃原因K-means 是逐像素独立分类不考虑邻域一致性所以边界处像素容易被分到不同类形成锯齿。坐标权重wxy设得太小也会加剧这个问题。解决加大wxy让空间坐标参与距离计算相邻像素更容易同簇。但wxy太大会导致过平滑小目标被吞掉。另一个办法是后处理对标签图做中值滤波或形态学开闭运算bwareaopen去掉小面积区域imclose连接断裂边界。如果对边界精度要求高可以在 K-means 之后接一个图割或者条件随机场做边界优化但那是另一个量级的复杂度了。4.4 现象K 值调来调去都不满意目标区域总是被切碎或合并原因K 值的选择没有唯一标准取决于图像内容和分割目的。同一张图分 3 类可能把叶片和背景分开分 5 类可能把叶片的高光、阴影、健康组织、病斑都分开。没有正确的 K只有适合当前任务的 K。解决先用肘部法则或轮廓系数缩小范围。肘部法则看聚类总距离随 K 变化的曲线拐点处的 K 通常是候选。轮廓系数更直接但计算量大建议在降采样后的图上算。然后结合目视检查把 K 候选值各跑一遍看哪个 K 下目标区域完整且内部一致。如果颜色特征不够加纹理特征如果空间连续性不够加坐标特征。K 值不是孤立的它和特征构造、权重设置一起决定结果。4.5 现象Lab 转换后数值范围不对聚类结果全黑或全白原因rgb2lab对 uint8 输入按 sRGB 规则转换L 范围 0 到 100a 和 b 范围约 -128 到 127。如果输入是 double 且范围 0 到 1rgb2lab会按线性 RGB 处理结果范围不同。很多人读图后直接im2double再转 Lab数值范围就乱了。解决读图后先确认数据类型。imread返回 uint8 的话直接传给rgb2lab不要先转 double。如果已经是 double确保范围是 0 到 1并且知道rgb2lab会按线性 RGB 处理结果和 sRGB 有差异。最稳的做法是统一转成 uint8 再转 Lab或者手动按公式归一化。转换后检查 L、a、b 的最小最大值L 应该在 0 到 100 附近a 和 b 在 -128 到 127 附近偏离太多说明转换有问题。5. 让 K-means 分割真正可用K 值自动搜索与分块大图处理5.1 用轮廓系数自动选 K代码与计算量控制轮廓系数在像素级算太慢我的做法是先对图像做超像素分割比如用superpixels函数生成 500 到 1000 个超像素每个超像素用平均颜色作为特征再在超像素级别算轮廓系数。这样计算量从百万级降到千级几秒钟就能出结果。% auto_select_k.m % 用超像素 轮廓系数自动搜索 K img imread(leaf_sample.jpg); img imresize(img, 0.3); lab rgb2lab(img); lab double(lab); % 生成超像素N 为超像素个数 N 800; [Lsp, Nsp] superpixels(lab, N); % 计算每个超像素的平均 Lab 颜色 spFeatures zeros(Nsp, 3); for i 1:Nsp mask Lsp i; spFeatures(i, :) squeeze(mean(mean(lab .* mask, 1), 2)) ./ sum(mask(:)); end spFeatures normalize(spFeatures); % 归一化 % 在 K 2 到 8 之间搜索 Krange 2:8; silScore zeros(size(Krange)); for ki 1:length(Krange) K Krange(ki); rng(42); idx kmeans(spFeatures, K, Replicates, 3, Display, off); s silhouette(spFeatures, idx, sqeuclidean); silScore(ki) mean(s); end % 找轮廓系数最大的 K [~, bestIdx] max(silScore); bestK Krange(bestIdx); fprintf(推荐 K %d轮廓系数 %.3f\n, bestK, silScore(bestIdx)); % 画曲线 figure; plot(Krange, silScore, -o); xlabel(K); ylabel(平均轮廓系数); title(K 值选择曲线); grid on;这段代码的关键是超像素降维。superpixels把颜色相近且空间相邻的像素合并成小块每个小块用平均颜色代表既保留了空间信息又把样本数从百万降到几百。silhouette函数在超像素特征上算速度快结果稳定。注意normalize默认按列归一化对特征矩阵正好。5.2 分块处理大幅图像重叠对齐与中心点合并大图分块聚类的难点是块与块之间的标签不一致。块 A 的第 1 类可能是绿色块 B 的第 1 类可能是白色直接拼起来边界会错乱。我的做法是每块聚类后保存聚类中心相邻块的重叠区域里用中心点距离做匹配把块 B 的中心点映射到块 A 的标签体系然后按映射关系重编号块 B 的标签。% block_kmeans.m % 分块 K-means重叠区域对齐中心点 img imread(large_image.jpg); img imresize(img, 0.5); [h, w, ~] size(img); lab double(rgb2lab(img)); blockSize 512; overlap 64; step blockSize - overlap; K 4; % 初始化全局标签图和中心点存储 globalLabel zeros(h, w); allCenters []; for r 1:step:h for c 1:step:w r2 min(r blockSize - 1, h); c2 min(c blockSize - 1, w); block lab(r:r2, c:c2, :); [bh, bw, ~] size(block); feats reshape(block, [], 3); feats normalize(feats); rng(42); [idx, C] kmeans(feats, K, Replicates, 3, Display, off); labelBlock reshape(idx, [bh, bw]); if isempty(allCenters) % 第一块直接写入 globalLabel(r:r2, c:c2) labelBlock; allCenters C; else % 后续块用重叠区域匹配中心点 if r 1 overlapRegion globalLabel(r:roverlap-1, c:c2); blockOverlap labelBlock(1:overlap, :); % 对每个新中心找重叠区里最匹配的旧标签 mapping zeros(K, 1); for k 1:K mask blockOverlap k; if any(mask(:)) oldLabels overlapRegion(mask); mapping(k) mode(oldLabels(:)); end end % 重编号 newLabelBlock zeros(size(labelBlock)); for k 1:K if mapping(k) 0 newLabelBlock(labelBlock k) mapping(k); else newLabelBlock(labelBlock k) max(globalLabel(:)) 1; end end labelBlock newLabelBlock; end globalLabel(r:r2, c:c2) labelBlock; allCenters [allCenters; C]; end end end % 显示结果 figure; imshow(label2rgb(globalLabel, jet, k, shuffle)); title(分块 K-means 分割结果);这段代码的核心是重叠区域的标签映射。mode取众数把新块里某一类在重叠区对应的旧标签找出来然后整块重编号。如果某一类在重叠区没有出现就分配一个新标签。这样拼起来的大图标签是连续的。注意normalize对每块单独做可能导致块间特征尺度不一致更稳的做法是用全局均值和标准差归一化但需要先扫一遍全图统计内存换精度。5.3 从分割到测量把标签图变成面积、周长和计数分割不是终点测量才是。拿到标签图后对每个标签做二值化用regionprops算面积、周长、质心、偏心度。面积用Area周长用Perimeter注意Perimeter对锯齿边界会偏大可以用bwperim后算像素数再乘像素物理尺寸。计数用bwlabel对二值图连通域编号max就是个数。% measure_regions.m % 对分割结果做区域测量 K 4; stats struct(label, {}, area, {}, perimeter, {}, centroid, {}); for k 1:K bw globalLabel k; bw bwareaopen(bw, 50); % 去掉小于 50 像素的噪点 bw imfill(bw, holes); props regionprops(bw, Area, Perimeter, Centroid); for j 1:length(props) stats(end1).label k; stats(j).area props(j).Area; stats(j).perimeter props(j).Perimeter; stats(j).centroid props(j).Centroid; end end % 输出统计表 T struct2table(stats); disp(T);bwareaopen去掉小噪点imfill填洞这两个后处理在测量前几乎必做。regionprops的Area是像素数要换算成物理面积需要知道像素分辨率。如果是显微镜图像通常有标定信息用im2double和空间参考对象处理。没有标定的话面积只能做相对比较。5.4 一个我常犯的错误忘了归一化导致坐标项吞掉颜色早期做分割时我把像素坐标直接拼到 RGB 后面坐标范围 1 到 1000RGB 范围 0 到 255结果距离计算完全被坐标主导分割出来是按位置切的方块颜色信息几乎没用。后来改成坐标归一化到 0 到 1颜色也归一化再乘权重才正常。这个坑很隐蔽因为代码不报错结果看起来也分了但分的是位置不是内容。现在我的习惯是任何多源特征拼接前先normalize或者手动缩放到 [0, 1]再乘权重。权重本身也要试wxy从 0 开始每次加 0.1看分割结果什么时候开始出现明显的空间块状那个点就是上限。希望帮到你。本文还有配套的精品资源点击获取
