-
摘要: 针对压缩感知（Compressive sensing，CS）中未知稀疏度信号的重建问题，本文提出一种变步长稀疏度自适应子空间追踪算法.首先，采用一种匹配测试的方法确定固定步长，然后以该固定步长与变步长方式相结合，通过不同支撑集原子个数下的重建残差变化确定信号稀疏度，算法采用子空间追踪方法确定相应支撑集原子，并完成原始信号准确重建.实验结果表明，与同类算法相比，该算法可以更准确重建原始信号，且信号稀疏度值较高时，运算量低于同类算法.Abstract: A novel variable step size sparsity adaptive subspace pursuit algorithm is proposed to rebuild the sparse signals with unknown sparsity in compressive sensing. Firstly, the initial fixed step size is obtained by matching test, which is combined with the variable step size method. Then, the sparsity is accurately estimated according to the change of signal rebuilding residual error under vary support set. Subspace pursuit algorithm is used to determine the support set and exactly rebuild the sparse signal. Simulation results show that the proposed algorithm is competitive in recovering accuracy and running speed, compared to other similar algorithms, when sparsity is large.
-
压缩感知(Compressive sensing, CS)是Donoho、Candès及Tao等提出的一种新的信息获取理论[1?4].对于具有稀疏性或者在特定域上具有稀疏性的原始信号, 根据压缩感知理论, 我们能够采用远低于Nyquist采样的速率来对该信号同时进行采样和压缩, 得到低维观测信号, 并可通过观测信号准确完成高维原始信号的重建, 从而大大降低了原始信号获取、存储及传输的代价, 缓解了高速采样对硬件系统造成的压力.
CS理论的一个重要方面是信号重建算法, 常用的重建方法主要有贪婪追踪算法、凸优化算法和组合算法三大类.其中贪婪追踪算法由于结构简单、计算复杂度较低而备受关注.传统的贪婪算法有匹配追踪(Matching pursuit, MP)[5]以及正交匹配追踪(Orthogonal matching pursuit, OMP)[6]算法等, 这些算法在每一次迭代时对支撑集都只选择增加一个原子, 所以其运算复杂度较大; 而后提出的正则化正交匹配追踪(Regularized orthogonal matching pursuit, ROMP)算法[7]及分段正交匹配追踪(Stagewise orthogonal matching pursuit, StOMP)算法[8]对此进行了改进, 在每次迭代时可一次选择增加多个支撑集原子, 所以提高了重建速度; 后来提出的压缩采样匹配追踪(Compressive sampling matching pursuit, CoSaMP)算法[9]和子空间追踪(Subspace pursuit, SP)算法[10]在每次迭代时均一次选择多个原子, 同时采用回溯检验方法, 先求得一个原子个数较多的候选支撑集合, 再每次迭代中从该集合中增加与残差相关度高的原子、淘汰权值系数小的原子, 直到算法收敛, 该方法可提高重建精度.
上述算法实现信号重建的必要条件是必须获知原始信号稀疏度.但是在现实环境中, 这个条件难以满足, 因为大多数实际信号在其稀疏域仅是近似稀疏的, 例如图像信号, 不同的图像其稀疏度各不相同, 所以这些算法往往难以处理实际问题.为解决该问题, 文献[11]提出了稀疏度自适应匹配追踪(Sparsity adaptive matching pursuit, SAMP)算法, 该算法采用StOMP算法中的分段的思想来逐段扩充真实支撑集, 同时也利用了SP/CoSaMP算法中的回溯的思想, 算法不需要知道稀疏度的先验信息.还有其他一些算法[12?15]也可以实现稀疏度自适应信号重建.
这些算法均可以实现在未知稀疏度时的信号重建, 但是算法在迭代确定支撑集原子个数时, 一般采用固定步长方式逐步增加原子个数.当步长较小时, 迭代次数就多, 算法运算量大; 步长较大时, 可降低运算量, 但对稀疏度的估计误差就会加大.算法在迭代时, 对支撑集原子个数是单调递增的, 因此很容易产生过度估计问题, 导致重建误差变大.
为解决上述问题, 本文根据重建残差与支撑集原子个数及信号稀疏度之间的关系, 针对由于SAMP算法采用固定增加步长所带来的精度较低和可能由此产生的过度估计问题[12], 提出了一种新的变步长稀疏度自适应子空间追踪算法(Variable step size sparsity adaptive subspace pursuit, VssSASP)算法.该算法不需要知道信号稀疏度信息, 先采用固定步长确定信号稀疏度的范围, 再用可变步长逐步缩小该范围并确定稀疏度, 对支撑集原子选择采用子空间追踪算法的回溯方法, 仿真实验结果表明该算法信号重建运算量较小, 重建精度高.
1. 压缩感知基本理论框架
将要测量的原始信号x作为实空间xN内N维列向量.那么由于在x内的所有信号都可以用N维的基向量Ψ={Ψi|i=1, 2, · · ·, N}的线性组合来表示, 这里规定基向量均为规范且正交的, 那么x可以表示为
$ \boldsymbol{x} = \sum\limits_{i = 1}^N {{\alpha _i}} {\boldsymbol{\Psi} _i} = \Psi \alpha $
(1) 其中, x和α是同一信号的等价表示, α是信号x在Ψ域表示.如果α仅有K « N个非零系数(或远大于零的系数, 且其他系数均接近于零)时, 那么可以判定信号x为K稀疏的, 这里Ψ为N × N维变换基.
如果x是K稀疏的, 我们用一个观测矩阵Φ ∈ RM×N (M < N)对x进行观测, 得到观测向量y ∈ RM.
$ \boldsymbol{y} = \Phi \boldsymbol{x} $
(2) 这里要求Φ与Ψ是不相关的, 结合式(1)有
$ \boldsymbol{y} = \Phi \Psi \alpha = \Theta \boldsymbol{\alpha} $
(3) 其中, Θ=ΦΨ为M × N维矩阵, 定义为传感矩阵, 针对具体应用, Φ与Ψ均为确定的, 所以传感矩阵Θ是固定的.压缩感知中常用的观测矩阵有伯努利矩阵、高斯随机矩阵等.
由于测量值y的维数M小于信号x的维数N, 从而实现了信号的压缩.式(2)中不能由y求解x, 而当α是K稀疏, K « M < N, 且Θ满足有限等距性质(Restricted isometry property, RIP)时, 我们可以根据已有的稀疏分解理论, 通过求解式(3)得到稀疏系数α, 然后再根据式(1)得到重建的原始信号${\boldsymbol{\hat x}} $.上述过程即为求解最小l0范数下式(3)的最优化问题:
$ \boldsymbol{\widehat \alpha} = \arg \min {\left\| \boldsymbol{\alpha} \right\|_0}\quad {\rm{s}}.{\mkern 1mu} {\rm{t}}.\quad \boldsymbol{y} = \Theta \boldsymbol{\alpha} $
(4) 由于l0范数下求解只能通过对所有可能的稀疏情况进行求解, 所以式(4)的求解是个NP-hard问题.为解决该问题, 研究人员根据信号稀疏分解的相关理论, 寻找到了更有效的求解途径.文献[4]表明, 在一定条件下, l1最小范数与l0最小范数具有等价性, 那么式(4)可转化为
$ \boldsymbol{\widehat \alpha} = \arg \min {\left\| \boldsymbol{\alpha} \right\|_1}\quad {\rm{s}}.{\mkern 1mu} {\rm{t}}.\quad \boldsymbol{y} = \Theta \boldsymbol{\alpha} $
(5) l1最小范数可用凸优化算法求解式(5)中的α, 如内点法、BP算法等, 这类算法重构精度较高, 复杂度也高, 不适合处理维数较高数据.相对而言, 基于l0最小范数的贪婪追踪算法计算复杂度低, 算法结构简单, 但重建精度稍低, 研究快速有效的贪婪追踪算法, 是本文要研究的主要内容.
2. 变步长稀疏度自适应子空间追踪算法
CS的信号重建可以看成信号稀疏分解问题, 如果把Θ的各列向量θi (1≤i≤N)看作原子并将其作归一化处理, 那么这些列向量就组成超完备字典D.根据设定, 在α中只有K个非零元素, 那么观测值y就能用D中K个原子线性组合来表示.我们把这K个原子的索引组成集合Γ, y可以表示成ΘΓαΓ, 其中ΘΓ、αΓ分别表示根据Γ索引在对应Θi、α中的子集.
y在字典D上具有稀疏的表示, 我们的目标就是寻找最少的一组原子{νi|i ∈ Γ}, 使得残差r=y -ΘΓαΓ信号能量最小.
已有的SAMP算法及其他改进算法可以在未知信号的稀疏度K的情况下, 采用逐步扩充支撑集数量并回溯验证的方法来选择支撑集, 逐步逼近原始信号.但这些算法在每次迭代中采用增加固定数目(步长)原子的方法估计信号稀疏度, 存在计算量过大及过度估计的问题.
2.1 算法描述
本文提出的VssSASP算法在K未知的条件下, 首先通过一种测试方法确定固定步长值, 再将固定步长与可变步长相结合, 通过残差变化精确估计K值, 在每个阶段的迭代均采用SP算法迭代选择支撑集原子.算法分为以下三个阶段:
1)确定固定迭代步长.在文献[16]中给出了一种稀疏度估计方法, 通过测试找到满足
$ {\left\| {\Phi _{{{\bf{\Gamma }}_0}}^{\rm{T}}\boldsymbol{y}} \right\|_2} < \frac{{1 - {\delta _k}}}{{1 - {\delta _k}}}{\left\| \boldsymbol{y} \right\|_2} $
(6) 的最大K0值.其中, δk为传感矩阵Θ的K阶有限等距常数(Restricted isometry constant, RIC), Γ0为$ {\left\| {\Phi _{{\boldsymbol{\Gamma} _0}}^T\boldsymbol{y}} \right\|_2} < \frac{{1 - {\delta _k}}}{{1 - {\delta _k}}}{\left\| \boldsymbol{y} \right\|_2}$中前K0个最大值的索引. K0为接近且小于信号稀疏度值K, 我们可以把它作为固定迭代步长.
2)确定K的范围.以K0作为固定步长, 分步增加信号支撑集的个数L=nK0, 并以SP算法迭代计算出相应信号支撑集及残差, 当L < K时, 增加的支撑集含有效原子, 所以残差会迅速减小; 而当L≥K后, 新增加的支撑集一般不含有效原子, 残差基本不变甚至变大.所以根据残差在L=nK0各阶段的变化, 可以粗略确定K的范围(在一个固定步长K0区间内).
3) K值逐步逼近.采用二分法, 取信号支撑集个数为上个阶段确定K的范围中间值, 再利用SP算法求出当前支撑个数下的残差值, 将该残差值与K的范围两端对应的残差值对比, 根据若残差基本不变, 则对应范围的信号支撑集个数大于K的原则, 逐步缩小对K值的估计, 直到最终确定K值, 并求出对应支撑集及其系数, 进而实现原信号的重建.
2.2 算法步骤
VssSASP算法具体步骤如下.流程中, rt表示残差, t表示迭代次数, Ø表示空集, Γt表示t次迭代的索引(列序号)集合(元素个数为L), θj表示矩阵Θ的第j列, Θt={θj} (∀j ∈ Ct)表示按索引集合Ct选出的矩阵Θ的列集合(列数为Lt), αt为Lt × 1的列向量, 符号∪表示集合并运算, $\left\langle { \cdot, \cdot } \right\rangle $表示求向量内积, abs[·]表示求模值(绝对值).
输入. M维测量向量y, M × N维传感矩阵Θ, 判定阈值ε.
输出.信号x的稀疏表示系数估计${\boldsymbol{\hat \alpha} }$.
步骤1.确定固定迭代步长.
步骤1.1.初始化: K0=1, 计算u=abs[ΘTy] (即计算$ \left\langle {\boldsymbol{y}, {\boldsymbol{\theta} _j}} \right\rangle $, 1≤j≤N);
步骤1.2.选择u中K0个最大值, 将这些值对应Θ的列序号j构成(列序号集合) Γ0;
步骤1.3.如果${\left\| {\Phi _{{\boldsymbol{\Gamma} _0}}^{\rm{T}}\boldsymbol{y}} \right\|_2} < \frac{{1 - {\delta _k}}}{{1 - {\delta _k}}}{\left\| \boldsymbol{y} \right\|_2} $, 则K0=K0 + 1, 转至步骤1.2;
步骤1.4.得到稀疏度估计K0=K0 -1, 固定迭代步长step0=K0.
步骤2.确定K的范围.
步骤2.1.初始化: r0=rs0=y, Γ0=Ø, L=S=step0, k=1, t=1;
步骤2.2.计算u=abs[ΘTrt-1], 选择u中L个最大值, 将这些值对应Θ的列序号j构成(列序号集合) St;
步骤2.3.令ct=Γt-1 ∪ St, Θt={θj} (∀ j ∈ Ct);
步骤2.4.求y=Θtαt的最小二乘解:
$ \left\| r \right\|_2^2 = \left\| {y - {\Theta _{{{\bf{\Gamma }}_{_L}}\alpha {{\bf{\Gamma }}_{_L}}}}} \right\|_2^2 = \\ \;\;\;\;\;\;\left\| {{\Theta _{{\bf{\Gamma }}\alpha {\bf{\Gamma }}}} - {\Theta _{{{\bf{\Gamma }}_{_L}}\alpha {{\bf{\Gamma }}_{_L}}}}} \right\|_2^2 \approx \left\| {{\Theta _{{{\bf{\Gamma }}_R}\alpha {{\bf{\Gamma }}_R}}}} \right\|_2^2 $
步骤2.5.从${{\boldsymbol{\hat \alpha} }_t} $中选出绝对值最大的L项记作$ {{\boldsymbol{\hat \alpha} }_{tL}}$, 对应的Θt中的L列记为ΘtL, 对应的Θ的列序号为ΓtL;
步骤2.6.更新残差rnew=y -ΘtL(ΘT tLΘtL)-1ΘT tLy;
步骤2.7.如果${\left\| {{\boldsymbol{r}_{{\rm{new}}}}} \right\|_2} < {\left\| {{\boldsymbol{r}_{t - 1}}} \right\|_2} $, Γt=ΓtL, rt=rnew, t=t + 1, 返回步骤2.2;
步骤2.8.否则rsk=rt-1, pak=Γt, 如果${\left\| {\boldsymbol{r}{\boldsymbol{s}_{k - 1}}} \right\|_2} - {\left\| {\boldsymbol{r}{\boldsymbol{s}_k}} \right\|_2} \ge \varepsilon $, L=L + S, k=k + 1, 返回步骤2.2, 继续迭代; 否则完成该阶段计算, 确定K值估计区间为[L-2S, L-S].
步骤3.变步长逐步逼近.
步骤3.1.初始化: rsl=rsk-2, rsu=rsk-1, Kl=L -2S, Ku=L -S, r0=rsl, Γ0=pak-2, L=round((Kl + Ku)/2), t=1;
步骤3.2.同步骤2.2～2.7, 利用SP算法计算出在支撑集为L最优估计时的残差rsm;
步骤3.3.如果${\left\| {\boldsymbol{rsm}} \right\|_2} - {\left\| {\boldsymbol{rsu}} \right\|_2} < \varepsilon $, 则Ku=L, rsu=rsm, L=round((Kl + Ku)/2), 返回步骤3.2;否则Kl=L, rsl=rsm, L=round((Kl + Ku)/2), 返回步骤3.2.
步骤3.4.重复步骤3.2～3.3, 不断缩小K值范围[Kl, Ku], 直到最终确定K值和相应的支撑集及其系数以及残差.
算法中阈值ε用于增大支撑集个数时, 根据残差判断恢复精度是否提高.没有环境噪声, 当信号为理想K稀疏(α仅有K个非零系数)时, ε可取任意小于α中最小系数的值; 当信号为近似稀疏(α的系数按大小满足指数衰减)时, ε取期望的重建精度, ε越小, 重建所需的支撑集越多, 重建精度越高.如果存在环境噪声, 则阈值ε与信号噪声水平、能量等有密切关系, 根据大量实验结果, ε取${K_0}{10^{ - SNR/20}}{\left\| \boldsymbol{y} \right\|_2} $, K0取0.12时重建效果较好.
2.3 算法分析
对于稀疏度为K的原始信号, 进行压缩采样, 信号重建时, 若选择的支撑集原子个数L < K, 设α的K个稀疏系数中较大的L个对应的索引集为ΓL, 其余较小的K -L个稀疏系数对应的索引集为ΓR, 则各种重建算法的最优结果是选择集合ΓL对应的Θ的各列向量作为支撑集, 且残差的能量最小为
$ \begin{array}{l} \left\| {\bf{r}} \right\|_2^2 = \left\| {{\bf{y}} - {\Theta _{{{\bf{\Gamma }}_L}}}{{\bf{\alpha }}_{{{\bf{\Gamma }}_L}}}} \right\|_2^2 = \\ \;\;\;\;\;\;\;\;\;\left\| {{\Theta _{\bf{\Gamma }}}{{\bf{\alpha }}_{\bf{\Gamma }}} - {\Theta _{{{\bf{\Gamma }}_L}}}{{\bf{\alpha }}_{{{\bf{\Gamma }}_L}}}} \right\|_2^2 \approx \left\| {{\Theta _{{{\bf{\Gamma }}_R}}}{{\bf{\alpha }}_{{{\bf{\Gamma }}_{\bf{R}}}}}} \right\|_2^2 \end{array} $
(7) 由于Θ的各列向量作为支撑集满足RIP条件, 则
$ \left\| {\left( {1 - {\delta _k}} \right){{\bf{\alpha }}_{{{\bf{\Gamma }}_R}}}} \right\|_2^2 \le \left\| {{\Theta _{{{\bf{\Gamma }}_R}}}{{\bf{\alpha }}_{{{\bf{\Gamma }}_R}}}} \right\|_2^2 \le \left\| {\left( {1 + {\delta _k}} \right){{\bf{\alpha }}_{{{\bf{\Gamma }}_R}}}} \right\|_2^2 $
一般δk « 1, 残差的能量主要取决于$\left\| {{\boldsymbol{\alpha} _{{\boldsymbol{\Gamma} _R}}}} \right\|_2^2 $, 而αΓR由K -L个较小的稀疏系数组成, 所以随着L的增加, 重建残差迅速减小.
类似分析可知, 当支撑集原子个数L≥K时, 则成功重建时K个稀疏系数都会被选中, 那么残差能量为
$ \left\| \boldsymbol{r} \right\|_2^2 = \left\| {{\Theta _{{{\bf{\Gamma }}_0}}}{\boldsymbol{\alpha} _{{{\bf{\Gamma }}_0}}}} \right\|_2^2 $
(8) 其中, Γ0为α中除了K个稀疏系数外, 其他为零或接近于零系数的索引集, 共有N -L个元素.很明显这时随着L增加, 重建残差将保持为一个很小的值, 且基本稳定在较小范围内.不同支撑集原子个数L与重建残差大小的变化趋势如图 1所示, 此即为本文算法用来估计信号稀疏度的基本依据.
图 1 支撑集原子个数与残差对应关系示意图Fig. 1 The relationship of residual error and the number of support set atom在本文算法的步骤2, 确定残差变化拐点的大致位置, 在步骤3用二分法确定拐点准确位置, 该位置即对应稀疏度的估计.我们计算得到残差接近最小的情况下L的最小值, 计算过程中, 可知rsk能量单调递减, 算法至少收敛到一个局部最小点.
在SAMP及其他稀疏度自适应算法中, 通过递增方式估计K值, 步长和运算量与K的准确估计存在矛盾, 特别是有噪声存在时, K值估计误差大, 恢复精度不高.而VssSASP算法通过恢复残差的变化估计稀疏度K值大小范围, 可以得到比较准确的K值估计, 具有较高的恢复精度.
关于运算量方面, 在确定固定迭代步长阶段主要运算为一次内积运算及小于K次范数运算, 计算量较小; 确定K的范围阶段主要运算为K/K0 + 1次SP运算(K0略小于K), 变步长阶段逐步逼近采用二分法, 主要运算为log2 K0次SP运算, 后两个阶段计算量较大, 总的运算复杂度上限为O(MNK log2 K).除了第1次, 后面每次SP运算可以利用上次的SP运算结果, 所以单次SP运算迭代次数均很少; 与SAMP算法相比, 当K很大时, VssSASP算法的运算复杂度小于SAMP算法运算复杂度O(MNK2).
3. 仿真实验结果及分析
为验证本文提出的VssSASP算法对各种原始信号的重建效果及性能, 对算法进行了一系列仿真验证实验, 将典型贪婪追踪算法OMP、ROMP、SP以及SAMP算法与本文算法进行了性能仿真比较.其中SAMP算法由于不同步长相应的算法性能不同, 所以这里取步长值s为1, 5和10三个值分别进行测试.
这些实验均是在惠普g14笔记本(4 GB DDR3内存, i5-4200U)上运行, 采用的仿真软件版本为Matlab R2009a, 若非特殊说明, 对于不同的数据点均为运行500次的平均值.
3.1 稀疏度估计实验
为验证本文算法对稀疏度K估计的准确度, 通过仿真实验加以测试.实验中, 取N=256, M=128, 观测矩阵Φ为M × N阶独立分布、零均值、单位方差的高斯随机矩阵, 原始信号x为直接构造的K稀疏信号(从x中随机取K个元素, 每一项值为独立分布、零均值、单位方差高斯随机变量, 其他元素值为零), 所以稀疏基Ψ为单位阵, 通过y=Φx得到观测向量y.实验中对原始信号x叠加了不同强度的高斯白噪声.
图 2为不同噪声环境下对稀疏度K值的估计误差, 图中横坐标表示原始信号x的实际稀疏度K, 纵坐标表示K的估计误差.在信噪比较小时, 由于加入的噪声信号较大, 算法会将能量较大的部分噪声当成实际信号, 所以估计误差较大; 而随着信噪比增加, 噪声减小, 估计误差也逐渐变小, 当没有噪声时, 本文算法可以准确估计信号的实际稀疏度K.
图 2 不同信噪比下稀疏度K估计误差Fig. 2 Estimating error of signal sparsity under different noise level3.2 信号重建性能实验
将本文的VssSASP算法与OMP、ROMP、SP、SAMP其他典型贪婪追踪算法进行了性能仿真比较, 检验各算法在准确重建概率、重建精度及重建运行时间方面的对比.实验仿真中基本设置与前面相同, 并根据仿真需要改变其中部分参数.
首先对信号的准确重建概率进行比较, 这里信号准确重建定义为在无噪声的理想情况下, 实际信号x与恢复信号$ {\boldsymbol{\hat x}}$中非零元素的位置相同, 且误差的能量小于某一个阈值, 这里阈值取10-6.
图 3 (a)中给出了K=20时, 不同采样点M的信号准确重建率.从图 3 (a)可以看出, 对于所有重建算法, 原始信号的准确重建概率均随着采样点数M的增加而增大; 对于本文算法, 当M > 80时, VssSASP算法重建概率接近于1;当M > 60时, VssSASP算法重建概率明显超过OMP、ROMP和SP算法, 与SAMP算法重建概率相当.相对而言, 对于同样的信号, VssSASP算法能稳定重建信号所需的采样点数较少.
图 3 不同算法信号准确重建率Fig. 3 Exact recovery ratio of recovery algorithms by different reconstruction algorithms图 3 (b)中给出了M=128时, 各重建算法对不同稀疏度K信号的准确重建率.同样可以看出, 信号准确重建概率随着信号稀疏度K的增大而减小, 当K < 45时, VssSASP算法重建概率接近于1;当K > 45时, 重建概率开始明显下降, 但重建概率明显高于OMP、ROMP和SP算法, 与SAMP算法重建概率相当.相对而言, 对于同样的采样点数, VssSASP算法能稳定重建稀疏度更高的信号.
值得注意的是, 在采样点较少或稀疏度较大时, VssSASP算法重建概率要低于SAMP算法.这两种算法都是以SP算法为基础, SAMP算法以线性增加稀疏度K估计的方式进行搜索, 尽管此时单次SP算法重建成功概率较低, 但在搜索过程中只要有一次重建成功, 即认为SAMP算法重建完成, 所以这时SAMP算法重建成功率远高于SP算法, 且步长s越小成功率越高, 但实际上此时算法估计稀疏度一般要明显高于原信号实际稀疏度K. VssSASP算法采用残差判定稀疏度K的区间, 然后再进行二分法搜索, 搜索次数较少, 且此时残差变化的不规律也可能导致区间定位错误, 所以算法重建成功率提高有限, 低于SAMP算法.
综上所述, VssSASP算法在无噪声环境下重建概率超过OMP、ROMP和SP算法, 而与SAMP算法性能相近, 信号重建概率较高.
图 4给出了在原始信号x叠加了不同强度的高斯白噪声情况下, 各算法对信号的重建精度.从图 4可以看出, 各算法的重建精度均随着信噪比的增大而增大, 其中VssSASP算法重建精度略优于OMP、SP算法, 而明显优于ROMP及各种不同步长的SAMP算法.对比SAMP算法, 在低信噪比(小于6dB)时, VssSASP算法重建精度约高3dB, 在高信噪比(大于6dB)时, 重建精度约高1.8dB. VssSASP算法在不同噪声环境下重建精度高于稀疏度自适应算法SAMP, 而与非自适应算法OMP、SP算法相当; 相比原信号, 重建后的信号信噪比要高约3～5dB, 说明重建过程具有一定的去噪功能.
图 4 噪声环境下的重建精度(M=80, N=256, K=20)Fig. 4 Reconstruction accuracy in noise environment (M=80, N=256, K=20)图 5是各种算法运行时间的对比, 为体现对高K值信号的重建性能, 这里N和M值均取较大值.从图 5可以看出, 各算法的运算时间均随着K的增大而增加.总体而言, 由于ROMP和SP算法已知稀疏度K, 每次迭代时可选择支撑集原子个数多, 迭代次数少, 运算时间最少; OMP算法每次迭代只增加支撑集中的一个原子, 运算时间较多; 对于SAMP和本文的稀疏度自适应算法, 在K较小时, VssSASP算法的运行时间与SAMP算法相差不大, 但随着K的增加, SAMP算法运算时间迅速增加(步长越小, 增加越多), 而VssSASP算法变步长优势逐渐体现, 运算速度很快超过SAMP算法.
图 5 算法运行时间对比(M=1 024, N=2 048)Fig. 5 Running time by different reconstruction algorithms (M=1 024, N=2 048)3.3 图像重建实验
在实验中采用的是大小为256像素× 256像素的Lena灰度图像, 因为图像各列(行)的稀疏性较差, 为提高重建效果, 这里先将其用双正交小波(bior3.7)进行变换, 然后将变换后矩阵的各行作为一维列信号进行压缩采样, 再用不同算法进行重建, 重建完成后再通过同样的小波反变换得到重建图像.其中观测矩阵Φ采用正交观测矩阵, 稀疏基Ψ仍为单位阵.
定义图像的采样率为M/N, 重建图像的质量用峰值信噪比(Peak signal to noise ratio, PSNR)来表示.图 6给出了在采样率为0.5的条件下, Lena图像的原图像及OMP, ROMP, SP, SAMP, VssSASP算法重建图像和PSNR值.从图 6可以看出, VssSASP算法的重建质量略优于ROMP算法, 优于OMP、SP及各种步长的SAMP算法, 经VssSASP算法重建后的图像与原始图像比较接近, 且细节保持较好, 具有良好的视觉效果.
图 6 样率0.5时, Lena原图像及各算法重建图像Fig. 6 Original image and reconstructed image of different algorithms for Lena when sampling rate is 0.5表 1给出了不同算法在不同压缩比条件下重建图像的PSNR和运算时间对比, 每个数据均为50次实验的平均结果.
表 1 各算法的重建质量及运行时间对比Table 1 Comparison of the qualities of images reconstructed and running time by different algorithms重建算法 M=0.3 × N M=0.4 × N M=0.5 × N PSNR(dB) 运算时间 PSNR(dB) 运算时间 PSNR(dB) 运算时间 OMP 23.07 0.69 27.48 1.06 30.96 1.24 ROMP 24.72 0.32 27.89 0.68 31.52 1.05 SP 23.15 0.76 26.51 0.96 30.84 1.32 SAMP1 24.80 11.57 27.79 18.25 30.96 26.96 SAMP5 24.14 2.71 27.49 4.55 30.85 7.46 SAMP10 23.65 1.55 27.42 2.84 30.72 3.91 VssSASP 25.26 3.21 28.13 4.92 31.63 7.63 从表 1可以看出, 随着采样率的增大, 各种重建算法重建图像的PSNR均显著增大, 通过增加观测数目可以提高重建图像的重建质量; 在相同压缩比的条件下, VssSASP算法重建图像的PSNR值均为最大, 说明该算法对图像重建精度较高.
另外, 根据表 1中时间数据可以看出, 随着采样率的增大, 各种重建算法运算时间均增大, 且稀疏度自适应算法SAMP和VssSASP运算时间明显高于其他非稀疏度自适应算法.其中ROMP算法在每次迭代时会选取多个原子, 其运算时间少于OMP算法; SAMP算法要进行多次类SP运算, 其运行时间高于SP算法, 且步长s越小, 运算时间越长; VssVSASP算法也是进行多次类SP运算, 由于我们进行实验时, 分别对每一行进行压缩重建, N值较小, 算法变步长优势不明显, 算法运算时间低于步长为1的SAMP算法, 而高于步长为10的SAMP算法.
4. 结论
本文提出了一种变步长稀疏度自适应子空间追踪算法VssSASP, 可在未知原始信号稀疏度的先验信息的情况下准确重建信号.算法采用固定步长与可变步长相结合的方式, 采用子空间追踪通过迭代估计的更新, 根据各阶段残差变化确定信号稀疏度, 可精确重建原信号.实验结果表明, VssSASP算法可以高概率重建稀疏信号, 在噪声环境下具有高重建精度, 当信号稀疏度较大时, 算法的运算速度也高于SAMP稀疏度自适应算法.
-
图 1 支撑集原子个数与残差对应关系示意图
Fig. 1 The relationship of residual error and the number of support set atom
图 2 不同信噪比下稀疏度K估计误差
Fig. 2 Estimating error of signal sparsity under different noise level
图 3 不同算法信号准确重建率
Fig. 3 Exact recovery ratio of recovery algorithms by different reconstruction algorithms
图 4 噪声环境下的重建精度(M=80, N=256, K=20)
Fig. 4 Reconstruction accuracy in noise environment (M=80, N=256, K=20)
图 5 算法运行时间对比(M=1 024, N=2 048)
Fig. 5 Running time by different reconstruction algorithms (M=1 024, N=2 048)
图 6 样率0.5时, Lena原图像及各算法重建图像
Fig. 6 Original image and reconstructed image of different algorithms for Lena when sampling rate is 0.5
表 1 各算法的重建质量及运行时间对比
Table 1 Comparison of the qualities of images reconstructed and running time by different algorithms
重建算法 M=0.3 × N M=0.4 × N M=0.5 × N PSNR(dB) 运算时间 PSNR(dB) 运算时间 PSNR(dB) 运算时间 OMP 23.07 0.69 27.48 1.06 30.96 1.24 ROMP 24.72 0.32 27.89 0.68 31.52 1.05 SP 23.15 0.76 26.51 0.96 30.84 1.32 SAMP1 24.80 11.57 27.79 18.25 30.96 26.96 SAMP5 24.14 2.71 27.49 4.55 30.85 7.46 SAMP10 23.65 1.55 27.42 2.84 30.72 3.91 VssSASP 25.26 3.21 28.13 4.92 31.63 7.63 
下载: 导出CSV
-
[1] Donoho D L. Compressed sensing. IEEE Transactions on Information Theory, 2006, 52(4):1289-1306 doi: 10.1109/TIT.2006.871582 [2] Candés E J. Compressive sampling. In:Proceedings of the 2006 International Congress of Mathematicians. Madrid, Spain:EMS Publishing House, 2006. 1433-1452 http://www.researchgate.net/publication/41537648_Compressive_Sampling [3] Candés E J, Tao T. Near-optimal signal recovery from random projections:universal encoding strategies. IEEE Transactions on Information Theory, 2006, 52(12):5406-5425 doi: 10.1109/TIT.2006.885507 [4] 荆楠, 毕卫红, 胡正平, 王林.动态压缩感知综述.自动化学报, 2015, 41(1):22-37 http://www.aas.net.cn/CN/abstract/abstract18580.shtmlJing Nan, Bi Wei-Hong, Hu Zheng-Ping, Wang Lin. A survey on dynamic compressed sensing. Acta Automatica Sinica, 2015, 41(1):22-37 http://www.aas.net.cn/CN/abstract/abstract18580.shtml [5] Mallat S G, Zhang Z F. Matching pursuits with time-frequency dictionaries. IEEE Transactions on Signal Processing, 1993, 41(12):3397-3415 doi: 10.1109/78.258082 [6] Tropp J A, Gilbert A C. Signal recovery from random measurements via orthogonal matching pursuit. IEEE Transactions on Information Theory, 2007, 53(12):4655-4666 doi: 10.1109/TIT.2007.909108 [7] Needell D, Vershynin R. Uniform uncertainty principle and signal recovery via regularized orthogonal matching pursuit. Foundations of Computational Mathematics, 2009, 9(3):317-334 doi: 10.1007/s10208-008-9031-3 [8] Donoho D L, Tsaig Y, Drori I, Starck J L. Sparse solution of underdetermined systems of linear equations by stagewise orthogonal matching pursuit. IEEE Transactions on Information Theory, 2012, 58(2):1094-1121 doi: 10.1109/TIT.2011.2173241 [9] Needell D, Tropp J A. CoSaMP:iterative signal recovery from incomplete and inaccurate samples. Communications of the ACM, 2010, 53(12):93-100 [10] Dai W, Milenkovic O. Subspace pursuit for compressive sensing signal reconstruction. IEEE Transactions on Information Theory, 2009, 55(5):2230-2249 doi: 10.1109/TIT.2009.2016006 [11] Do T T, Gan L, Nguyen N, Tran T D. Sparsity adaptive matching pursuit algorithm for practical compressed sensing. In:Proceedings of the 42nd Asilomar Conference on Signals, Systems, and Computers. Pacific Grove, California:IEEE Press, 2008. 581-587 [12] 沈燕飞, 李锦涛, 朱珍民, 张勇东, 代锋.基于非局部相似模型的压缩感知图像恢复算法.自动化学报, 2015, 41(2):261-272 http://www.aas.net.cn/CN/abstract/abstract18605.shtmlShen Yan-Fei, Li Jin-Tao, Zhu Zhen-Min, Zhang Yong-Dong, Dai Feng. Image reconstruction algorithm of compressed sensing based on nonlocal similarity model. Acta Automatica Sinica, 2015, 41(2):261-272 http://www.aas.net.cn/CN/abstract/abstract18605.shtml [13] 伍飞云, 周跃海, 童峰.基于似零范数和混合优化的压缩感知信号快速重构算法.自动化学报, 2014, 40(10):2145-2150 http://www.aas.net.cn/CN/abstract/abstract18489.shtmlWu Fei-Yun, Zhou Yue-Hai, Tong Feng. A fast sparse signal recovery algorithm based on approximate l0 norm and hybrid optimization. Acta Automatica Sinica, 2014, 40(10):2145-2150 http://www.aas.net.cn/CN/abstract/abstract18489.shtml [14] 徐泽芳, 刘顺兰.一种自适应正则化子空间追踪算法.计算机工程与应用, 2015, 51(3):208-211 http://youxian.cnki.com.cn/yxdetail.aspx?filename=JSGG20140123002&dbname=CAPJ2014Xu Ze-Fang, Liu Shun-Lan. Adaptive regularized subspace pursuit algorithm. Computer Engineering and Applications, 2015, 51(3):208-211 http://youxian.cnki.com.cn/yxdetail.aspx?filename=JSGG20140123002&dbname=CAPJ2014 [15] 张成, 沈川, 程鸿, 章权兵, 陈岚, 韦穗.彩色全息压缩重构.自动化学报, 2015, 41(2):419-428 http://www.aas.net.cn/CN/abstract/abstract18620.shtmlZhang Cheng, Shen Chuan, Cheng Hong, Zhang Quan-Bing, Chen Lan, Wei Sui. Compressed reconstruction of color holography. Acta Automatica Sinica, 2015, 41(2):419-428 http://www.aas.net.cn/CN/abstract/abstract18620.shtml [16] 杨成, 冯巍, 冯辉, 杨涛, 胡波.一种压缩采样中的稀疏度自适应子空间追踪算法.电子学报, 2010, 38(8):1914-1917 http://www.cnki.com.cn/Article/CJFDTOTAL-DZXU201008031.htmYang Cheng, Feng Wei, Feng Hui, Yang Tao, Hu Bo. A sparsity adaptive subspace pursuit algorithm for compressive sampling. Acta Electronica Sinica, 2010, 38(8):1914-1917 http://www.cnki.com.cn/Article/CJFDTOTAL-DZXU201008031.htm 期刊类型引用(8)
1. 杜秀丽 ，张文龙 ，邱少明 ，刘庆利 . 基于支撑集先验的脑电信号正则化子空间重构. 计算机应用与软件. 2022(04): 105-109+148 . 
百度学术2. 田金鹏，闵天，薛莹，郑国莘. 自适应线性预测卡尔曼滤波压缩感知算法. 控制与决策. 2020(01): 83-90 . 百度学术3. 王可心，韩太林，高杨，王啸. 基于改进子空间追踪算法的冲击波信号采集. 长春理工大学学报(自然科学版). 2019(01): 89-94 . 百度学术4. 孙润润. 基于稀疏度区间的变步长最优子空间追踪算法. 计算机技术与发展. 2019(02): 48-53 . 百度学术5. 朱华，岳峻，李振波，张志旺，寇光杰. 基于l_(2, 1)范数原子选择的图像分块稀疏重构. 计算机应用研究. 2019(05): 1560-1563 . 百度学术6. 谢瑶，李思敏，唐智灵. 基于压缩感知和循环谱理论的信号载频估计方法. 计算机应用研究. 2019(09): 2760-2763+2768 . 百度学术7. 姚万业，姚吉行. 一种稀疏度自适应广义正交匹配追踪算法. 仪器仪表用户. 2018(08): 16-20 . 百度学术8. 孙娜，刘继文，肖东亮. 基于BFGS拟牛顿法的压缩感知SL0重构算法. 电子与信息学报. 2018(10): 2408-2414 . 百度学术其他类型引用(6)
-
下载:
计量
- 文章访问数: 2259
- HTML全文浏览量: 552
- PDF下载量: 832
- 被引次数: 14
