迭代硬阈值压缩感知重构算法――IIHT

收稿日期:2011-01-26;修回日期:2011-03-07。基金项目:广东省自然科学基金资助项目（9151170003000017）。

作者简介:张宗念（1963-），男，河北深州人,副教授，博士，主要研究方向：压缩感知、图像分析与处理； 李金徽（1980-）,男,辽宁沈阳人,工程师,主要研究方向:分布式计算机网络; 黄仁泰（1964-），男，广东东莞人,副教授，主要研究方向:嵌入式系统设计。

文章编号:1001-9081（2011）08-02123-03doi:10.3724/SP.J.1087.2011.02123

（1.东莞理工学院 电子工程学院,广东 东莞523808； 2.东莞理工学院 网络中心,广东 东莞523808；

3.东莞理工学院 计算机学院,广东 东莞523808）

（）

摘 要:研究了压缩感知信号重构算法的理论,针对迭代硬阈值（IHT）重构算法对测量矩阵的过分依赖、计算复杂度高、运算时间长的缺点，通过修订迭代硬阈值重构算法的代价函数和自适应地调整迭代步长的选取原则,设计了一种迭代硬阈值重构算法――IIHT。IIHT算法显著提高了信号精确重构的概率,降低了算法的计算复杂度,进一步减少了算法的运算时间,加快了算法的收敛速度。

关键词:迭代；硬阈值；压缩感知

中图分类号: TN911.73文献标志码:A

IIHT: New improved iterative hard thresholding algorithm for compressive sensing

ZHANG Zong-nian1, LI Jin-hui2, HUANG Ren-tai3

（1. School of Electronic Engineering, Dongguan University of Technology, Dongguan Guangdong 523808, China;

2. Network Center, Dongguan University of Technology, Dongguan Guangdong 523808, China;

3. School of Computer Science, Dongguan University of Technology, Dongguan Guangdong 523808, China）

Abstract: To overcome the shortcomings of the overdependence on the measurement matrix, the high computation complexity, the long computation time of the Iterative Hard Thresholding （IHT） algorithm, a new improved iterative hard thresholding （IIHT） algorithm was proposed by studying the theory of signal reconstruction for compressive sensing. It improved the cost function and the selection method of step size for the IHT algorithm. The simulation results show that the proposed algorithm increases the probability of recovery and the speed of convergence and reduces the computational complexity and time.

Key words: iteration; hard thresholding; compressive sensing

0 引言

由Donoho［1-3］和Candès等人［4-7］提出了一种新颖的信号分析理论――压缩感知理论,其基本思想是:只要信号本身或在某个变换域上是稀疏的或可压缩的，那么就可用一个与变换矩阵不相关的测量矩阵将变换域的高维信号投影到一个低维空间上，然后通过求解一个优化问题就能从少量的投影中以高概率或精确地重构出原信号。压缩感知理论的主要研究内容包括信号的稀疏变换、测量矩阵设计和信号重构算法设计等。压缩感知信号重构就是研究从M个测量值中恢复出原信号x的过程。不幸的是，该过程的求解是一个组合的复杂问题,没有一个通用高效的方法。通常都是寻找其次优解。在满足某些条件下,已经有许多重构算法可以保证能够找到次优解,如正则化正交匹配跟踪（Regularized Orthogonal Matching Pursuit, ROMP）法［8］、压缩抽样匹配跟踪（Compressived Sensing Matching Pursuit, CoSaMP）法［9］和子空间跟踪（Space Pursuit, SP）法［10］，这些算法与基跟踪（Basis Pursuit, BP）算法［11］一样可以保证收敛，计算复杂度与正交匹配跟踪（Orthogonal Matching Pursuit, OMP）法［12］相当。互补匹配跟踪［13］和梯度跟踪［14］也取得了不错的性能。迭代硬阈值（Iterative Hard Thresholding, IHT）算法［15-17］在理论上推导出了IHT算法收敛的理论依据,但是为了保证算法的稳定性,要求测量矩阵除了满足有限等距特性（Restricted Isometry Property, RIP）条件以外［18-21］,还要满足其最大奇异值小于1;这样就造成了算法对测量矩阵的尺度缩放十分敏感,依赖性很强,而且运算复杂度较高,计算时间较长。为了克服算法的这些不足,本文通过改进IHT算法的代价函数和自适应调整迭代步长的选取原则,设计了一种迭代硬阈值压缩感知重构算法――IIHT（Improved Iterative Hard Thresholding）。

