基于精确步长的贪婪随机Kaczmarz图像重建方法

本发明涉及信号处理中的信号压缩采样和重建方法，具体讲的是一种基于精确步长的贪婪随机kaczmarz图像重建方法。

背景技术：

1、压缩感知(compressed sensing,cs)是一种近些年来兴起的信息获取、信息传输以及信息恢复的理论和方法，其基本思想和人类视觉系统的原理类似，利用人眼对图像的特定特征的敏感度来压缩图像，例如，人眼一般对于图像中的细节变化比较敏感，对图像中的平滑部分的变化不太敏感，因此cs会保留图像中的细节，对平滑部分进行压缩。压缩感知算法通常使用稀疏表示来实现压缩，稀疏表示是指对图像进行基于特征的分解，只保留重要的少数特征，其他特征则可以忽略。因此压缩后的图像通常是高度稀疏，重建的时候需要通过恢复这些重要的特征来重建图像。

2、压缩感知图像重建算法的有点在于能保证图像质量的前提下实现较高的压缩比。这使其在图像存储和传输等领域具有较好的应用前景。然而，传统压缩感知图像重建方法主要是基于凸优化理论进行重建，计算复杂度较高，因此在实际应用中需要进行改进以提高效率。稀疏图像的恢复问题可以看作一个线性系统ax＝b的求解问题。基于迭代投影的kaczmarz方法由于其信息存储量不大和运算速度较快等优点，在压缩感知图像重建领域受到人们的关注。

3、kaczmarz算法是求解线性系统ax＝b的一类重要方法，传统kaczmarz方法是一种求解相容线性系统的方法，按照行的顺序，把初始的解向量依次投影到由系数矩阵每一个行向量和所对应的观测值决定的超平面上。很明显，kaczmarz的收敛性质依赖于行顺序，如果系数矩阵a的行顺序不好的时候，kaczmarz的收敛结果非常慢，而且传统的kaczmarz方法的收敛速度很难估计。

4、另一方面kaczmarz算法求解的是最小二乘解，通常是稠密的。而压缩感知图像重建方法中需要的解向量是稀疏的，为此需要引入l1范数，从而求解线性系统的稀疏解可以转化为正则化的最小二乘问题，因此kaczmarz方法可以延伸为稀疏kaczmarz算法。后续出现的一系列kaczmarz类算法的改进主要在通过选择合适的行，使得收敛速度得到提高，例如随机稀疏kaczmarz算法的特点是，当系数矩阵的行范数相差不大的时候，随机kaczmarz将等概率选取系数矩阵的行，此时随机kaczmarz算法的收敛速率相比于传统kaczmarz有所提高，但是仍然较慢；再例如贪婪随机kaczmarz算法提出采用一种新的概率准则来选取工作行，这种概率准则结合了贪婪的思想和随机的策略，使得其收敛速度相比于随机kaczmarz算法又有所有所提高。但是对于需要获得稀疏解向量的稀疏kaczmarz算法，贪婪随机稀疏kaczmarz算法仍然不够高效，主要因为它采用正交投影，正交投影的步长是固定的，但不一定是最优的，这使得迭代不够高效。确定每一次投影的最优步长对于贪婪随机稀疏kaczmarz是一个有效改进方法，本发明提出一种基于精确步长的贪婪随机稀疏kaczmarz算法用于图像重建。

技术实现思路

1、本发明要解决的技术问题是提供一种基于精确步长的贪婪随机kaczmarz图像重建方法。为解决以上技术问题，本发明采用以下技术方案：

2、一种基于精确步长的贪婪随机kaczmarz图像重建方法，该方法包括以下步骤：

3、步骤1、将待重建的图像按列拆分；

4、步骤2、设置求解过程的迭代终止条件，软阈值函数系数，初始迭代向量；

5、步骤3、采用贪婪随机的方式选择迭代行向量确定投影方向；

6、步骤4、使用线性搜索的方法，得到精确步长，并且更新解向量，迭代次数加1，判断是否达到迭代终止条件，若是则输出最后一次迭代得到的稀疏解向量，否则转到步骤3；

7、步骤5、将计算得到的所有最终稀疏解向量分别进行反稀疏变换得到最终解向量，再按列合并，即得到重建的图像。

8、进一步的，所述步骤4具体包括以下步骤：

9、步骤4.1、确定最优步长的问题实质是解决如下优化问题：

