简介这份MATLAB源码包面向信号处理、无线通信与图像处理方向的学习者和工程师聚焦压缩感知CS理论中从低采样率测量值恢复稀疏信号的核心问题。包内共18个文件以17个.m脚本和1个.fig图形文件为主压缩包约20KB脚本分别对应OMP、CoSaMP、IST等经典重构算法的实现与对比以及信道估计演示、MMSE与LS误差计算等辅助模块fig文件可用于查看算法性能对比曲线。资源围绕稀疏表示与重构算法展开涵盖正交匹配追踪的贪婪迭代、压缩采样匹配追踪的多原子更新、迭代软阈值的阈值收缩等思路并延伸至无线通信信道估计场景便于读者运行脚本、复现实验并比较不同算法在重构精度与抗噪性上的差异。目前已有1477人学习下载适合希望快速上手压缩感知重构实验、对照源码理解算法流程并开展性能验证的读者。1. 压缩感知重构信号为什么采样率降了一个数量级信号还能还原如果你做过数据采集大概率遇到过这种局面传感器采样率被硬件卡死奈奎斯特率要求你每秒采几千甚至几万个点但存储、带宽、功耗三座大山压着根本存不下也传不走。压缩感知Compressed Sensing, CS给出的答案有点反直觉——只要信号本身在某组基下是稀疏的你就可以用远低于奈奎斯特率的采样点数把它压着采再靠重构算法把原始信号还原出来。这不是玄学是有数学保证的。这套东西真正落地的场景很具体稀疏信道估计、单像素成像、核磁共振加速、宽带频谱感知、无线传感器网络里的低功耗节点。适合谁看手里有 MATLAB想跑通「稀疏表示 → 观测矩阵 → 重构算法 → 误差评估」这条完整链路的人。下面我按自己实际做项目的顺序把压缩感知重构信号这件事从选型讲到代码再讲到踩过的坑。2. 压缩感知重构算法的三条技术路线从凸松弛到贪婪迭代2.1 重构问题的数学骨架与稀疏前提先把问题写清楚。设原始信号为 $x \in \mathbb{R}^N$观测矩阵 $\Phi \in \mathbb{R}^{M \times N}$且 $M \ll N$观测向量 $y \Phi x$。这是一个欠定方程组解有无穷多个。压缩感知能唯一重构的前提是 $x$ 在某个稀疏基 $\Psi$ 下是稀疏的即 $x \Psi s$而 $s$ 中非零元素个数 $K \ll N$。此时观测变成 $y \Phi \Psi s A s$$A$ 叫感知矩阵。重构的本质是求解$$\min |s|_0 \quad \text{s.t.} \quad y A s$$$\ell_0$ 范数最小化是 NP 难的所以实际算法都在做两件事之一要么把 $\ell_0$ 松弛成 $\ell_1$ 变成凸问题要么用贪婪策略一步步逼近稀疏解。这就是后面所有算法的分水岭。观测矩阵 $\Phi$ 不能随便选。它必须满足 RIP有限等距性质或者至少和稀疏基 $\Psi$ 不相关。工程上最省事的做法是取高斯随机矩阵或伯努利随机矩阵因为它们以极高概率满足 RIP而且和绝大多数稀疏基都不相关。我一般先用高斯矩阵跑通再考虑能不能换成部分傅里叶矩阵、托普利兹矩阵这类有物理可实现性的结构。2.2 凸松弛类算法BP、ISTA、FISTA 的取舍凸松弛的核心是把 $\ell_0$ 换成 $\ell_1$$$\min |s|_1 \quad \text{s.t.} \quad y A s$$这就是基追踪Basis Pursuit, BP。它重构精度高、需要的观测数少但计算量大本质是一个线性规划问题。MATLAB 里可以用linprog硬解但维度一上去就慢得让人想砸键盘。更实用的是迭代软阈值类算法。ISTAIterative Shrinkage-Thresholding Algorithm的迭代格式是$$s_{k1} \text{soft}\left(s_k \frac{1}{L} A^T (y - A s_k),\ \frac{\lambda}{L}\right)$$其中 $L$ 是 $A^T A$ 的最大特征值 Lipschitz 常数$\text{soft}$ 是软阈值算子。FISTA 在 ISTA 基础上加了一个 Nesterov 动量项收敛速度从 $O(1/k)$ 提升到 $O(1/k^2)$代价只是多存一个中间变量。我实测下来同样精度下 FISTA 的迭代次数通常只有 ISTA 的一半左右。下面是我常用的 FISTA 最小实现直接可跑function [s_hat, err_hist] fista_cs(y, A, lambda, max_iter, tol) % FISTA 求解 min 0.5*||y-A*s||_2^2 lambda*||s||_1 % y: 观测向量 Mx1, A: 感知矩阵 MxN, lambda: 正则化参数 % max_iter: 最大迭代次数, tol: 收敛阈值 [~, N] size(A); L norm(A)^2; % Lipschitz 常数用谱范数平方估计 s zeros(N, 1); z s; % 动量变量 t 1; err_hist zeros(max_iter, 1); for k 1:max_iter grad A * (A * z - y); % 梯度 s_new soft_thresh(z - grad / L, lambda / L); t_new (1 sqrt(1 4 * t^2)) / 2; z s_new ((t - 1) / t_new) * (s_new - s); % Nesterov 动量 s s_new; t t_new; err_hist(k) norm(y - A * s); if err_hist(k) tol err_hist err_hist(1:k); break; end end s_hat s; end function y soft_thresh(x, T) y sign(x) .* max(abs(x) - T, 0); end逻辑说明L norm(A)^2是 Lipschitz 常数的保守估计用norm算谱范数在中小规模下够用大规模可以改用幂迭代。lambda控制稀疏度越大解越稀疏但可能丢真实分量。tol我一般设成观测噪声水平的量级比如1e-4 * norm(y)。动量项那两行是 FISTA 的精髓删掉就退化成 ISTA可以自己对比收敛曲线。2.3 贪婪类算法OMP、CoSaMP、SP 的适用边界贪婪算法的思路完全不同每次挑和当前残差最相关的原子把它加进支撑集再最小二乘更新系数。正交匹配追踪OMP是最基础的一个每轮只加一个原子迭代 $K$ 次$K$ 是稀疏度。OMP 的优点是快、实现简单、不需要调正则化参数缺点是稀疏度 $K$ 必须已知而且一旦某轮选错原子就无法回头。压缩采样匹配追踪CoSaMP每轮选 $2K$ 个候选原子再剪枝回 $K$ 个容错性更好。子空间追踪SP每轮选 $K$ 个再回溯精度介于两者之间。选型上我的经验是稀疏度已知且信噪比高用 OMP 最快稀疏度未知或噪声大用 CoSaMP 或 FISTA追求重构精度且不在乎时间用 BP。下面给一个 OMP 的实现function [s_hat, support] omp_cs(y, A, K) % 正交匹配追踪K 为已知稀疏度 [~, N] size(A); r y; % 初始残差 support []; % 支撑集 s_hat zeros(N, 1); for k 1:K proj abs(A * r); % 计算所有原子与残差的相关性 [~, idx] max(proj); support [support, idx]; % 加入支撑集 A_sel A(:, support); s_sel A_sel \ y; % 最小二乘更新 r y - A_sel * s_sel; % 更新残差 end s_hat(support) s_sel; end逻辑说明A * r是相关性度量选最大者加入支撑集。A_sel \ y用 MATLAB 反斜杠做最小二乘比显式求伪逆数值更稳。注意support可能重复选到同一个原子工程上要加去重判断否则A_sel会奇异。这个坑我在第一次写的时候踩得很实。3. 在 MATLAB 里跑通压缩感知重构从造数据到出误差曲线3.1 构造稀疏信号与观测矩阵的完整脚本光看公式没用得有一份能直接跑的实验脚本。下面这段从造稀疏信号开始到生成观测、调用重构、画误差曲线一条龙clear; clc; rng(42); % 固定随机种子保证可复现 N 512; % 信号长度 M 200; % 观测数压缩比约 0.39 K 20; % 稀疏度 lambda 0.01; % FISTA 正则化参数 % 1. 构造 K-稀疏信号 s_true zeros(N, 1); pos randperm(N, K); s_true(pos) randn(K, 1); % 2. 高斯随机观测矩阵 Phi randn(M, N) / sqrt(M); % 3. 观测无噪声 y Phi * s_true; % 4. FISTA 重构 [s_fista, err_hist] fista_cs(y, Phi, lambda, 2000, 1e-6); % 5. OMP 重构 [s_omp, ~] omp_cs(y, Phi, K); % 6. 误差评估 err_fista norm(s_fista - s_true) / norm(s_true); err_omp norm(s_omp - s_true) / norm(s_true); fprintf(FISTA 相对误差: %.4f\n, err_fista); fprintf(OMP 相对误差: %.4f\n, err_omp); % 7. 画图对比 figure; subplot(3,1,1); stem(s_true, .); title(原始稀疏信号); subplot(3,1,2); stem(s_fista, .); title([FISTA 重构, 误差, num2str(err_fista, %.4f)]); subplot(3,1,3); stem(s_omp, .); title([OMP 重构, 误差, num2str(err_omp, %.4f)]);逻辑说明rng(42)固定种子是必须的否则每次跑出来的误差都在跳没法对比算法。Phi randn(M,N)/sqrt(M)里的归一化让每列能量接近 1避免数值量级失控。lambda和K是两个最关键的参数下面单独说。3.2 观测数 M、稀疏度 K、正则化 lambda 三个参数的调法这三个参数决定了重构能不能成功我按重要性排。观测数 M理论上 M 要满足 $M \geq C \cdot K \log(N/K)$$C$ 通常取 2 到 4。以 N512、K20 为例$K \log(N/K) \approx 20 \times 3.2 64$所以 M 至少 128 到 256。我一般从 $M 4K$ 起步逐步往下压看误差什么时候开始崩。压到 $M 2K$ 基本就没救了。稀疏度 KOMP 和 CoSaMP 必须知道 K。如果实际信号稀疏度未知可以先跑 FISTA 看解的非零个数反推一个 K 的估计。或者用 CoSaMP它对 K 的过估计容忍度比 OMP 高。正则化 lambdaFISTA 里 lambda 控制稀疏度。lambda 太小解不稀疏、噪声被放大lambda 太大真实分量被阈值掉。经验值是lambda 0.01 * max(abs(A * y))到0.1 * max(abs(A * y))之间。有噪声时按噪声标准差调lambda ≈ sigma * sqrt(2 * log(N))是常用起点。参数作用推荐起点调过头会怎样M观测数4K太小误差爆炸K稀疏度真实非零数过估计引入伪影lambda稀疏惩罚0.05*max(abs(Ay))太大丢分量太小不稀疏3.3 用相对误差和支撑集恢复率评估重构质量只看相对误差不够还要看支撑集恢复得对不对。相对误差小但支撑集错位说明解在数值上接近但结构不对换到实际应用里可能完全不可用。% 支撑集恢复率 support_true find(abs(s_true) 1e-6); support_est find(abs(s_fista) 1e-3); recovery_rate numel(intersect(support_true, support_est)) / numel(support_true); fprintf(支撑集恢复率: %.2f%%\n, recovery_rate * 100);阈值1e-3是经验值信号量级不同要调。我一般同时看三个指标相对误差、支撑集恢复率、运行时间。三个都达标才算这个参数组合可用。4. 压缩感知重构的避坑与排查那些让误差曲线突然翘起来的细节4.1 观测矩阵没归一化导致 lambda 完全失效现象FISTA 跑出来全是零或者全是噪声误差比不重构还大。原因Phi randn(M, N)没除以sqrt(M)导致A * y的量级随 M 变化lambda 的绝对值失去意义。解决观测矩阵每列能量归一化或者用Phi randn(M,N); Phi Phi ./ repmat(sqrt(sum(Phi.^2,1)), M, 1)。归一化后 lambda 才有可比性。4.2 OMP 支撑集重复选原子导致矩阵奇异现象OMP 跑到某一轮报Matrix is singular或者结果突然发散。原因max(proj)选到的原子已经在支撑集里了A_sel出现重复列最小二乘无解。解决选原子前把已选位置的相关性置零proj(support) 0;再取 max。这个改动一行就够但不知道的话能卡半天。4.3 稀疏基和观测矩阵相关导致 RIP 不满足现象所有算法误差都下不去加大 M 也没用。原因稀疏基 $\Psi$ 和观测矩阵 $\Phi$ 相关性太高感知矩阵 $A \Phi \Psi$ 的 RIP 常数超标。解决换随机观测矩阵或者对 $\Psi$ 做随机化处理。如果 $\Psi$ 是傅里叶基观测矩阵就别用部分傅里叶矩阵改用高斯随机矩阵。我一般先用高斯矩阵验证算法本身没问题再换结构化矩阵。4.4 噪声水平估计错误导致 lambda 选偏现象有噪声时重构结果要么过拟合噪声要么把信号一起阈值掉。原因lambda 按无噪声场景设的没考虑噪声标准差。解决先估计噪声水平sigma median(abs(A*y))/0.6745这个 0.6745 是正态分布中位数绝对偏差的归一化常数再设lambda sigma * sqrt(2*log(N))。有噪声时这个起点比拍脑袋靠谱得多。4.5 用 linprog 解 BP 时维度上去直接卡死现象N 到几千linprog跑几个小时不出结果。原因BP 的线性规划规模随 N 立方增长MATLAB 的linprog默认内点法在大规模下内存和迭代都吃不消。解决N 超过 1000 就别用 BP 了换 FISTA 或 ADMM。真要解 BP用专门的 $\ell_1$ 求解器或者 ADMM 分裂把问题拆成小规模子问题。5. 进阶技巧用 ADMM 把重构问题拆开以及怎么验证你真的做对了5.1 ADMM 求解压缩感知的变量分裂写法FISTA 在中等规模下够用但遇到 $A$ 条件数很差或者要加额外约束比如非负、总变差时ADMM 更灵活。ADMM 的核心是把问题写成$$\min f(s) g(z) \quad \text{s.t.} \quad s z$$对压缩感知取 $f(s) 0.5|y - As|_2^2$$g(z) \lambda|z|_1$。迭代格式function s admm_cs(y, A, lambda, rho, max_iter) % ADMM 求解 min 0.5*||y-A*s||^2 lambda*||s||_1 [~, N] size(A); s zeros(N, 1); z s; u s; AtA A * A; Aty A * y; for k 1:max_iter % s 更新解线性方程组 s (AtA rho * eye(N)) \ (Aty rho * (z - u)); % z 更新软阈值 z_old z; z soft_thresh(s u, lambda / rho); % 对偶变量更新 u u s - z; if norm(s - z) 1e-6 norm(z - z_old) 1e-6 break; end end end逻辑说明rho是惩罚参数影响收敛速度但不影响最终解。rho太大收敛慢太小对偶变量震荡我一般从 1 开始试。(AtA rho*eye(N))可以预分解循环里只做回代这是 ADMM 比 FISTA 快的地方。soft_thresh函数和前面 FISTA 里共用。5.2 用已知真值做交叉验证别只看一条误差曲线重构算法最容易骗自己的地方是调参数调到某条误差曲线好看就以为成了。我的习惯是造多组不同 K、不同 M、不同噪声水平的信号每组跑 50 次蒙特卡洛看成功率的统计。K_list [10, 20, 30, 40]; M_list 100:20:300; success_rate zeros(numel(K_list), numel(M_list)); for i 1:numel(K_list) for j 1:numel(M_list) cnt 0; for trial 1:50 s_true zeros(N,1); s_true(randperm(N, K_list(i))) randn(K_list(i),1); Phi randn(M_list(j), N) / sqrt(M_list(j)); y Phi * s_true; s_hat fista_cs(y, Phi, 0.01, 2000, 1e-6); if norm(s_hat - s_true)/norm(s_true) 1e-3 cnt cnt 1; end end success_rate(i,j) cnt / 50; end end imagesc(M_list, K_list, success_rate); colorbar; xlabel(观测数 M); ylabel(稀疏度 K); title(重构成功率相图);逻辑说明成功率阈值1e-3是相对误差比这个松说明重构质量不够。imagesc画出来的相图能直观看到成功区和失败区的边界比单条曲线信息量大得多。这张图也是判断「我的参数组合到底行不行」最硬的证据。5.3 一个我反复用的验证习惯我做完任何压缩感知重构实验最后一定做一件事把重构出来的信号重新乘一遍观测矩阵看norm(Phi*s_hat - y)和norm(Phi*s_true - y)差多少。如果重构信号在观测域对不上那它在原始域再像也是假的。这个检查花不了几秒钟但帮我拦下过好几次「误差曲线好看、实际不可用」的情况。压缩感知重构信号这件事算法本身不复杂难的是参数和边界条件。我的习惯是先用高斯矩阵和无噪声数据把算法跑通再逐步加噪声、换结构化矩阵、压观测数每一步只改一个变量。这样出问题的时候能立刻定位到是哪一步引入的。希望帮到你。本文还有配套的精品资源点击获取 压缩感知重构信号：MATLAB实现FISTA与OMP算法及参数调优
简介这份MATLAB源码包面向信号处理、无线通信与图像处理方向的学习者和工程师聚焦压缩感知CS理论中从低采样率测量值恢复稀疏信号的核心问题。包内共18个文件以17个.m脚本和1个.fig图形文件为主压缩包约20KB脚本分别对应OMP、CoSaMP、IST等经典重构算法的实现与对比以及信道估计演示、MMSE与LS误差计算等辅助模块fig文件可用于查看算法性能对比曲线。资源围绕稀疏表示与重构算法展开涵盖正交匹配追踪的贪婪迭代、压缩采样匹配追踪的多原子更新、迭代软阈值的阈值收缩等思路并延伸至无线通信信道估计场景便于读者运行脚本、复现实验并比较不同算法在重构精度与抗噪性上的差异。目前已有1477人学习下载适合希望快速上手压缩感知重构实验、对照源码理解算法流程并开展性能验证的读者。1. 压缩感知重构信号为什么采样率降了一个数量级信号还能还原如果你做过数据采集大概率遇到过这种局面传感器采样率被硬件卡死奈奎斯特率要求你每秒采几千甚至几万个点但存储、带宽、功耗三座大山压着根本存不下也传不走。压缩感知Compressed Sensing, CS给出的答案有点反直觉——只要信号本身在某组基下是稀疏的你就可以用远低于奈奎斯特率的采样点数把它压着采再靠重构算法把原始信号还原出来。这不是玄学是有数学保证的。这套东西真正落地的场景很具体稀疏信道估计、单像素成像、核磁共振加速、宽带频谱感知、无线传感器网络里的低功耗节点。适合谁看手里有 MATLAB想跑通「稀疏表示 → 观测矩阵 → 重构算法 → 误差评估」这条完整链路的人。下面我按自己实际做项目的顺序把压缩感知重构信号这件事从选型讲到代码再讲到踩过的坑。2. 压缩感知重构算法的三条技术路线从凸松弛到贪婪迭代2.1 重构问题的数学骨架与稀疏前提先把问题写清楚。设原始信号为 $x \in \mathbb{R}^N$观测矩阵 $\Phi \in \mathbb{R}^{M \times N}$且 $M \ll N$观测向量 $y \Phi x$。这是一个欠定方程组解有无穷多个。压缩感知能唯一重构的前提是 $x$ 在某个稀疏基 $\Psi$ 下是稀疏的即 $x \Psi s$而 $s$ 中非零元素个数 $K \ll N$。此时观测变成 $y \Phi \Psi s A s$$A$ 叫感知矩阵。重构的本质是求解$$\min |s|_0 \quad \text{s.t.} \quad y A s$$$\ell_0$ 范数最小化是 NP 难的所以实际算法都在做两件事之一要么把 $\ell_0$ 松弛成 $\ell_1$ 变成凸问题要么用贪婪策略一步步逼近稀疏解。这就是后面所有算法的分水岭。观测矩阵 $\Phi$ 不能随便选。它必须满足 RIP有限等距性质或者至少和稀疏基 $\Psi$ 不相关。工程上最省事的做法是取高斯随机矩阵或伯努利随机矩阵因为它们以极高概率满足 RIP而且和绝大多数稀疏基都不相关。我一般先用高斯矩阵跑通再考虑能不能换成部分傅里叶矩阵、托普利兹矩阵这类有物理可实现性的结构。2.2 凸松弛类算法BP、ISTA、FISTA 的取舍凸松弛的核心是把 $\ell_0$ 换成 $\ell_1$$$\min |s|_1 \quad \text{s.t.} \quad y A s$$这就是基追踪Basis Pursuit, BP。它重构精度高、需要的观测数少但计算量大本质是一个线性规划问题。MATLAB 里可以用linprog硬解但维度一上去就慢得让人想砸键盘。更实用的是迭代软阈值类算法。ISTAIterative Shrinkage-Thresholding Algorithm的迭代格式是$$s_{k1} \text{soft}\left(s_k \frac{1}{L} A^T (y - A s_k),\ \frac{\lambda}{L}\right)$$其中 $L$ 是 $A^T A$ 的最大特征值 Lipschitz 常数$\text{soft}$ 是软阈值算子。FISTA 在 ISTA 基础上加了一个 Nesterov 动量项收敛速度从 $O(1/k)$ 提升到 $O(1/k^2)$代价只是多存一个中间变量。我实测下来同样精度下 FISTA 的迭代次数通常只有 ISTA 的一半左右。下面是我常用的 FISTA 最小实现直接可跑function [s_hat, err_hist] fista_cs(y, A, lambda, max_iter, tol) % FISTA 求解 min 0.5*||y-A*s||_2^2 lambda*||s||_1 % y: 观测向量 Mx1, A: 感知矩阵 MxN, lambda: 正则化参数 % max_iter: 最大迭代次数, tol: 收敛阈值 [~, N] size(A); L norm(A)^2; % Lipschitz 常数用谱范数平方估计 s zeros(N, 1); z s; % 动量变量 t 1; err_hist zeros(max_iter, 1); for k 1:max_iter grad A * (A * z - y); % 梯度 s_new soft_thresh(z - grad / L, lambda / L); t_new (1 sqrt(1 4 * t^2)) / 2; z s_new ((t - 1) / t_new) * (s_new - s); % Nesterov 动量 s s_new; t t_new; err_hist(k) norm(y - A * s); if err_hist(k) tol err_hist err_hist(1:k); break; end end s_hat s; end function y soft_thresh(x, T) y sign(x) .* max(abs(x) - T, 0); end逻辑说明L norm(A)^2是 Lipschitz 常数的保守估计用norm算谱范数在中小规模下够用大规模可以改用幂迭代。lambda控制稀疏度越大解越稀疏但可能丢真实分量。tol我一般设成观测噪声水平的量级比如1e-4 * norm(y)。动量项那两行是 FISTA 的精髓删掉就退化成 ISTA可以自己对比收敛曲线。2.3 贪婪类算法OMP、CoSaMP、SP 的适用边界贪婪算法的思路完全不同每次挑和当前残差最相关的原子把它加进支撑集再最小二乘更新系数。正交匹配追踪OMP是最基础的一个每轮只加一个原子迭代 $K$ 次$K$ 是稀疏度。OMP 的优点是快、实现简单、不需要调正则化参数缺点是稀疏度 $K$ 必须已知而且一旦某轮选错原子就无法回头。压缩采样匹配追踪CoSaMP每轮选 $2K$ 个候选原子再剪枝回 $K$ 个容错性更好。子空间追踪SP每轮选 $K$ 个再回溯精度介于两者之间。选型上我的经验是稀疏度已知且信噪比高用 OMP 最快稀疏度未知或噪声大用 CoSaMP 或 FISTA追求重构精度且不在乎时间用 BP。下面给一个 OMP 的实现function [s_hat, support] omp_cs(y, A, K) % 正交匹配追踪K 为已知稀疏度 [~, N] size(A); r y; % 初始残差 support []; % 支撑集 s_hat zeros(N, 1); for k 1:K proj abs(A * r); % 计算所有原子与残差的相关性 [~, idx] max(proj); support [support, idx]; % 加入支撑集 A_sel A(:, support); s_sel A_sel \ y; % 最小二乘更新 r y - A_sel * s_sel; % 更新残差 end s_hat(support) s_sel; end逻辑说明A * r是相关性度量选最大者加入支撑集。A_sel \ y用 MATLAB 反斜杠做最小二乘比显式求伪逆数值更稳。注意support可能重复选到同一个原子工程上要加去重判断否则A_sel会奇异。这个坑我在第一次写的时候踩得很实。3. 在 MATLAB 里跑通压缩感知重构从造数据到出误差曲线3.1 构造稀疏信号与观测矩阵的完整脚本光看公式没用得有一份能直接跑的实验脚本。下面这段从造稀疏信号开始到生成观测、调用重构、画误差曲线一条龙clear; clc; rng(42); % 固定随机种子保证可复现 N 512; % 信号长度 M 200; % 观测数压缩比约 0.39 K 20; % 稀疏度 lambda 0.01; % FISTA 正则化参数 % 1. 构造 K-稀疏信号 s_true zeros(N, 1); pos randperm(N, K); s_true(pos) randn(K, 1); % 2. 高斯随机观测矩阵 Phi randn(M, N) / sqrt(M); % 3. 观测无噪声 y Phi * s_true; % 4. FISTA 重构 [s_fista, err_hist] fista_cs(y, Phi, lambda, 2000, 1e-6); % 5. OMP 重构 [s_omp, ~] omp_cs(y, Phi, K); % 6. 误差评估 err_fista norm(s_fista - s_true) / norm(s_true); err_omp norm(s_omp - s_true) / norm(s_true); fprintf(FISTA 相对误差: %.4f\n, err_fista); fprintf(OMP 相对误差: %.4f\n, err_omp); % 7. 画图对比 figure; subplot(3,1,1); stem(s_true, .); title(原始稀疏信号); subplot(3,1,2); stem(s_fista, .); title([FISTA 重构, 误差, num2str(err_fista, %.4f)]); subplot(3,1,3); stem(s_omp, .); title([OMP 重构, 误差, num2str(err_omp, %.4f)]);逻辑说明rng(42)固定种子是必须的否则每次跑出来的误差都在跳没法对比算法。Phi randn(M,N)/sqrt(M)里的归一化让每列能量接近 1避免数值量级失控。lambda和K是两个最关键的参数下面单独说。3.2 观测数 M、稀疏度 K、正则化 lambda 三个参数的调法这三个参数决定了重构能不能成功我按重要性排。观测数 M理论上 M 要满足 $M \geq C \cdot K \log(N/K)$$C$ 通常取 2 到 4。以 N512、K20 为例$K \log(N/K) \approx 20 \times 3.2 64$所以 M 至少 128 到 256。我一般从 $M 4K$ 起步逐步往下压看误差什么时候开始崩。压到 $M 2K$ 基本就没救了。稀疏度 KOMP 和 CoSaMP 必须知道 K。如果实际信号稀疏度未知可以先跑 FISTA 看解的非零个数反推一个 K 的估计。或者用 CoSaMP它对 K 的过估计容忍度比 OMP 高。正则化 lambdaFISTA 里 lambda 控制稀疏度。lambda 太小解不稀疏、噪声被放大lambda 太大真实分量被阈值掉。经验值是lambda 0.01 * max(abs(A * y))到0.1 * max(abs(A * y))之间。有噪声时按噪声标准差调lambda ≈ sigma * sqrt(2 * log(N))是常用起点。参数作用推荐起点调过头会怎样M观测数4K太小误差爆炸K稀疏度真实非零数过估计引入伪影lambda稀疏惩罚0.05*max(abs(Ay))太大丢分量太小不稀疏3.3 用相对误差和支撑集恢复率评估重构质量只看相对误差不够还要看支撑集恢复得对不对。相对误差小但支撑集错位说明解在数值上接近但结构不对换到实际应用里可能完全不可用。% 支撑集恢复率 support_true find(abs(s_true) 1e-6); support_est find(abs(s_fista) 1e-3); recovery_rate numel(intersect(support_true, support_est)) / numel(support_true); fprintf(支撑集恢复率: %.2f%%\n, recovery_rate * 100);阈值1e-3是经验值信号量级不同要调。我一般同时看三个指标相对误差、支撑集恢复率、运行时间。三个都达标才算这个参数组合可用。4. 压缩感知重构的避坑与排查那些让误差曲线突然翘起来的细节4.1 观测矩阵没归一化导致 lambda 完全失效现象FISTA 跑出来全是零或者全是噪声误差比不重构还大。原因Phi randn(M, N)没除以sqrt(M)导致A * y的量级随 M 变化lambda 的绝对值失去意义。解决观测矩阵每列能量归一化或者用Phi randn(M,N); Phi Phi ./ repmat(sqrt(sum(Phi.^2,1)), M, 1)。归一化后 lambda 才有可比性。4.2 OMP 支撑集重复选原子导致矩阵奇异现象OMP 跑到某一轮报Matrix is singular或者结果突然发散。原因max(proj)选到的原子已经在支撑集里了A_sel出现重复列最小二乘无解。解决选原子前把已选位置的相关性置零proj(support) 0;再取 max。这个改动一行就够但不知道的话能卡半天。4.3 稀疏基和观测矩阵相关导致 RIP 不满足现象所有算法误差都下不去加大 M 也没用。原因稀疏基 $\Psi$ 和观测矩阵 $\Phi$ 相关性太高感知矩阵 $A \Phi \Psi$ 的 RIP 常数超标。解决换随机观测矩阵或者对 $\Psi$ 做随机化处理。如果 $\Psi$ 是傅里叶基观测矩阵就别用部分傅里叶矩阵改用高斯随机矩阵。我一般先用高斯矩阵验证算法本身没问题再换结构化矩阵。4.4 噪声水平估计错误导致 lambda 选偏现象有噪声时重构结果要么过拟合噪声要么把信号一起阈值掉。原因lambda 按无噪声场景设的没考虑噪声标准差。解决先估计噪声水平sigma median(abs(A*y))/0.6745这个 0.6745 是正态分布中位数绝对偏差的归一化常数再设lambda sigma * sqrt(2*log(N))。有噪声时这个起点比拍脑袋靠谱得多。4.5 用 linprog 解 BP 时维度上去直接卡死现象N 到几千linprog跑几个小时不出结果。原因BP 的线性规划规模随 N 立方增长MATLAB 的linprog默认内点法在大规模下内存和迭代都吃不消。解决N 超过 1000 就别用 BP 了换 FISTA 或 ADMM。真要解 BP用专门的 $\ell_1$ 求解器或者 ADMM 分裂把问题拆成小规模子问题。5. 进阶技巧用 ADMM 把重构问题拆开以及怎么验证你真的做对了5.1 ADMM 求解压缩感知的变量分裂写法FISTA 在中等规模下够用但遇到 $A$ 条件数很差或者要加额外约束比如非负、总变差时ADMM 更灵活。ADMM 的核心是把问题写成$$\min f(s) g(z) \quad \text{s.t.} \quad s z$$对压缩感知取 $f(s) 0.5|y - As|_2^2$$g(z) \lambda|z|_1$。迭代格式function s admm_cs(y, A, lambda, rho, max_iter) % ADMM 求解 min 0.5*||y-A*s||^2 lambda*||s||_1 [~, N] size(A); s zeros(N, 1); z s; u s; AtA A * A; Aty A * y; for k 1:max_iter % s 更新解线性方程组 s (AtA rho * eye(N)) \ (Aty rho * (z - u)); % z 更新软阈值 z_old z; z soft_thresh(s u, lambda / rho); % 对偶变量更新 u u s - z; if norm(s - z) 1e-6 norm(z - z_old) 1e-6 break; end end end逻辑说明rho是惩罚参数影响收敛速度但不影响最终解。rho太大收敛慢太小对偶变量震荡我一般从 1 开始试。(AtA rho*eye(N))可以预分解循环里只做回代这是 ADMM 比 FISTA 快的地方。soft_thresh函数和前面 FISTA 里共用。5.2 用已知真值做交叉验证别只看一条误差曲线重构算法最容易骗自己的地方是调参数调到某条误差曲线好看就以为成了。我的习惯是造多组不同 K、不同 M、不同噪声水平的信号每组跑 50 次蒙特卡洛看成功率的统计。K_list [10, 20, 30, 40]; M_list 100:20:300; success_rate zeros(numel(K_list), numel(M_list)); for i 1:numel(K_list) for j 1:numel(M_list) cnt 0; for trial 1:50 s_true zeros(N,1); s_true(randperm(N, K_list(i))) randn(K_list(i),1); Phi randn(M_list(j), N) / sqrt(M_list(j)); y Phi * s_true; s_hat fista_cs(y, Phi, 0.01, 2000, 1e-6); if norm(s_hat - s_true)/norm(s_true) 1e-3 cnt cnt 1; end end success_rate(i,j) cnt / 50; end end imagesc(M_list, K_list, success_rate); colorbar; xlabel(观测数 M); ylabel(稀疏度 K); title(重构成功率相图);逻辑说明成功率阈值1e-3是相对误差比这个松说明重构质量不够。imagesc画出来的相图能直观看到成功区和失败区的边界比单条曲线信息量大得多。这张图也是判断「我的参数组合到底行不行」最硬的证据。5.3 一个我反复用的验证习惯我做完任何压缩感知重构实验最后一定做一件事把重构出来的信号重新乘一遍观测矩阵看norm(Phi*s_hat - y)和norm(Phi*s_true - y)差多少。如果重构信号在观测域对不上那它在原始域再像也是假的。这个检查花不了几秒钟但帮我拦下过好几次「误差曲线好看、实际不可用」的情况。压缩感知重构信号这件事算法本身不复杂难的是参数和边界条件。我的习惯是先用高斯矩阵和无噪声数据把算法跑通再逐步加噪声、换结构化矩阵、压观测数每一步只改一个变量。这样出问题的时候能立刻定位到是哪一步引入的。希望帮到你。本文还有配套的精品资源点击获取