1 IHT算法

对于N维任意矢量x∈RN，定义y∈RM是x的M个测量值，即yΦx，其中Φ∈RM×N称为测量矩阵。IHT算法就是求解满足约束条件x0≤K的miny－Φx22优化问题。文献［15-17］推导出了迭代硬阈值算法理论,为该领域的创造性研究奠定了基础。IHT算法十分简洁,采用迭代公式（1）：

xn+1HK（xn+ΦT（y－Φxn））（1）

其中HK（a）是一个非线性算子,它将矢量a中幅度最大的前K个元素以外的所有元素设置为零。

IHT算法的计算复杂度估计。若测量矩阵是满足RIP的一般矩阵,则计算复杂度为O（MN）,很高。若测量矩阵是通过傅里叶变换矩阵或小波变换矩阵后构成的结构化矩阵,则计算复杂度为O（N log M）。

测量矩阵选取。若定义测量矩阵Φ∈RM×N和集合I{1,2,…,N}，矩阵ΦI由I中下标为i的列向量构成,i∈I,I≤K和任意矢量q∈RI,ΚcK ln （N/K）,而且Φ的元素是从方差为1/M的高斯独立同分布中选取的;或者Φ的元素是等概率取1/M或－1/M,那么测量矩阵Φ就能够以高概率具有小的有限等距常数,也可以计算出常数c和δK的值。

IHT算法的收敛性讨论。

定理1 如果测量矩阵满足Φ2

此定理确保了每一次的迭代误差会逐步减小,第n次迭代估计值xn最终收敛到一个点,即保证了IHT算法的稳定性。然而,它没有给出该估计值与原始信号的近似程度有多大。

定理2 若给定含有噪声e的测量值yΦx+e,其中x是任意矢量。定义xK是x的最佳K项近似值。如果Φ满足有限等距特性而且δ3K

尽管该定理保证了算法能够收敛到任意矢量x的最佳K项近似值的近邻,但是它还不能保证收敛于一个点。x的近似值距离其全局最优值有多远取决于K个稀疏矢量近似x的程度有多好。

2 迭代硬阈值算法――IIHT

算法的简要实现过程如下:

1）有关矢量初始化,令x10,Γ 1supp（HK（ΦTy））。

2）设定迭代顺序n1,每次迭代结束后令nn+1,直到满足终止条件,迭代结束。

3）迭代过程。首先求出gnΦT（y－Φxn）,n+1HK（xn+sngn）,Γ n+1supp（ n+1）,sn（gTΓ ngΓ n）/（gTΓ nΦTΓ nΦΓ ngΓ n）,如果Γ n+1Γ n,则有xn+1n+1；否则,如果Γ n+1≠Γ n,进一步判断sn≤（1-c）（n+1－xn22）/（Φ（n+1－xn）22））是否成立,如果成立,则xn+1n+1;否则,继续迭代直到sn≤（1-c）（n+1－xn22）/（Φ（n+1－xn）22））为止。然后修改参数snsn/（k（1-c））,n+1HK（xn+sngn）;Γ nsupp（n+1）,xn+1n+1。

IIHT算法主要作了两个方面的修改。

1）修正了代价函数。IHT算法实际上就是求约束条件为0≤K的代价函数y－Φn+122的最小值。此优化问题等价于函数：

y－Φ22+－Φxn22－Φ（－xn）22（2）

当Φ2

xn+1HK（xn+sΦT（y－Φxn））（3）

2）在迭代过程中自适应地不断改变和调整步长s的值。首先定义Γ n是xn的支撑集,gnΦT（y－Φxn）是由当前估计值xn决定的y－Φxn22的负梯度。假定已经确定了支撑集,Γ n是y的K项最佳近似值的支集,目标是优化y－ΦΓ nxΓ n22。采用梯度下降算法,通过不断迭代xn+1Γ nxnΓ n+sΦTΓ n（y－ΦΓ nxnΓ n）来实现。在支撑集固定的情况,可以计算出最佳步长,即每次迭代时最大限度地减小误差所对应的步长,每次迭代的最佳步长计算公式如式（4）：

