简介这份资源是Normalized Cuts规范化切割算法的Matlab实现代码包面向计算机视觉、图像处理方向的学习者与研究者帮助理解并复现这一经典图论分割方法。包内共77个文件以50个m脚本、12个cpp源码、5个mexw64编译文件为主另有少量jpg示例图、txt说明与db数据文件压缩包约103KB涵盖相似度矩阵构建、归一化拉普拉斯矩阵计算、特征分解、特征向量选择与K-means聚类等完整流程并附有图像分割与聚类演示脚本。已有344人学习下载。读者可借助这些代码深入理解Normalized Cuts从图模型构建到谱聚类求解的实现细节对照示例图像验证分割效果并在此基础上扩展到视频分割、社区发现等场景是入门谱聚类与图像分割的实用参考。1. 从一次图像分割翻车说起这个 Matlab 实现到底能干什么去年帮一个做医学影像的朋友处理一批细胞切片图他丢过来一句“把前景背景分开就行”我心想阈值分割三行代码搞定。结果光照不均加上细胞粘连阈值法直接翻车边缘碎得像被狗啃过。后来翻出这份 Normalized Cuts 的 Matlab 实现才把粘连区域干净地切开。Normalized Cuts简称 Ncut是谱聚类在图像分割里的经典落地核心思路是把图像看成一张带权图像素是节点像素间的相似度是边权然后通过求解广义特征值问题把图切成若干块让“块内相似度高、块间相似度低”。这份资源就是它在 Matlab 环境下的完整实现适配 Matlab2014a 及以上版本适合做图像分割、区域提取、聚类实验的开发者尤其是需要快速验证算法效果又不想从零推导公式的人。它不依赖深度学习框架纯矩阵运算拿过来改改相似度函数就能跑自己的数据。2. 算法骨架与工程结构先看清 Ncut 是怎么被拆成代码的2.1 从图论到矩阵Ncut 的数学骨架Ncut 的目标函数长这样给定图 G(V,E)把节点分成 A、B 两部分最小化Ncut(A,B) cut(A,B)/assoc(A,V) cut(A,B)/assoc(B,V)其中 cut(A,B) 是 A 和 B 之间所有边的权重和assoc(A,V) 是 A 中节点与全图所有边的权重和。这个归一化项是关键它避免了最小割偏向切掉孤立小区域。直接优化是 NP 难的所以工程上把它松弛成连续问题求解 (D-W)x λDx 的广义特征值问题W 是相似度矩阵D 是度矩阵。取第二小特征值对应的特征向量用它做阈值分割或 k-means 聚类就得到最终分割结果。这份 Matlab 实现就是围绕这个流程展开的。它没有用任何工具箱里的谱聚类函数而是手写相似度矩阵构建、稀疏矩阵求解和特征向量后处理方便你逐行看懂每一步在干什么。2.2 目录结构与模块职责拿到源码包目录大致是这样组织的不同版本可能略有差异但核心文件跑不掉Ncut/ ├── main.m % 入口脚本读图、调参、显示结果 ├── computeW.m % 构建相似度矩阵 W ├── computeD.m % 由 W 计算度矩阵 D ├── ncut.m % 核心求解广义特征值返回分割向量 ├── discretize.m % 把连续特征向量离散化成标签 ├── imgsegment.m % 封装完整流程输入图像输出分割图 ├── utils/ % 辅助函数 │ ├── imread_color.m % 统一读图接口 │ └── showmask.m % 叠加显示分割掩膜 └── data/ % 示例图片main.m是给你直接跑的里面通常预设了一张示例图。computeW.m和computeD.m负责图构建ncut.m是计算最重的部分discretize.m决定最终标签怎么定。这种拆分方式很“教学向”每个文件对应一个数学步骤改起来定位快。2.3 相似度矩阵构建W 的三种常见写法相似度矩阵 W 直接决定分割质量。这份实现里computeW.m一般提供两种模式全连接和稀疏连接。全连接就是每个像素和所有像素算相似度公式是W(i,j) exp(-||F(i)-F(j)||^2 / sigma_F) * exp(-||X(i)-X(j)||^2 / sigma_X)F 是像素特征灰度或颜色X 是空间坐标。两个 sigma 分别控制特征相似度和空间距离的敏感度。全连接矩阵是 N×N 稠密矩阵N 是像素数一张 200×200 的图就是 40000×40000内存直接爆炸。所以实际代码里通常用稀疏化只保留每个像素周围 r 邻域内的连接或者用 kd-tree 找 k 近邻。% computeW.m 片段稀疏邻域相似度构建 function W computeW(img, sigmaF, sigmaX, radius) [h, w, ~] size(img); N h * w; % 把图像拉成 N x 3 的特征矩阵RGB feat reshape(img, N, 3); % 生成像素坐标 [xx, yy] meshgrid(1:w, 1:h); coord [xx(:), yy(:)]; % 用 kd-tree 找每个点的邻域避免全连接 % 常见做法是 rangesearch 或 knnsearch idx rangesearch(coord, coord, radius); % 逐点计算权重组装稀疏矩阵 I []; J []; S []; for i 1:N nb idx{i}; for j nb if i j, continue; end df sum((feat(i,:) - feat(j,:)).^2); dx sum((coord(i,:) - coord(j,:)).^2); w exp(-df / sigmaF) * exp(-dx / sigmaX); I(end1) i; J(end1) j; S(end1) w; end end W sparse(I, J, S, N, N); W (W W) / 2; % 保证对称 end这段代码里sigmaF和sigmaX是调参重点。sigmaF越大颜色差异被容忍得越多分割区域更倾向于按空间连通性走sigmaX越大空间距离影响越弱容易出现远距离同色区域被合并。radius控制邻域半径一般取 5 到 15 像素。注意最后做了对称化因为逐点计算可能因为浮点误差导致 W 不对称而后续特征值求解要求对称矩阵。2.4 广义特征值求解ncut.m 里的计算核心ncut.m是整个项目最“重”的文件。它要解的是(D - W) x λ D xMatlab 里用eigs求解前 k 个最小特征值。因为 D 是对角矩阵可以转成标准特征值问题令 y D^(1/2) x则D^(-1/2) (D - W) D^(-1/2) y λ y这样就能用eigs直接算对称矩阵的特征值。% ncut.m 片段求解第二小特征向量 function [eigvec, eigval] ncut(W, D, k) N size(W, 1); % 构造归一化拉普拉斯矩阵 Dinv_sqrt spdiags(1./sqrt(diag(D)), 0, N, N); L Dinv_sqrt * (D - W) * Dinv_sqrt; L (L L) / 2; % 强制对称避免数值误差 % 求前 k 个最小特征值对应的特征向量 % 注意最小特征值恒为 0对应平凡解要跳过 [V, E] eigs(L, k 1, sm); eigval diag(E); eigvec V(:, 2:k1); % 丢掉第一个平凡解 % 把特征向量还原回原空间 eigvec Dinv_sqrt * eigvec; end这里有几个坑点。第一eigs的sm选项在小矩阵上可能不稳定如果报错可以换成sa或者直接eig(full(L))但后者对大规模矩阵不现实。第二求出来的特征向量需要归一化常见做法是按行做 L2 归一化让每个像素的特征向量长度为 1这样后续聚类对尺度不敏感。第三k是你要分割的区域数二分割取 k1多分割取 k1。代码里eigs多求一个特征值就是为了丢掉那个恒为 0 的平凡解。2.5 离散化从连续特征向量到最终标签特征向量是连续的不能直接当标签用。discretize.m负责把它变成离散区域编号。二分割最简单取第二小特征向量找一个阈值常用 0 或者中位数把像素分成两类。多分割则要对前 k 个特征向量做 k-means 聚类。% discretize.m 片段二分割与多分割 function labels discretize(eigvec, k) if k 1 % 二分割用零阈值或中位数 v eigvec(:, 1); thresh median(v); % 比 0 更鲁棒 labels double(v thresh) 1; else % 多分割对特征向量矩阵做 k-means % 常见做法是先按行归一化 X eigvec ./ repmat(sqrt(sum(eigvec.^2, 2)), 1, k); labels kmeans(X, k, Replicates, 5, MaxIter, 200); end end二分割里用中位数而不是 0 做阈值是因为特征向量的分布不一定以 0 为中心中位数能保证两类像素数量大致均衡。多分割里kmeans的Replicates参数设成 5 是为了避免陷入局部最优代价是计算时间翻几倍。如果数据量大可以降到 2 或 3。3. 跑起来从零到一张分割图的完整操作3.1 环境准备与版本适配这份实现针对 Matlab2014a 及以上版本。2014a 是个分水岭之前的版本eigs对稀疏矩阵支持不够好kmeans的参数名也有差异。如果你用的是 2016b 之后的版本基本可以直接跑但要注意两点一是rangesearch在旧版本里叫rangesearch新版本没变但knnsearch的返回参数格式在 2019a 之后有调整二是kmeans的Replicates参数在 2014a 里是支持的但默认距离度量是sqeuclidean新版本默认也是这个一般不用改。没有额外工具箱依赖Image Processing Toolbox 只在读图和显示时用到核心算法纯基础矩阵运算。如果你连 IPT 都不想装把imread换成imread原生函数、显示部分用imagesc替代即可。3.2 最小可运行示例跑通 main.m打开main.m通常长这样% main.m最小可运行入口 clear; clc; close all; % 1. 读图 img imread(data/sample.jpg); img imresize(img, [128, 128]); % 降采样控制矩阵规模 figure; imshow(img); title(原图); % 2. 参数设置 sigmaF 0.1; % 特征相似度尺度 sigmaX 5; % 空间距离尺度 radius 8; % 邻域半径 k 2; % 分割区域数 % 3. 构建 W 和 D W computeW(double(img)/255, sigmaF, sigmaX, radius); D computeD(W); % 4. 求解 Ncut [eigvec, eigval] ncut(W, D, k-1); % 5. 离散化 labels discretize(eigvec, k-1); % 6. 显示结果 mask reshape(labels, size(img,1), size(img,2)); figure; imshow(mask, []); title(分割结果);跑之前先把data/sample.jpg换成你自己的图或者把imread路径改成绝对路径。imresize到 128×128 是为了控制计算量N16384稀疏矩阵还能扛住。如果你直接上 512×512N262144eigs会跑很久内存也吃紧。3.3 参数怎么调sigmaF、sigmaX、radius 的联动这三个参数没有理论最优值只能按图调。我一般按这个顺序试参数作用调大效果调小效果常用范围sigmaF控制颜色相似度敏感度颜色差异被忽略区域更连通颜色差异被放大分割更碎0.05~0.5sigmaX控制空间距离敏感度远距离同色区域可能合并只保留局部连接区域更紧凑2~20radius邻域半径连接更多计算变慢连接稀疏可能断连5~15调参时先固定 radius8sigmaX5然后只动 sigmaF。如果分割结果把该分开的区域粘在一起说明 sigmaF 太大往下调如果同一区域被切碎说明 sigmaF 太小往上调。sigmaX 一般不用大动除非你的图像分辨率特别高或特别低。3.4 结果可视化与掩膜叠加showmask.m通常做一件事把分割标签转成彩色掩膜叠加到原图上。核心代码就几行% showmask.m叠加显示 function showmask(img, labels) % 生成随机颜色表 nLabels max(labels); colors rand(nLabels, 3); % 把标签映射成 RGB maskRGB reshape(colors(labels, :), size(img)); % 叠加原图 0.6 掩膜 0.4 overlay 0.6 * double(img)/255 0.4 * maskRGB; figure; imshow(overlay); title(分割叠加); end这样能直观看到分割边界落在哪里。如果边界和实际物体边缘对不齐回去调 sigmaF 和 radius。4. 避坑与排查跑 Ncut 最容易翻车的五个地方4.1 现象eigs 报错“矩阵太大”或直接卡死原因全连接 W 是 N×N 稠密矩阵N 超过 5000 就基本不可行。很多人第一次跑忘了稀疏化或者 radius 设得太大导致稀疏矩阵退化成稠密。解决确认computeW.m里用的是sparse而不是zeros。把 radius 降到 5~8或者先用imresize把图像降到 128×128 以下。如果还卡改用eigs(L, k1, sm, Tolerance, 1e-6, MaxIter, 300)限制迭代次数。4.2 现象分割结果全黑或全白没有分界原因特征向量阈值选错了。二分割时如果直接用 0 做阈值而特征向量整体偏移就会全分到一类。解决把阈值改成median(v)或mean(v)。更稳的做法是对特征向量做归一化后再取中位数。如果还是不行检查ncut.m里是不是忘了丢掉第一个平凡解导致拿到的特征向量恒为常数。4.3 现象kmeans 报错“X 包含 NaN 或 Inf”原因特征向量里有 NaN通常是因为 D 的对角线有 0求逆时除了零。D 的对角线是 W 的行和如果某个像素和所有像素的相似度都下溢到 0行和就是 0。解决在computeD.m里加一个极小值保护D spdiags(max(sum(W,2), 1e-10), 0, N, N)。或者在ncut.m里对 Dinv_sqrt 做同样处理。这个坑在图像有大片纯黑或纯白区域时特别常见。4.4 现象分割边界呈块状锯齿不平滑原因Ncut 本身是像素级分割没有做任何平滑或后处理。加上降采样后边界自然粗糙。解决如果只是可视化需求可以在discretize之后加一个中值滤波或形态学闭运算。如果是精度需求别降采样直接在全分辨率上跑但要做好等几分钟的准备。另一个办法是增大 radius让邻域连接更密边界会更连续代价是计算量上升。4.5 现象换一张图就完全不能用参数全要重调原因sigmaF 和 sigmaX 是绝对尺度对图像亮度和分辨率敏感。一张图调好的参数换一张图就失效。解决把 sigmaF 改成自适应先算全图特征差值的标准差令 sigmaF 0.1 * std(df)。sigmaX 同理用坐标差值的标准差做基准。这样参数变成相对值跨图迁移性好很多。这个改动我一般直接写在computeW.m开头省得每次手动调。5. 二次开发把 Ncut 改成你自己的分割工具5.1 替换相似度函数从 RGB 到多特征融合默认的 W 只用 RGB 和坐标。如果你做医学影像可以把灰度、纹理、梯度方向都拼进特征向量。改法很简单在computeW.m里把feat从 N×3 扩成 N×d然后调整 sigmaF 的维度。比如加一个 LBP 纹理特征% 扩展特征RGB LBP gray rgb2gray(img); lbp extractLBPFeatures(gray); % 需要 Computer Vision Toolbox % 或者手写简化版 LBP feat [reshape(img, N, 3), reshape(lbp, N, 1)]; % sigmaF 要相应调大因为维度变了注意特征维度增加后欧氏距离的尺度会变sigmaF 要重新标定。我一般先用pdist算一下特征距离的分布取中位数的 0.1 倍作为 sigmaF 初值。5.2 加速策略从 eigs 到 Lobpcg 或子采样eigs对大规模稀疏矩阵还是慢。两个加速方向一是用 LOBPCG 算法Matlab 没有内置但可以自己实现或者找第三方二是子采样先在小图上算特征向量再用插值或最近邻把标签映射回全分辨率。第二种更实用代码改动也小% 子采样加速先降采样求标签再上采样 img_small imresize(img, 0.25); labels_small imgsegment(img_small, sigmaF, sigmaX, radius, k); labels_full imresize(labels_small, size(img(:,:,1)), nearest);这样计算量降到 1/16边界会粗一点但大多数场景够用。如果边界精度要求高可以在全分辨率上只对边界附近的像素重新计算。5.3 多区域分割的标签合并技巧k2 时kmeans 出来的标签顺序是随机的同一张图跑两次可能标签编号对不上。如果你要做批量处理需要固定标签顺序。常见做法是按区域面积排序或者按区域质心的空间位置排序。我一般用面积降序% 固定标签顺序按面积从大到小重排 areas histcounts(labels, 1:k1); [~, order] sort(areas, descend); map zeros(1, k); map(order) 1:k; labels map(labels);这样每次跑出来的标签 1 都是最大区域方便后续统计和对比。5.4 一个我踩过的坑稀疏矩阵对称化的代价早期我为了省事在computeW.m里只算上三角然后W W W。结果发现分割结果左右不对称。排查半天才意识到sparse组装时如果 I、J 有重复索引Matlab 会把对应值相加而不是覆盖。上三角加转置后对角线元素被加了两次。后来改成先组装完整矩阵再对称化或者组装时跳过对角线。从那以后我每次写完稀疏矩阵都强制走一遍issymmetric(W)检查不通过就不往下跑。希望这个习惯也能帮你省掉几个小时的排查时间。本文还有配套的精品资源点击获取