简介这是一份面向压缩感知研究的稀疏贝叶斯算法代码包包含SBL、TSBL与TMSBL三种主流实现适合需要从少量测量值中重构稀疏信号的科研人员和工程师可应用于图像重建、医学成像、无线通信信道估计等任务。作者已亲自测试代码可顺利运行减少了调试成本。资源共15个文件主体为11个MATLAB脚本实现了SBL、TSBL、TMSBL及MSBL、MFOCUSS等算法并附带演示示例另有2份PDF和1份Word快速上手指南讲解算法使用要点txt文件说明整体流程。压缩包仅479KB获取方便、部署简单。目前已有929人学习该资源。通过阅读源码和文档读者不仅能复现时间变化、不同信噪比等实验场景还能深入对比SBL与TSBL、TMSBL在时间相关性和多尺度建模上的差异掌握贝叶斯稀疏重构的核心思路便于后续改造和拓展应用是一份高价值的压缩感知工具资源。1. 压缩感知稀疏贝叶斯SBL 家族为什么值得亲测一遍压缩感知稀疏贝叶斯学习SBL这个词做过信号重构的人多半不陌生。和 L1 范数类算法相比SBL 的迭代往往显得慢但在多快拍场景下它的支撑集恢复精度和噪声鲁棒性通常会反超。这份资源把 SBL、TSBL、TMSBL 三种稀疏贝叶斯算法的 MATLAB 实现打包在一起我实测过 demo.m、demo_time_varying.m 等脚本能直接跑通出图。适合正在做 DOA 估计、MIMO 信道估计、脑电或心电信号重构的工程师也适合刚接触多测量向量MMV压缩感知的学生拿去做基线对比。下文先说清楚三种算法的差别和选型边界再给出完整复现步骤、参数说明和踩坑记录。2. 三种算法的原理与选型从稀疏先验到时间相关性建模2.1 SBL 的核心思想把稀疏约束放进超参数SBL 和 L1 正则化的本质差别在于L1 把稀疏约束写成目标函数里的惩罚项而 SBL 把稀疏建模到参数的概率分布里。以包里的 MSBL.m 为例观测模型是 Y A X NX 是 N×T 的多快拍信号T 是时间快拍数。核心假设是行稀疏所有时间快拍共享同一个支撑集也就是 X 的非零行都落在同样的行索引上。这个假设在 DOA 估计里很自然——同一个声源在连续多个快拍里都位于相同方位只是幅度在变。实现上算法给每一行 X(i,:) 赋予一个方差超参数 gamma_i然后做 Type-II 最大似然估计。迭代过程中大部分 gamma_i 会被压到极小值对应行被自动剪掉剩下的大 gamma 行就是恢复出的支撑集。相比 L1 需要手调正则化系数SBL 顺手把噪声方差 sigma 也当成可学习参数这是它鲁棒性好的一个直接原因。我对这个机制的理解是L1 是在罚稀疏SBL 是在学稀疏后者对噪声水平的适应能力强得多。核心迭代可以写成下面的示意代码理解它对后续调参加大有帮助% SBL 类算法的核心迭代示意对应 MSBL.m 的内部逻辑 % Sigma inv(diag(gamma) A*A / sigma^2) % 每轮迭代: % Sigma_old inv(diag(gamma) A*A / sigma^2); % mu Sigma_old * A * Y; % 后验均值 % gamma_new mean(mu.^2, 2) diag(Sigma_old) / T; % sigma_new norm(Y - A*mu, fro)^2 / (M*T) ... % sum(diag(Sigma_old) .* sum(A.^2, 1)) / (M*T);逻辑说明mu 是后验均值Sigma_old 是后验协方差。gamma_new 的更新公式里第一项 mean(mu.^2, 2) 衡量该行在数据拟合里的能量第二项 diag(Sigma_old)/T 是超参数自身的惯性惩罚正是这一项让无关行的 gamma 迅速衰减。sigma_new 的更新拆成残差项和协方差迹项避免噪声估计被过拟合带偏。参数说明gamma 是 N×1 列向量对应 X 每一行的方差sigma^2 是标量噪声方差T 是快拍数。gamma 初始化通常全设为 1sigma 初始化可以先用观测残差估一个粗值后面算法会自己迭代修正。需要注意这段代码是理解用不是让你替换包内实现——MSBL.m 里还有数值稳定性处理直接改容易翻车。2.2 TSBL把时间连续性写进先验TSBLTemporal Sparse Bayesian Learning解决的是支撑集随时间不变、但幅度随时间变化的信号。典型场景同一个 DOA 源在连续多个快拍里一直存在只是幅度有起伏或者生理信号里某个频点在整段记录中都活跃。这种时间结构如果不用相当于浪费了宝贵的先验信息。TSBL 和 MSBL 的差别在于MSBL 假设不同时间快拍之间独立TSBL 则显式估计时间相关矩阵 B把 X 每一行在时间轴上的相关性建模进去。这样在快拍数不充裕时TSBL 能靠时间结构把支撑集估得更稳。包里的 demo_time_varying.m 就是这种场景的复现脚本信号支撑集固定、幅度随时间平滑变化。% TSBL 调用以包内 TSBL.m 实际签名为准 % 常见签名: [X, gamma, sigma] TSBL(A, Y, option); option.prune_gamma 1e-4; % gamma 剪枝阈值 option.maxit 500; % 最大迭代次数 option.learn_b 1; % 是否学习时间相关矩阵 B [X, gamma, sigma] TSBL(A, Y, option);逻辑说明learn_b 置 1 时算法会在每轮迭代里用当前残差估计时间相关矩阵再把它代回先验更新置 0 时退化成时间独立的版本。对平滑变化的信号learn_b1 能显著提高支撑集恢复率但对时间上本来就不相关的信号learn_b1 反而是负资产。参数说明prune_gamma 是剪枝阈值gamma 低于该值的行会在迭代中被直接剔除既省时间也防止数值震荡。maxit 是最大迭代次数实际多数场景在 100300 轮内收敛设 500 是为了防止算法在边界条件下没跑完就退出真跑满 500 轮还没收敛问题通常不在迭代次数上。2.3 TMSBL多尺度分块解决了什么问题TMSBL 是在 TSBL 基础上的进一步扩展引入了分块结构按先验知识把 N 行分成若干组每一组内部共享同一个方差超参数。这就是Multiple-Scale的含义——不同块可以用不同的稀疏力度去约束而不是所有行一刀切。这在信号天然有分组结构时特别有用比如一组传感器对应一组行、一组频点对应一组行。当信号的行天然存在分组结构时TMSBL 能把块稀疏先验用上样本效率更高。代价是必须提供合理的 group 划分group 给错了效果反而比 MSBL 差。包里的 TMSBL.m 和 TMSBL_code 目录就是干这个的配套的 Master the Usage of TMSBL in 3 Minutes.pdf 快速上手文档写得很清楚建议跑代码之前先花三分钟过一遍。% TMSBL 需要额外传入 group分块索引 % 常见签名: [X, gamma, sigma] TMSBL(A, Y, group, option); N 100; % 行数 block_size 10; group cell(1, N / block_size); for g 1:(N / block_size) group{g} (g-1)*block_size 1 : g*block_size; end option.prune_gamma 1e-4; [X, gamma, sigma] TMSBL(A, Y, group, option);逻辑说明group 是元胞数组每个元胞是一组行索引算法会为每个块学习一个独立的超参数同时对块内的时间相关性做联合建模。分块越符合信号的真实结构收敛越快、支撑集恢复越准。参数说明block_size 不是算法参数是信号本身的先验。如果你并不知道真实块大小宁可把块设小一点比如 510 行一组也不要设成整块一块后者会让 TMSBL 失去区分不同尺度的能力退化成近似 TSBL。2.4 选型边界什么时候该换算法选型判断我总结成三条剩下的交给实验说话只有一个观测向量T1用 SBL 类最朴素的版本就行TSBL/TMSBL 的时间建模没有用武之地。多快拍且支撑集固定信号在时间上平滑变化直接用 TMSBL时间相关性一般不明显的TSBL 足够。完全不确定信号结构先用 MSBL 跑一遍当基线再拿 TSBL/TMSBL 对比。如果信号本身没有时间相关性时间建模会帮倒忙结果会比 MSBL 更差这是正常现象不是代码坏了。我自己的习惯是永远先跑 MSBL 基线。没有基线就上 TMSBL出了问题你分不清是参数没调好、group 给错了还是算法本身不适合这个信号。有了基线三行对比就能定位问题。3. 跑通代码包目录结构、demo 与核心函数3.1 压缩包里的文件怎么组织解压后能看到 TMSBL_code 目录、MSBL.m、TSBL.m、TMSBL.m、MFOCUSS.m 以及一堆 demo 脚本。ReadMe.txt 是总的说明两个 PDF/DOCX 是 TSBL 和 TMSBL 的快速上手文档。我建议按这个顺序看ReadMe.txt 先读确认 MATLAB 版本要求以及函数依赖关系。跑 demo.m看 MSBL 的基本重构流程。跑 demo_time_varying.m 和 demo_identicalVector.m对比 TSBL/TMSBL 对不同时间结构的处理。需要复现论文图的再看 demo_fig3.m、demo_fig6_SNR10.m、demo_fig8.m。perfSupp.m 是支撑集恢复性能的计算函数后面验证指标会用到。MFOCUSS.m 是 FOCUSS 的 MMV 版本通常被当作对比基线跑对比实验时拿它和三兄弟放一起表格会很好看。3.2 demo.m 拆解从生成信号到调用 MSBLdemo.m 是最短的完整流程我把它的关键步骤拆出来% demo.m 的核心流程示例 clear; clc; close all; rng(42); % 固定随机种子保证可复现 N 100; % 原子数 / 信号维度 M 50; % 观测维度 T 10; % 时间快拍数 K 6; % 非零行数稀疏度 A randn(M, N); % 感知矩阵 for col 1:N A(:, col) A(:, col) / norm(A(:, col)); % 列归一化 end X_true zeros(N, T); supp randperm(N, K); X_true(supp, :) randn(K, T); Y A * X_true 0.01 * randn(M, T); % 观测加噪 [X, gamma, sigma] MSBL(A, Y); % 调用 MSBL support_est find(gamma 1e-3); hit length(intersect(support_est, supp)) / K; fprintf(support recovery: %.2f\n, hit);逻辑说明先固定随机种子保证每次跑出来的支撑集、噪声都一样这是复现和对比的前提。感知矩阵逐列归一化这是 SBL 类算法的隐藏前置条件后面避坑章节会专门讲。调用 MSBL 后算法返回后验均值 X、超参数 gamma 和估计的噪声方差 sigma。参数说明N、M、T、K 四个数决定问题规模。压缩比 M/N 越低重构难度越大T 太小比如小于 3TSBL/TMSBL 的时间协方差估计就不稳。K 是真实稀疏度对算法是未知的只能从 gamma 的分布间接观察所以 demo 里算 support recovery 前先用 gamma 阈值提取估计支撑集。3.3 TSBL/TMSBL 的调用方式与 group 参数TSBL 和 TMSBL 的调用和 MSBL 基本一致区别在选项和时间建模% TSBL 与 TMSBL 的调用对照 % TSBL: [X, gamma, sigma] TSBL(A, Y, option); % TMSBL: [X, gamma, sigma] TMSBL(A, Y, group, option); common_opt.prune_gamma 1e-4; common_opt.maxit 500; % 用同一份观测同时跑三个算法 [X_msbl, g_msbl, s_msbl] MSBL(A, Y); [X_tsbl, g_tsbl, s_tsbl] TSBL(A, Y, common_opt); [X_tmsbl, g_tmsbl, s_tmsbl] TMSBL(A, Y, group, common_opt); % 对比三者的支撑集 support_msbl find(g_msbl 1e-3); support_tsbl find(g_tsbl 1e-3); support_tmsbl find(g_tmsbl 1e-3);逻辑说明把三个算法放在同一份观测上跑是判断算法是否适合当前信号的最快方式。支撑集用 gamma 阈值提取阈值取 1e-3 还是 1e-4 要看 gamma 的量级量级偏小时取 1e-4 更安全。注意三个算法返回的 gamma 量级可能不同统一用一个阈值时要先各自看一眼分布再定。参数说明common_opt 在两个算法间复用保证对比公平。注意 TMSBL 多一个 group 参数如果这里传了错误的 groupTMSBL 的结果会明显劣于 TSBL不能据此判断TMSBL 不行只能说明 group 给错了。这里也是新手最容易误判的地方。4. 参数设置与结果验证怎么判断恢复得好不好4.1 关键参数表把三个算法的参数汇总成一张表方便对照参数所属算法建议取值作用gamma 初始值SBL/TSBL/TMSBL全 1每行方差超参数的起点prune_gammaSBL/TSBL/TMSBL1e-3 ~ 1e-4剪枝阈值低于此值的行被剔除maxitSBL/TSBL/TMSBL300 ~ 500最大迭代次数sigma 初始值SBL/TSBL/TMSBL残差粗估计噪声方差起点迭代中更新learn_bTSBL/TMSBL1 或 0是否学习时间相关矩阵groupTMSBL按先验分块块稀疏结构定义需要多说一句的是 sigma 初始值。很多人直接取 0.01 或 0.1这在 demo 里能跑换数据就崩。我一般先用最小二乘残差估一个初值量级对了算法自己会迭代修正量级错了就是灾难。prune_gamma 也别拍脑袋可以先跑一轮不剪枝看 gamma 的直方图分布选择两个量级之间的位置作为阈值。4.2 验证指标支持集恢复率与 NMSE跑通不是终点支撑集恢复率和归一化均方误差NMSE才是判断算法效果的两个硬指标% 支撑集恢复率 support_true supp; support_est find(gamma 1e-3); hit length(intersect(support_true, support_est)) / length(support_true); % 归一化均方误差 nmse norm(X_est - X_true, fro)^2 / norm(X_true, fro)^2; fprintf(hit %.2f, NMSE %.4f\n, hit, nmse);逻辑说明支撑集恢复率算的是真实非零行有多少被找到这个指标比 NMSE 对稀疏重构更敏感——NMSE 好看但支撑集找错位置对后续定位问题毫无意义。两个指标搭配看一个管位置准不准一个管幅度像不像缺一不可。参数说明gamma 阈值取 1e-3 时略激进可能把微弱行也剪掉取 1e-4 相对保守。对比不同算法时阈值必须一致否则对比结论是无效的。demo_fig6_SNR10.m 和 demo_fig8.m 分别对应不同 SNR 下的实验跑完用 perfSupp.m 汇总就能看到支撑集恢复率随信噪比的曲线。4.3 收敛判断看 gamma 的变化轨迹SBL 类算法没有显式的残差下降曲线可看我一般直接观察 gamma 的迭代轨迹。最实用的做法是记录每一轮 gamma 的最大值和非零个数画出来看是否进入平稳区。如果 gamma 最大值连续 20 轮变化小于 1%基本可以认为收敛了如果 gamma 一直在震荡优先怀疑 sigma 初始化过小或感知矩阵条件数太差。另外两个容易被忽略的点一是 prune_gamma 设得过大比如 1e-2会让算法过早剪枝把微弱但真实的非零行误杀恢复率会掉得很明显二是 demo_identicalVector.m 这个脚本信号各时间快拍完全相同时算法收敛极快别拿它的收敛轮数去预估自己数据的收敛速度二者没有可比性。提示对比算法必须用同一个随机种子、同一份观测 Y。换了观测谈对比等于把变量搞混结论不可信。5. 避坑记录SBL/TSBL/TMSBL 实战中我踩过的五个坑5.1 换了自己的数据后输出全是 NaN现象demo 跑得挺好把代码里的随机感知矩阵换成自己的测量矩阵MSBL 直接输出 NaN 或者全零矩阵。原因感知矩阵没做列归一化。SBL 类算法的更新公式里隐式假设 A 的列能量一致列范数差异过大时gamma 的迭代会数值崩溃。demo 里的 randn 矩阵列范数天然接近所以问题不容易暴露换真实数据才会翻车。解决在调用任何算法之前先对 A 做逐列归一化并同步检查 Y 的尺度。这个前置处理我放在封装函数的第一行比让每个调用者自己记得要靠谱得多。5.2 TSBL/TMSBL 报维度不匹配现象报错 Subscripted assignment dimension mismatch 或 Index exceeds array bounds。原因group 分块没有覆盖全部 N 行。比如 N100group 只包了 90 行或者块之间有重叠、有空洞。解决用一行代码校验 group 覆盖assert(isequal(sort([group{:}]), 1:N), group must cover 1:N exactly);这个断言我每次都会写进脚本里比看报错猜原因快得多。TMSBL 对 group 的完整性要求极高缺一行都会导致索引错乱。5.3 gamma 迟迟不收敛支撑集来回跳现象算法迭代到 maxit 也没停最后输出的支撑集和真实支撑集偏差很大而且换一次随机种子结果就变。原因sigma 初始化太小。噪声方差起点过低时算法会把噪声当成信号去拟合gamma 更新被残差主导收敛不到正确位置。解决先用最小二乘残差估一个粗噪声方差作为初始值sigma_init norm(Y - A * (A \ Y), fro)^2 / (M * T);这样一个粗估计通常就在正确量级附近算法迭代几轮就能收敛到稳定值。从那以后我每次跑新数据第一件事就是算这个值不再手动猜噪声水平。5.4 T1 时 TMSBL 反而比 MSBL 差现象只有一个观测向量的时候TMSBL 的时间协方差估计完全不 work结果比 MSBL 差一截。原因T1 时时间相关矩阵 B 根本没有足够的样本去估计TMSBL 等于多了一个无效约束自然拖后腿。这属于模型假设和应用场景不匹配。解决T 小于 5 时直接回退到 MSBL不要硬上 TSBL/TMSBL。这个规律在文档里没有明说但代码跑出来的趋势非常明显。我遇到单快拍任务时就直接跳过 TSBL/TMSBL 的对比环节了。5.5 时间不相关的信号用 TSBL 精度下降现象信号在时间上其实是独立的或者每快拍支撑集都在变用 TSBL 去重构支撑集恢复率反而低于 MSBL。原因TSBL 的先验假设是支撑集固定、幅度相关。信号不满足这个假设时时间相关矩阵 B 学到的是虚假相关性等于给模型强加了错误结构属于先验和事实拧着来。解决先用 demo_time_varying.m 和 demo_identicalVector.m 两个脚本对照判断自己的信号属于哪种形态再决定选型。支撑集随时间变化的信号本质上不是 TSBL/TMSBL 的目标场景这时应该回到 MSBL 或者换用别的动态稀疏算法别在不合适的算法上浪费时间调参。6. 进阶把 TMSBL 封装成自己的重构入口跑通 demo 之后最值得做的一件事是把这三个算法封装成一个统一的入口函数这样以后换数据只需要改参数不需要重读文档function [X, gamma, info] my_sbl_recovery(A, Y, alg, group, opts) % alg: msbl / tsbl / tmsbl % group: 仅 tmsbl 需要其余传 [] if nargin 5, opts struct(); end if ~isfield(opts, prune_gamma), opts.prune_gamma 1e-4; end if ~isfield(opts, maxit), opts.maxit 500; end switch lower(alg) case msbl [X, gamma, sigma] MSBL(A, Y); case tsbl [X, gamma, sigma] TSBL(A, Y, opts); case tmsbl assert(~isempty(group), TMSBL requires group); assert(isequal(sort([group{:}]), 1:size(A,2)), group mismatch); [X, gamma, sigma] TMSBL(A, Y, group, opts); otherwise error(unknown algorithm: %s, alg); end info.gamma gamma; info.sigma sigma; info.support find(gamma opts.prune_gamma); end逻辑说明封装的核心是统一输出接口把 gamma、sigma、支撑集一次性返回。info 结构体里带上支撑集提取结果后面无论是画图还是接入下游定位流程都不用再重复写 find(gamma thr) 这类代码。A 的列归一化也可以放进这个函数里一步到位。参数说明opts 用 isfield 判断缺省值自动补上这样新来的同事只传 A、Y、alg 也能跑不至于被参数项吓退。三个算法的返回值结构一致switch 只是分派调用真正的差异被挡在函数内部。封装完之后我建议你把 demo_fig3.m、demo_fig6_SNR10.m、demo_fig8.m 三张图分别在不同 SNR 下重跑一遍用 perfSupp.m 统计支撑集恢复率形成一张算法×SNR的对比表。这张表就是你以后跟别人讨论压缩感知稀疏贝叶斯算法时的硬证据比口头说效果不错有用得多。从那以后我每次拿到新的测量数据都会强制先跑一遍 MSBL 基线和两个时间序列 demo 的对照确认信号形态属于哪一类再决定上 TSBL 还是 TMSBL。这个习惯帮我少走了很多弯路希望帮到你。本文还有配套的精品资源点击获取