sn（gTΓ ngΓ n）/（gTΓ nΦTΓ nΦΓ ngΓ n）（4）

其中：gΓ n是舍弃掉g中对应于Γ n的元素以外的元素构成的矢量;ΦΓ n是舍弃Φ中对应于Γ n的元素以外的对应列矢量构成的矩阵。这里，在每一次迭代时首先按照式（4）计算步长s和按照式（1）计算出新的待定值n+1。如果n+1和xn的支撑集相等,那么代价函数保证有最大降幅,所以选择xn+1n+1。然而,如果n+1的支撑集与xn的不同,步长的最优值就不能保证。此时,算法收敛的充分条件是s≤ω,其中：

ωn（1-c）（n+1－xn22）/（Φ（n+1－xn）22））（5）

其中c是一个小的固定常数。首先计算ω,然后检查s

关于步长s的取值范围。步长s的取值范围由测量矩阵Φ的特性决定。对于非零元素小于2K的任意矢量x有0≤α2K≤Φx2/x2≤β2K成立,其中α2K和β2K是两个常数。由于gΓ n只有K个非零元素,所以s的取值范围由不等式（1/β22K）≤（1/β2K）≤sn≤（1/α2K）≤（1/α22K）确定,详细证明见文献［17］。

3 仿真实验

实验1 研究IIHT算法与IHT以及其他流行的压缩感知重构算法CoSaMP和l1范数最小化算法的精确重构概率,运算时间与稀疏程度之间的关系。从独立同高斯分布中选取测量矩阵Φ使得Φ21,稀疏度为K的矢量x也是从高斯独立同分布（Independent Identical Distribution, IID）中随机选取。用这几种不同的算法分别进行重复实验400次,然后比较精确重构x的概率与信号稀疏程度之间的关系，结果见图1、2。从图中可以看出,在一定范围（K/M小于40%）内,IIHT算法精确重构信号的概率比IHT算法有显著提高,也比CoSaMP算法精确重构信号的概率高。

实验2 研究噪声对几种重构算法性能的影响。由于现实世界中,精确稀疏的信号是很少的,信号的测量值通常包含噪声,所有研究有噪声的情况下算法的性能十分必要。一般情形下,测量矢量yΦx+e,其中x是稀疏度为K的待测矢量,Φ为测量矩阵,e是高斯噪声。本文选择Φx和e归一化后的比值分别为10dB和30dB两种情况,用几种不同的重构算法分别进行实验。l1范数最小化算法是求在约束条件y－Φx2≤e2下,minx1最小化对应的解。观察信号重构后的信噪比随着稀疏度的变化规律,结果见图3、4。从图中可以看出,在低信噪比10dB,信号稀疏度与测量矢量数目的比值小于20%的情形下,IIHT算法的性能超过了其他算法的性能;当稀疏度与测量矢量数目的比值增加到超过20%时,IIHT算法的性能比l1范数最小化算法的性能稍差,CoSaMP算法的性能较差,而且信号的稀疏度越小,其性能就越差。在中等信噪比30dB,信号稀疏度与测量矢量数目的比值小于33%的情况下,IIHT算法的性能超过了其他算法的性能;但是当稀疏度与测量矢量数目的比值增加到超过33%时,即测量矢量数目较少的情况下,IIHT算法的性能比l1范数最小化算法的性能略差。

图1 重构概率与稀疏度的关系

图2 运算时间与稀疏度的关系

图3 测量噪声的影响（Φx和e归一化后的比值分别为10dB）

4 结语

尽管IHT算法具有简洁及容易实现的特点,但是它对测量矩阵的尺度缩放十分敏感,依赖性很强,而且运算复杂度较高,计算时间较长,特别不适合应用到如图像信号等复杂的大尺度信号领域。通过改进代价函数和自适应调整迭代步长的选取原则,设计并实现了一种迭代硬阈值压缩感知重构算法。理论推导与实验数据证明:IIHT算法与其他算法相比,信号精确重构概率有了显著提高;计算复杂度进一步降低,运算时间得到了进一步减少,收敛速度更快;特别适宜于对大尺度信号的重构。对于含义噪声的大尺度图像信号快速重构算法及其应用（如医学图像和合成孔径雷达图像）是进一步研究的方向。

