Matlab 2014a 实现 Normalized Cuts 图像分割：原理、代码与避坑指南

简介：这份资源是Normalized Cuts（规范化切割）图像分割算法的Matlab实现代码包，面向计算机视觉方向的学生、研究人员及需要复现经典分割算法的开发者。算法由Shi和Malik提出，通过构建图模型、定义归一化切分并转化为谱聚类问题来求解最优分割，资源完整覆盖从相似度矩阵构建、归一化拉普拉斯矩阵计算到特征分解与K-means聚类的实现链路。压缩包共77个文件，约103KB，以50个m脚本为主体，辅以12个cpp源码及5个mexw64预编译文件用于加速矩阵运算，另有jpg示例图像、txt说明与db数据文件，便于直接运行demo并对照理解各模块职责。目前已有344人学习下载。读者可借此掌握Ncut算法的Matlab工程实现细节，理解谱分析、特征向量选择与聚类后处理流程，并可将思路迁移至视频分割、社区发现等场景，是入门图割分割与谱聚类实践的实用参考。

1. 从一次图像分割翻车说起：这个 Matlab 实现到底能干什么

去年帮一个做医学影像的朋友处理一批细胞切片图，他想要把粘连在一起的细胞核自动分开。我一开始用阈值加形态学操作，调了一下午参数，结果稍微换一批图就崩了。后来换成谱聚类思路，用 Normalized Cuts 做图划分，效果稳得多。但网上的 Python 实现要么依赖太重，要么跑一张 512×512 的图要等好几分钟。直到我翻出这份基于 Matlab 2014a 的 Normalized Cuts 实现，才真正把流程跑顺了。

这份资源解决的就是「把图像分割当成图划分问题」这件事。它把像素当成图的节点，像素之间的相似度当成边的权重，然后通过求解广义特征值问题，找到让子图之间「割」最小、同时子图内部尽量紧密的划分方式。适合谁用？做图像分割、医学影像处理、计算机视觉课程作业的开发者，尤其是实验室环境还在用 Matlab 2014a 的老项目。如果你手头有大量图片要批量分割，又不想折腾深度学习环境，这个实现值得拆开看看。

2. 核心原理与代码结构：从相似度矩阵到特征向量

2.1 Normalized Cuts 的数学骨架

Normalized Cuts 的核心目标函数是：

Ncut(A,B) = cut(A,B)/assoc(A,V) + cut(A,B)/assoc(B,V)

其中cut(A,B)是子图 A 和 B 之间所有边的权重和，assoc(A,V)是 A 中所有节点到全图节点的权重和。这个归一化项是关键——它避免了把孤立点单独切出来的退化解。

求解过程分三步：先构造相似度矩阵 W 和度矩阵 D，然后求解广义特征值问题(D-W)x = λDx，最后对第二小特征值对应的特征向量做离散化，得到二值分割结果。

我一般会提醒自己：这一步的玄学在于相似度函数的选取。用高斯核exp(-||Fi-Fj||²/σ²)时，σ 的取值直接决定分割粒度。σ 太大，所有像素都相似，割不出来；σ 太小，图变成一堆孤岛，计算量爆炸。

2.2 目录结构与文件职责

这份实现的文件组织很紧凑，常见做法是：

NormalizedCuts/ ├── main.m # 入口脚本，读图、调参、显示结果 ├── computeW.m # 构造相似度矩阵 W ├── computeD.m # 构造度矩阵 D ├── ncut.m # 核心求解：特征值分解 + 离散化 ├── discretize.m # 将连续特征向量转为二值分割 ├── image2graph.m # 图像像素到图节点的映射 └── utils/ ├── imgload.m # 图像读取与预处理 └── visualize.m # 分割结果可视化

main.m是唯一需要手动改参数的入口。computeW.m和computeD.m负责把图像转成图结构。ncut.m是黑匣子，里面调用eigs做稀疏特征值分解。discretize.m决定最终分割边界落在哪里。

2.3 相似度矩阵构造：computeW.m 逐段拆解