10、

11、其中是上一次迭代中未经过软阈值函数作用之前的解向量，tk为待求的本次迭代最优步长，a是选取的工作行向量的转置，β是选取的工作行向量对应的观测值，sλ(·)是软阈值函数，||·||2是l2范数；

12、使用线性搜索的方法解决上述优化问题得到精确步长tk；

13、步骤4.2、通过如下公式更新本次迭代解向量xk：

14、

15、其中a为投影方向，tk为步长。

16、进一步的，所述步骤4.1中，使用线性搜索的方法得到精确步长的方法包括以下步骤：

17、步骤4.11、通过如下公式初始化截距b的估计值：

18、b＝atx-β

19、其中at表示向量的转置，s＝sign(b),其中sign(·)表示符号函数，即当b>0时，sign(b)＝1；当b＝0时，sign(b)＝0；当b<0时，sign(b)＝-1。

20、若s的值等于0，则中止搜索，最终的步长t＝0；若s的值为-1，则将a和b取相反数，即a＝-a,b＝-b；

21、步骤4.11中，初始化截距b的目的是让b＝-g′(0),即让b的初始值为待优化问题的目标函数的一阶导函数在t＝0处的函数值的相反数，在后续的步骤中用于确定t是否跨过折点。同时如果s的值为-1，对a和b取相反数，是为了保证最终找到的步长t>0。

22、步骤4.12、记录a中非零项的索引，并且利用这个索引将a和z中的非零元素保留，为零的元素移除；

23、步骤4.12中，只保留a和z中非零的元素进行后续计算，这样可以避免对a中为0的元素进行不必要的搜索，以提高搜索效率。

24、步骤4.13、计算出目标函数g(t)所有的折点，并将所有折点根据如下计算公式分为左折点l和右折点r：

25、

26、

27、随后将l和r合并，记作kink，并将其中的元素按照升序排列；

28、步骤4.13中，由于软阈值算子sλ(·)的存在，目标函数g(t)的导函数g′(t)无法直接表示出来，但是g′(t)是一个单调递增函数，且具有分段线性函数的性质，即斜率m和截距-b在不同的区间内会发生变化，因此若将导函数记作g′(t)＝m*t-b,需要找出每一个间断点，在每两个间断点构成的区间进行搜索，因此计算出左折点和右折点即为所有的间断点，另一方面，在后续的搜索过程中根据左右折点的不同，更新斜率的公式不同。

29、步骤4.14、初始化g(t)的斜率m，初始的m的值为：a中所有满足l>0和r<0的元素的平方的和；

30、步骤4.15、共进行n次搜索，其中n等于折点的个数，在第i次搜索时，先检查t是否跨过了一个折点，判断的依据为：

31、m*kink(i)<b

32、其中kink(i)表示kink中的第i个值；满足上述条件则代表跨过了折点，否则进行下一次搜索；如果t跨过了一个左折点，则通过如下公式更新b和m：

33、b＝b-kink(i)*a2(i)

34、m＝m-a2(i)

35、其中a2(i)表示数组a的平方后的第i个值；若步长跨过了一个右折点，则通过如下公式更新b和m：

36、b＝b+kink(i)*a2(i)

37、m＝m+a2(i)

38、完成n次搜索后即可确定最终的步长tk，如果更新后m的值为0，则步长为：

39、tk＝s*kink(n)

40、即把数组kink中最大的折点赋值给tk；如果更新后m的值不为0，则步长为：

41、

42、当t穿过拐点时，需要更新m和b来确定新的区间，这样才能遍历不同的区间，并找到最优解。遍历所有的区间后最终的精确步长为：若m＝0，则令t＝s*kink(n)为近似最优步长。

43、进一步的，所述终止条件为：前后两次迭代对应的解向量的改变量充分小或者达到最大迭代次数。

44、本发明采用以上技术方案后，与现有技术相比，具有以下优点：

45、本发明的方法适用于任何已知测量矩阵和稀疏基矩阵的已压缩的测量图像的图像重建，本发明使用基于精确步长的贪婪随机稀疏kaczmarz图像重建方法进行图像重建，在贪婪随机稀疏kaczmarz方法的基础上使用线性搜索的方法得到每一次迭代的精确步长，减少了迭代所需次数和计算量，提高了重建效率和重建精度。

46、下面结合附图和实例对本发明进行详细说明。