图4 测量噪声的影响（Φx和e归一化后的比值分别为30dB）

参考文献:

[1] DONOHO D L. Compressed sensing ［J］. IEEE Transactions on Information Theory, 2006, 52（4）: 1289-1306.

［2］ DONOHO D L, TSAIG Y. Extensions of compressed sensing ［J］. Signal Processing, 2006, 86（3）: 533-548.

［3］ DONOHO D L. For most large underdetermined systems of linear equations the minimal l1-norm solution is also the sparsest solution ［J］. Communications on Pure and Applied Mathematics, 2006, 59（6）: 797-829.

［4］ CANDKS E, TAO T. Decoding by linear programming ［J］. IEEE Transactions on Information Theory, 2005,51（12）: 4203-4215.

［5］ CANDKS E, ROMBERG J. Quantitative robust uncertainty principles and optimally sparse decompositions ［J］. Foundations of Computational Mathematics, 2006, 6（2）: 227-254.

［6］ CANDKS E, ROMBERG J, TAO T. Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information ［J］. IEEE Transactions on Information Theory, 2006, 52（2）: 489-509.

［7］ CANDKS E, TAO T. Near-optimal signal recovery from random projections: Universal encoding strategies？ ［J］. IEEE Transactions on Information Theory, 2006, 52（12）: 5406-5425.

［8］ NEEDELL D, VERSHYNIN R. Uniform uncertainty principle and signal recovery via regularized orthogonal matching pursuit ［J］. Foundations of Computational Mathematics, 2009, 9（3）: 317-334.

［9］ NEEDELL D, TROPP J A. CoSaMP: Iterative signal recovery from incomplete and inaccurate samples ［J］. Applied and Computational Harmonic Analysis, 2009, 26（3）: 301-321.

［10］ DAI W, MILENKOVIC O. Subspace pursuit for compressive sensing signal reconstruction ［J］. IEEE Transactions on Information Theory, 2009, 55（5）: 2230-2249.

［11］ CHEN S S, DONOHO D L, SAUNDERS M A. Atomic decomposition by basis pursuit ［J］. SIAM Journal of Science Computing, 1988, 20（1）: 33-61.

［12］ TROPP J A, GILBERT A C. Signal recovery from random measurements via orthogonal matching pursuit ［J］. IEEE Transactions on Information Theory, 2007, 53（12）: 4655-4666.

［13］ RATH G, GUILLEMOT C. Complementary matching pursuit algorithms for sparse approximation ［J］. IEEE Transactions on Signal Processing, 2008, 56（4）: 2370-2382.

［14］ BLUMENSATH T, DAVIES M E. Gradient pursuits ［J］. IEEE Transactions on Signal Processing, 2008, 56（6）: 2370-2382.

［15］ BLUMENSATH T, DAVIES M E. Iterative thresholding for sparse approximations ［J］. Journal of Fourier Analysis and Applications, 2008, 14（5）: 629-654.

［16］ BLUMENSATH T, DAVIES M E. Iterative hard thresholding for compressed sensing ［J］. Applied and Computational Harmonic Analysis, 2009, 27（3）: 265-274.

［17］ BLUMENSATH T, DAVIES M E. How to use the iterative hard thresholding algorithm ［C］// SPARS'09: Proceeding of Signal Processing with Adaptive Sparse Structure Representations. Saint-Malo, France: ［s.n.］, 2009: 33-38.

［18］ 石光明,刘丹华,高大化,等.压缩感知理论及其研究进展［J］.电子学报,2009,37（5）:1070-1081.

［19］ 吴宗亮,窦衡.一种新的最小二乘支持向量机稀疏化算法［J］.计算机应用,2009,29（6）:1559-1562.

［20］ 曾小波,魏祖宽,金在弘.协同过滤系统的矩阵稀疏性问题的研究［J］.计算机应用,2010,30（4）:1079-1082.

［21］ 姚远,刘鹏,王辉,等.基于稀疏矩阵存储的状态表压缩算法［J］.计算机应用,2010,30（8）:2157-2160.