function W = computeW(img, sigma, radius) % img: 输入图像，灰度或 RGB % sigma: 高斯核带宽，控制相似度衰减速度 % radius: 邻域半径，只计算每个像素周围 radius 范围内的边 [H, W_img, ~] = size(img); N = H * W_img; img = double(img); % 将图像展平成 N x 3 的特征矩阵（RGB）或 N x 1（灰度） features = reshape(img, N, []); % 初始化稀疏矩阵，预估非零元素个数 % 每个像素最多连接 (2*radius+1)^2 个邻居 max_edges = N * (2*radius+1)^2; W = spalloc(N, N, max_edges); for i = 1:N % 计算当前像素的行列坐标 [r, c] = ind2sub([H, W_img], i); % 确定邻域范围，注意边界裁剪 r_min = max(1, r - radius); r_max = min(H, r + radius); c_min = max(1, c - radius); c_max = min(W_img, c + radius); for rr = r_min:r_max for cc = c_min:c_max j = sub2ind([H, W_img], rr, cc); if i == j continue; % 对角线不计算 end % 高斯核计算相似度 diff = features(i,:) - features(j,:); dist_sq = sum(diff.^2); W(i,j) = exp(-dist_sq / (2 * sigma^2)); end end end % 对称化，消除浮点误差 W = (W + W') / 2; end

这段代码的逻辑很直白：遍历每个像素，只和邻域内的像素建立边。sigma控制相似度衰减——两个像素颜色差越大，权重越小。radius控制图的稀疏程度，半径越大，边越多，计算越慢但分割越平滑。

参数说明：sigma常见取值 0.1 到 0.5 之间，具体看图像动态范围。如果是 0-255 的灰度图，先归一化到 0-1 再算。radius一般取 3 到 10，太大内存扛不住，太小分割结果碎。

注意：这里用双重循环在 Matlab 里跑会很慢。实际项目里我一般会向量化，或者用pdist2配合sparse一次性构造。但这份实现保留了循环写法，方便理解每一步在干什么。

2.4 特征值求解与离散化：ncut.m 的关键几行

function [labels, eigvecs, eigvals] = ncut(W, D, k) % W: 相似度矩阵 % D: 度矩阵，对角线上是每行权重和 % k: 期望的分割区域数 % 求解广义特征值问题 (D-W)x = lambda*D*x % 取前 k 个最小特征值对应的特征向量 opts.issym = 1; opts.isreal = 1; [eigvecs, eigvals] = eigs(D - W, D, k, 'sm', opts); % 忽略第一个特征向量（对应特征值 0，全 1 向量，无分割信息） eigvecs = eigvecs(:, 2:k); eigvals = diag(eigvals); eigvals = eigvals(2:k); % 对特征向量做行归一化，消除尺度影响 for i = 1:size(eigvecs, 1) eigvecs(i,:) = eigvecs(i,:) / norm(eigvecs(i,:)); end % 离散化：把连续特征向量转成离散标签 labels = discretize(eigvecs, k); end

eigs是 Matlab 的稀疏特征值求解器，'sm'表示取最小特征值。这里有个坑：D-W是半正定矩阵，最小特征值接近 0，数值稳定性差。常见做法是加一个小的正则项，或者用sigma参数做位移反演。

discretize.m通常用 k-means 对特征向量聚类，或者用简单的阈值分割。如果是二分割，直接看第二小特征向量的正负号就行。

3. 跑起来：Matlab 2014a 环境配置与完整执行流程

3.1 环境准备与依赖检查

Matlab 2014a 是个老版本，但这份实现只依赖基础工具箱和eigs函数，不需要额外安装。检查清单：

• Matlab 2014a 或更高版本（2016b 之后eigs行为有变化，但接口兼容）
• 图像处理工具箱（用于imread、imresize、rgb2gray）
• 统计工具箱（如果discretize.m里用了kmeans）

我一般会在脚本开头加一行ver确认工具箱版本。如果eigs报错说找不到，检查是不是装了 Symbolic Math Toolbox 冲突了——这种情况少见但遇到过。

3.2 入口脚本 main.m 的完整走读

% main.m clear; clc; close all; % 1. 读取图像并预处理 img = imread('test.jpg'); if size(img, 3) == 3 img = rgb2gray(img); % 转灰度，减少计算量 end img = imresize(img, [128, 128]); % 降采样，控制图规模 img = double(img) / 255; % 归一化到 0-1 % 2. 构造图结构 sigma = 0.15; % 高斯核带宽 radius = 5; % 邻域半径 W = computeW(img, sigma, radius); D = computeD(W); % 3. 求解 Normalized Cuts k = 2; % 二分割 [labels, eigvecs, eigvals] = ncut(W, D, k); % 4. 可视化 visualize(img, labels);

这段脚本的每一步都有讲究。imresize到 128×128 是血泪经验——原图 512×512 的话，W 矩阵是 262144×262144，稀疏存储也扛不住。降到 128×128 后，节点数 16384，稀疏矩阵非零元素约 16384×121≈200 万，内存占用可控。

sigma=0.15和radius=5是调参起点。如果分割结果太碎，增大sigma或radius；如果分割不出来，减小sigma。

3.3 运行结果解读与参数调优

跑完main.m后，visualize.m会显示原图和分割掩码。常见输出：

• 二分割：前景和背景用不同颜色标注
• 多分割：k>2时，每个区域一种颜色

调参时我一般按这个顺序：

1. 先固定radius=5，调sigma。从 0.05 开始，每次加 0.05，看分割边界是否合理。
2. sigma定下来后，调radius。增大半径会让分割更平滑，但计算时间线性增长。
3. 如果k>2，检查discretize.m里的聚类数是否和k一致。

提示：Matlab 2014a 的eigs在矩阵规模超过 10000 时可能报内存不足。解决办法是进一步降采样，或者用eigs的'lm'选项配合位移反演。

4. 避坑与排查：跑 Normalized Cuts 时最容易翻车的 5 个地方

4.1 特征值求解报错「Matrix is singular」

现象：eigs抛出警告或错误，说矩阵接近奇异。

原因：D-W的最小特征值接近 0，数值精度不够时求解器认为矩阵奇异。

解决：在ncut.m里给D-W加一个小的正则项，比如D - W + 1e-6 * speye(N)。或者改用eigs(D-W, D, k, 'sm', opts)时把opts.tol设小一点，比如1e-8。

4.2 分割结果全黑或全白

现象：labels全是 1 或全是 0，没有分割效果。

原因：sigma太大，所有像素相似度都接近 1，图退化成完全图，第二小特征向量没有区分度。

解决：减小sigma。如果是 0-1 归一化图像，sigma从 0.05 开始试。另外检查computeW.m里有没有把对角线元素也算了——对角线应该是 0。

4.3 内存溢出「Out of memory」

现象：跑 256×256 以上图像时 Matlab 卡死或报内存不足。

原因：W 矩阵是 N×N 稀疏矩阵，N=65536 时，即使稀疏存储，非零元素也可能上千万。

解决：降采样到 128×128 以下。或者用radius=3减小邻域。如果必须处理大图，考虑分块处理，每块单独分割后再合并。

4.4 分割边界锯齿严重

现象：分割结果边缘不平滑，像马赛克。

原因：radius太小，图连接稀疏，特征向量在边界处跳变。

解决：增大radius到 7 或 9。或者对labels做后处理，用形态学闭运算平滑边界。

4.5 多分割时区域数不对

现象：设了k=3，结果只分出 2 个区域，或者分出 5 个。

原因：discretize.m里的聚类算法对初始值敏感，或者特征向量本身区分度不够。

解决：在discretize.m里用kmeans时设置'Replicates', 5，多跑几次取最优。或者检查eigs返回的特征向量个数是否真的是k-1个。

5. 二次开发：从二分割到多区域、从单图到批量处理

5.1 扩展到多区域分割

默认实现是二分割，要支持k>2，需要改两个地方：

% 在 ncut.m 里，取前 k 个最小特征值对应的特征向量 [eigvecs, eigvals] = eigs(D - W, D, k, 'sm', opts); % 注意：第一个特征向量对应特征值 0，要跳过 eigvecs = eigvecs(:, 2:k);

然后在discretize.m里用 k-means 聚类：

function labels = discretize(eigvecs, k) % 对特征向量矩阵的每一行做聚类 % 每一行是一个节点的 k-1 维嵌入 if k == 2 % 二分割直接看正负号 labels = double(eigvecs(:,1) > 0) + 1; else % 多分割用 k-means opts = statset('MaxIter', 500, 'Display', 'off'); labels = kmeans(eigvecs, k, 'Replicates', 5, 'Options', opts); end end

这里的关键是：特征向量矩阵的每一行代表一个节点在低维空间中的嵌入。k-means 在这个嵌入空间里聚类，比在原空间里直接聚类效果好得多。

5.2 批量处理与参数自动化

如果有一批图要处理，我一般会写一个批处理脚本：

% batch_process.m image_files = dir('images/*.jpg'); results = cell(length(image_files), 1); % 参数网格 sigma_list = [0.1, 0.15, 0.2]; radius_list = [3, 5, 7]; for i = 1:length(image_files) img = imread(fullfile('images', image_files(i).name)); img = rgb2gray(img); img = imresize(img, [128, 128]); img = double(img) / 255; best_score = -inf; best_labels = []; for sigma = sigma_list for radius = radius_list W = computeW(img, sigma, radius); D = computeD(W); [labels, ~, eigvals] = ncut(W, D, 2); % 用特征值间隙作为分割质量评分 score = eigvals(2) / eigvals(1); if score > best_score best_score = score; best_labels = labels; end end end results{i} = best_labels; fprintf('Processed %s, best score: %.4f\n', image_files(i).name, best_score); end

这个脚本用第二小特征值和第一小特征值的比值作为分割质量指标。比值越大，说明分割越明确。自动选参数比手动调省事，但计算量翻倍。

5.3 性能优化：向量化与 GPU 加速

Matlab 2014a 对 GPU 支持有限，但向量化还是可以做的。computeW.m里的双重循环可以改成：

% 向量化版本：一次性计算所有邻域像素的相似度 [rows, cols] = meshgrid(1:H, 1:W_img); indices = sub2ind([H, W_img], rows(:), cols(:)); % 构造邻域索引矩阵 neighbor_offsets = []; for dr = -radius:radius for dc = -radius:radius if dr == 0 && dc == 0 continue; end neighbor_offsets = [neighbor_offsets; dr, dc]; end end % 批量计算 for k = 1:size(neighbor_offsets, 1) dr = neighbor_offsets(k, 1); dc = neighbor_offsets(k, 2); r2 = rows + dr; c2 = cols + dc; valid = r2 >= 1 & r2 <= H & c2 >= 1 & c2 <= W_img; idx1 = indices(valid); idx2 = sub2ind([H, W_img], r2(valid), c2(valid)); diff = features(idx1,:) - features(idx2,:); dist_sq = sum(diff.^2, 2); weights = exp(-dist_sq / (2 * sigma^2)); W = W + sparse(idx1, idx2, weights, N, N); end W = (W + W') / 2;

向量化后速度能快 5 到 10 倍，具体看 Matlab 版本和 CPU。如果装了 Parallel Computing Toolbox，把最外层的for k改成parfor还能再快一点。

5.4 一个我常犯的错误

有次赶项目，直接把sigma设成 1.0，想着「让相似度平滑一点总没错」。结果跑出来全图一个颜色，分割完全失效。后来才想明白：sigma是高斯核的带宽，图像归一化到 0-1 后，像素差最大也就 1，sigma=1.0意味着所有像素相似度都接近exp(-0.5)≈0.6，图变成均匀权重，Ncut 退化成随机切。

从那以后我每次调sigma都强制走一遍「从 0.05 开始，每次加 0.05，观察分割边界变化」的流程，再也不敢拍脑袋设参数了。希望这个习惯也能帮到你。

本文还有配套的精品资源，点击获取