豐 祥， 萬(wàn)旺根

1.上海大學(xué)通信與信息工程學(xué)院，上海200444

2.上海大學(xué)智慧城市研究院，上海200444

信號(hào)處理的依據(jù)是奈奎斯特采樣定理，即如果要精確地重構(gòu)原始信號(hào)，必須使用大于信號(hào)最高頻率兩倍的頻率進(jìn)行信號(hào)采樣，這就導(dǎo)致信號(hào)采樣過(guò)程中產(chǎn)生大量冗余數(shù)據(jù)，而最終只有部分重要的信息被應(yīng)用.文獻(xiàn)[1]在可壓縮信號(hào)的稀疏分解與信號(hào)恢復(fù)基礎(chǔ)上提出了壓縮感知理論框架，隨后該理論迅速發(fā)展，為解決瓶頸問(wèn)題提供了理論基礎(chǔ).在此基礎(chǔ)上，Candes等人正式提出了壓縮感知這一術(shù)語(yǔ).此后，人們?cè)谙∈璞硎镜南∈栊訹2]與不相干性、壓縮采樣的穩(wěn)健性、感知測(cè)量等[3]方面開展了許多研究工作，并將壓縮感知作為測(cè)量技術(shù)應(yīng)用于天文學(xué)、核磁共振、模式識(shí)別等領(lǐng)域，取得了較好的研究成果[4].

壓縮感知理論是對(duì)經(jīng)典信號(hào)處理領(lǐng)域的補(bǔ)充和完善，引起了眾多學(xué)者和組織機(jī)構(gòu)的關(guān)注，其應(yīng)用領(lǐng)域包含圖像處理與重建、數(shù)據(jù)融合和醫(yī)學(xué)信號(hào)檢測(cè)等.文獻(xiàn)[5]研究了壓縮感知稀疏表示算法；文獻(xiàn)[6-7]詳細(xì)探討了壓縮感知理論中約束等距性（restricted isometry property,RIP）矩陣的構(gòu)造問(wèn)題；文獻(xiàn)[8]研究了壓縮感知理論在探地雷達(dá)三維成像中的應(yīng)用；文獻(xiàn)[9]探索了感知測(cè)量矩陣在實(shí)際應(yīng)用中的設(shè)計(jì)問(wèn)題；文獻(xiàn)[10-11]將壓縮感知理論用于解決醫(yī)學(xué)圖像重建和圖像去噪等問(wèn)題.

文獻(xiàn)[12]提出了基于多觀測(cè)向量和稀疏貝葉斯學(xué)習(xí)的重建算法，并將壓縮感知理論用于圖像重建，但沒有考慮圖像分塊和稀疏表示方法對(duì)圖像重建的影響.文獻(xiàn)[13]基于多種稀疏變換基和觀測(cè)矩陣研究圖像重建，但沒有考慮圖像分塊和采樣率對(duì)圖像稀疏重建的影響.本文分析比較了離散余弦變換和離散小波變換兩種稀疏表示算法，通過(guò)調(diào)節(jié)分塊大小和采樣率大小對(duì)標(biāo)準(zhǔn)測(cè)試圖像進(jìn)行重建實(shí)驗(yàn)，并在運(yùn)行時(shí)間、重構(gòu)誤差和視覺效果方面綜合比較兩種算法對(duì)稀疏重構(gòu)效果和效率的影響，指出了各自算法的優(yōu)缺點(diǎn)，分析了分塊大小和采樣率大小對(duì)圖像重建效果的影響.

1 壓縮感知的研究?jī)?nèi)容

壓縮感知理論采用低于奈奎斯特采樣定理要求的頻率采樣、編碼和重構(gòu)，適用于更廣泛的信號(hào)類型.該理論主要涉及稀疏變換、觀測(cè)測(cè)量和恢復(fù)重建3個(gè)基本核心問(wèn)題[14]，其流程圖見圖1.壓縮感知要求原始信號(hào)為可壓縮信號(hào)，即首先可以通過(guò)稀疏變換把原始信號(hào)轉(zhuǎn)化為稀疏信號(hào)，然后使用觀測(cè)矩陣對(duì)稀疏變換進(jìn)行觀測(cè)測(cè)量，得到低維的觀測(cè)信號(hào)，經(jīng)傳輸和存儲(chǔ)等過(guò)程后，根據(jù)恢復(fù)重建算法將觀測(cè)信號(hào)重構(gòu)出原始信號(hào).

稀疏變換是在假設(shè)自然信號(hào)可以被壓縮表示時(shí)對(duì)多維數(shù)據(jù)進(jìn)行線性分解的一種表示方法，它具有過(guò)完備性和稀疏性兩個(gè)基本特征.過(guò)完備性表示字典的行數(shù)大于列數(shù)，稀疏性表示信號(hào)非零元數(shù)目足夠少.本文針對(duì)離散余弦變換和離散小波變換這兩種稀疏變換進(jìn)行研究.觀測(cè)測(cè)量指的是：對(duì)信號(hào)的線性投影采用一個(gè)與稀疏變換矩陣不相關(guān)的測(cè)量矩陣得到感知測(cè)量值，以確保精確地重構(gòu)原始信號(hào).測(cè)量矩陣應(yīng)滿足非相干性和限制等容性兩個(gè)基本原則.常用的隨機(jī)性測(cè)量矩陣包括隨機(jī)高斯測(cè)量矩陣、隨機(jī)伯努利矩陣、部分正交矩陣、稀疏隨機(jī)矩陣等，本文則采用正交歸一化的隨機(jī)高斯矩陣進(jìn)行觀測(cè)測(cè)量.信號(hào)重構(gòu)算法是指用壓縮測(cè)量的低維數(shù)據(jù)精確地重構(gòu)高維的原始信號(hào)或圖像，即利用M維測(cè)量值重建N維信號(hào)（M遠(yuǎn)小于N）的過(guò)程.信號(hào)重構(gòu)是壓縮感知研究中最重要而關(guān)鍵的部分.目前主流的重構(gòu)算法分為3類：基于最小l范數(shù)的凸優(yōu)化算法、貪婪算法、非凸優(yōu)化算法，其中貪婪算法在保證重建精度的同時(shí)具有較快的重建速度，于是本文采用貪婪算法中的正交匹配追蹤算法進(jìn)行恢復(fù)重建.

圖1 壓縮感知理論流程圖Figure 1 Flowchart of compressed sensing theory

2 圖像的稀疏表示與恢復(fù)重建

為了對(duì)圖像信號(hào)進(jìn)行稀疏表示，應(yīng)先把信號(hào)變換到稀疏域.本文在圖像稀疏表示階段采用離散余弦變換和離散傅里葉變換進(jìn)行對(duì)比性研究，在圖像恢復(fù)重建階段采用正交匹配追蹤算法進(jìn)行重建.

2.1 圖像的稀疏表示

2.1.1 離散小波變換

若ψ(t)∈ L2(R)，（L2(R)為ψ(t)的矢量 空間，R為實(shí)數(shù)集），a0＞1為擴(kuò)張步長(zhǎng)，ψj,k(t)=則定義一維信號(hào)x(t)的離散小波變換（DWT變換）為

當(dāng)a0=2且τ0=1時(shí)，所有的采樣點(diǎn)所形成的圖形被稱為二進(jìn)網(wǎng)格，此類離散小波變換被稱為二進(jìn)小波變換.它只選擇部分縮放因子和平移參數(shù)進(jìn)行離散小波變換，大大減少了分析的數(shù)據(jù)量.

假設(shè)一幅圖像大小為N×N，該圖像可以用長(zhǎng)度為N的二維方陣X進(jìn)行表示，w w表示長(zhǎng)度為N的一維離散小波變換所對(duì)應(yīng)的正交規(guī)范化變換矩陣，則矩陣X的二維離散小波變換為

該變換可以分解為以下兩步來(lái)執(zhí)行：1）對(duì)矩陣X的每一行數(shù)據(jù)序列進(jìn)行一維離散小波變換Xi=X·w w′；2）對(duì)變換后矩陣Xi的每一列進(jìn)行離散一維小波小波Xo=ww·Xi.Xo往往是稀疏的，它屬于稀疏矩陣；X屬于可壓縮信號(hào).在實(shí)際過(guò)程中，為了提高稀疏重構(gòu)效率，通常將整體圖像分成若干個(gè)大小為Nb×Nb的小塊進(jìn)行處理，每一分塊在水平方向和垂直方向的索引號(hào)分別為x和y，對(duì)應(yīng)的正交小波變換為=w w·Xx,y·w w′.記M為測(cè)量數(shù)據(jù)量大小，η=M/N0.測(cè)量矩陣取為行正交規(guī)范化的隨機(jī)高斯矩陣R，各個(gè)分塊的測(cè)量值為Y=R·Xx,yo，采用正交匹配追蹤法求得重構(gòu)值^Xx,yo，經(jīng)過(guò)小波逆變換計(jì)算出圖像域的數(shù)據(jù)值=ww··w w′，最后根據(jù)分塊的順序合成整體圖像并進(jìn)行顯示.

圖像經(jīng)過(guò)DWT變換后得到低頻部分系數(shù)和高頻部分系數(shù)，低頻部分系數(shù)通常具有非稀疏性，而高頻部分系數(shù)具有良好的稀疏性.如果低頻部分系數(shù)與高頻部分系數(shù)都作為原始信號(hào)進(jìn)行感知測(cè)量，往往會(huì)導(dǎo)致重構(gòu)圖像質(zhì)量嚴(yán)重下降[15].尺度越小分辨率越高，則高頻成分越多，小波分解的尺度對(duì)重構(gòu)結(jié)果具有重要的影響.

2.1.2 離散余弦變換

離散余弦變換將數(shù)據(jù)線性地變換到頻域，并用一系列系數(shù)來(lái)表示數(shù)據(jù)，其優(yōu)點(diǎn)是可以把原始數(shù)據(jù)的能量集中在少數(shù)低頻成分，因?yàn)樵紨?shù)據(jù)的相關(guān)性越強(qiáng)，能量越集中.離散余弦變換（discrete Cosine transform,DCT）常用于圖像數(shù)據(jù)的壓縮，經(jīng)過(guò)DCT變換后的數(shù)據(jù)能量非常集中，主要集中于離散余弦變換后的直流部分和低頻部分.

對(duì)于N×N的二維圖像信號(hào)f(x,y),0≤x,y≤N-1，其二維離散余弦變換的正變換公式為

2018年9月10日習(xí)總書記在全國(guó)教育大會(huì)上強(qiáng)調(diào)：“教師是人類靈魂的工程師，是人類文明的傳承者，承載著傳播知識(shí)、傳播思想、傳播真理，塑造靈魂、塑造生命、塑造新人的時(shí)代重任”。由此可見，高校輔導(dǎo)員肩負(fù)著立德樹人的神圣使命，只有謹(jǐn)記教育之初心，率先垂范，為人師表，立德、立言、立行，不斷提高自身學(xué)識(shí)修養(yǎng)，才能做好學(xué)生健康成長(zhǎng)的人生導(dǎo)師。筆者認(rèn)為，高校輔導(dǎo)員具體可以從如下四個(gè)方面著手努力：

相應(yīng)的二維DCT變換的反變換公式為

式（3）和（4）中的系數(shù)為

2.2 信號(hào)的恢復(fù)重建

正交匹配追蹤（orthogonal matching pursuit,OMP）重構(gòu)算法的原理如下[16]：以貪婪迭代的方法選擇恢復(fù)矩陣R的列，使得每次迭代中所選擇的列與當(dāng)前的殘差向量最大程度地相關(guān)，再?gòu)臏y(cè)量向量中減去相關(guān)部分，反復(fù)迭代直到迭代次數(shù)達(dá)到稀疏度K時(shí)停止迭代.該算法的輸入?yún)?shù)為恢復(fù)矩陣R、測(cè)量向量y、稀疏度K，輸出參數(shù)為原始信號(hào)x的K稀疏逼近^x.算法的主要步驟如下：

步驟1 初始化殘差r0=y，索引集Λ0=?，迭代次數(shù)t=1，循環(huán)執(zhí)行步驟2～5.

步驟2 找出殘差r和感知矩陣的列φj積中最大值所對(duì)應(yīng)的角標(biāo)λt，即λt=argmaxj=1,···,N|〈rt-1,φj〉|，記錄最大值所對(duì)應(yīng)的列向量.

步驟3 更新索引集Λt=Λt-1∪{λt}，記錄找到的感知矩陣中的重構(gòu)原子集合Φt=[Φt-1,φλt].

步驟5 判斷迭代次數(shù)是否滿足t＞K，若滿足，則停止迭代并跳轉(zhuǎn)到步驟6；若不滿足，則執(zhí)行步驟2.

步驟6 將最終的^x值作為稀疏重構(gòu)向量.

3 實(shí)驗(yàn)仿真與分析

3.1 實(shí)驗(yàn)說(shuō)明

本文采用大小為N×N(N=512)的灰度圖像Lena、Cameraman和Peppers作為標(biāo)準(zhǔn)測(cè)試圖像，圖像分塊時(shí)的分塊邊長(zhǎng)Nb取值為32、64、128，采樣率η取值為0.3、0.4、0.5.本文將比較DCT變換和DWT變換在圖像稀疏表示時(shí)的性能，采用正交歸一化的隨機(jī)高斯測(cè)量矩陣進(jìn)行感知測(cè)量，根據(jù)正交匹配追蹤算法進(jìn)行圖像信號(hào)的恢復(fù)重建，以時(shí)間消耗t和重建誤差σ表征不同分塊大小Nb和采樣率η條件下的重建效率和性能.其中，時(shí)間消耗包含信號(hào)分解、感知測(cè)量和信號(hào)重建整個(gè)過(guò)程中所花的時(shí)間.令Xerror為整體重建圖像與整體原始圖像X之間各個(gè)像素差值的平方和，則重建誤差為

所有實(shí)驗(yàn)均在MATLAB R2011b環(huán)境下完成，計(jì)算機(jī)的CPU性能參數(shù)為i5-2400 3.10 GHz，內(nèi)存大小為2.00 GB.

3.2 實(shí)驗(yàn)仿真比較與數(shù)據(jù)分析

當(dāng)改變采樣率和分塊大小時(shí)，記錄時(shí)間消耗和重建誤差，采用DCT變換和DWT變換得到的實(shí)驗(yàn)數(shù)據(jù)如表1和2所示.由表1和2可以發(fā)現(xiàn)，無(wú)論是采用DCT變換還是采用DWT變換，都滿足以下結(jié)論：在分塊大小固定的條件下，時(shí)間消耗隨著采樣率的提高而增加，而重建誤差隨著采樣率的提高而降低；在采樣率固定的條件下，運(yùn)行時(shí)間隨著分塊大小的增加而增加，而重建誤差隨著分塊大小的增加而降低.本文以Lena圖像和DCT變換為例，在不同采樣率和分塊大小時(shí)分別比較重建誤差和時(shí)間消耗，結(jié)果如圖2所示.在實(shí)際應(yīng)用中，往往難以同時(shí)獲得最短的時(shí)間消耗和最小的重建誤差，這就需要根據(jù)實(shí)際需求合理地選擇分塊大小和采樣率，以實(shí)現(xiàn)最佳的重建效果.

表1 采用DCT變換時(shí)的重建性能Table 1 Reconstruction performance using DCT transform

表2 采用DWT變換時(shí)的重建性能Table 2 Reconstruction performance using DWT transform

對(duì)比表1和2的重建誤差和時(shí)間消耗數(shù)據(jù)容易發(fā)現(xiàn)：在相同的分塊大小和采樣率條件下，對(duì)相同的標(biāo)準(zhǔn)圖像進(jìn)行測(cè)試，采用DCT變換時(shí)重建誤差和時(shí)間消耗都比采用DWT變換時(shí)小，采用CT稀疏表示在圖像重建方面往往表現(xiàn)出更優(yōu)的性能.圖3為原始標(biāo)準(zhǔn)的測(cè)試圖像，圖4表示在Nb=128,η=0.5時(shí)分別采用DCT變換和DWT變換時(shí)所對(duì)應(yīng)的測(cè)試圖像重建效果圖.對(duì)比圖3和4可以發(fā)現(xiàn)，采用DWT變換時(shí)的重建圖像往往大部分區(qū)域較清晰，而在小部分區(qū)域會(huì)出現(xiàn)明顯的局部帶狀條紋；采用DCT變換重建圖像的整體清晰度比采用DCT變換時(shí)重建圖像的整體清晰度較差，但不會(huì)出現(xiàn)明顯的帶狀條紋.

圖2 以Lena圖像和DCT變換為例，重建誤差和時(shí)間消耗比較圖Figure 2 Comparison of reconstruction error and time consuming with Lena image and DCT for example

分析表1、表2、圖3、圖4可知，采樣率的適當(dāng)選取與測(cè)量信號(hào)的稀疏程度密切相關(guān)，而對(duì)于非稀疏信號(hào)，較低的采樣率則會(huì)嚴(yán)重影響重建圖像的質(zhì)量.相比于DWT變換，DCT變換因不受小波分解尺度的影響而在圖像的稀疏化表達(dá)方面表現(xiàn)得更加穩(wěn)定.本文采用sym8作為小波基，保持小波分解尺度不變，并且可以調(diào)整分塊大小，此時(shí)分塊大小越小，重構(gòu)圖像質(zhì)量則越差.在采樣率為0.5、分塊大小為128的情況下，小波分解尺度對(duì)稀疏度的影響相對(duì)較小，采用DWT變換時(shí)的重建誤差與采用DCT時(shí)的重建誤差相比差異也較??；而當(dāng)分塊大小減小或采樣率下降時(shí)，采用DWT變換時(shí)的重建效果則明顯降低.

圖3 原始標(biāo)準(zhǔn)測(cè)試圖像Figure 3 Original standard test images

圖4 分別為使用DCT變換和DWT變換時(shí)的稀疏重建圖像Figure 4 Reconstructed images using DCT transform and DWT transform respectively

4 結(jié)論

本文將壓縮感知理論應(yīng)用于圖像的稀疏重建，在圖像的感知測(cè)量和恢復(fù)重建過(guò)程中分別采用了正交歸一化的隨機(jī)高斯矩陣和正交匹配追蹤算法，在圖像的稀疏表示過(guò)程中采用DCT變換和DWT變換，并在不同的分塊大小和采樣率條件下，根據(jù)時(shí)間消耗和重建誤差兩個(gè)性能參數(shù)對(duì)各自的性能進(jìn)行了綜合比較.在分塊大小固定的條件下，時(shí)間消耗隨著采樣率的提高而增加，而重建誤差隨著采樣率的提高而降低；在采樣率固定的條件下，運(yùn)行時(shí)間隨著分塊大小的增加而增加，而重建誤差隨著分塊大小的增加而降低，在圖像的稀疏表示與重建方面，DCT變換整體上比DWT變換表現(xiàn)出相對(duì)更優(yōu)的性能.在實(shí)際應(yīng)用中需根據(jù)對(duì)重建圖像的要求，結(jié)合分塊大小和采樣率等關(guān)鍵因素對(duì)重建性能的影響，合理地選擇分塊大小和采樣率，以實(shí)現(xiàn)整體最優(yōu)的重建效果.下一步我們將針對(duì)不同的圖像進(jìn)行自適應(yīng)分塊，希望在時(shí)間消耗和重建質(zhì)量等方面取得更好的效果.

[1]CANDES E J,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.

[2]LIU J,MALLICK M,HAN C Z,YAO X H,LIAN F.Similar sensing matrix pursuit:an efficient reconstruction algorithm to cope with deterministic sensing matrix[J].Signal Processing,2014,95:101-110.

[3]LIU Xiaoman,ZHU Yonggui.A fast method for TVL1-MRI image reconstruction in compressive sensing[J].Journal of Computational Information Systems,2014,10(2):691-699.

[4]MICCHELLI C A,SHEN L X,XU Y S,ZENG X Y.Proximity algorithms for the L1/TV image denoising model[J].Advances in Computational Mathematics,2013,38(2):401-426.

[5]XU Zhiqiang.Deterministic sampling of sparse trigonometric polynomials[J].Journal of Complexity,2011,27(2):133-140.

[6]BOURGAIN J,DILw ORTH S J,FORD K,KONYAGIN S,KUTZAROVAD.Explicit constructions of RIP matrices and related problems[J].Duke Mathematical Journal,2011,159(1):145-185.

[7]ZHANG Tong.Sparse recovery with orthogonal matching pursuit under RIP[J].IEEE Transactions on Information Theory,2011,57(9):6215-6221.

[8]YU Huimin,FANG Guangyou.The application of compressive sensing in the three dimensional imaging of ground penetrating radar[J].Journal of Electronics and Information Technology,2010,32(1):12-16.

[9]FU Qiang,LI Qiong.The research of constructing the measurement matrix in compressive sensing[J].Computer and Telecommunication,2011,7(1):39-41.

[10]WANG Yan,LIAN Qiusheng,LI Kai.MRI image reconstruction based on joint regularization and compressed sensing[J].Optical Technology,2010,36(3):350-355.

[11]QU Xiaobo,GUO Di,NING Bende.Undersampled MRI reconstruction with the patch-based directional wavelets[J].Magnetic Resonance Imaging,2012,30(7):964-977.

[12]LI Zhilin,CHEN Houjin,LI Jupeng,YAO Chang,YANG Na.An efficient algorithm for compressed sensing image reconstruction[J].Acta Electronica Sinica,2011,39(12):2796-2800.

[13]HOU Jinman,HE Ning,LüKe.Image fast reconstruction method based on compressive sensing[J].Computer Engineering,2011,37(19),215-217.

[14]YINHongpeng,LIUZhaodong,CHAIYi,JIAOXuguo.Survey of compressed sensing[J].Control and Decision,2013,28(10):1441-1445.

[15]YUAN Quan,ZHANG Cheng,CHEN Jianjun,YAO Junxia,LIYingran,WANGHuan.Image compressed sensing based on wavelet transform[J].Computer Engineering,2012,38(20):209-218.

[16]CARMIA Y,MIHAYLOVAL S,GODSILLSJ.Introduction to compressed sensing and sparse f iltering[M].Berlin Heidelberg:Springer,2014:1-23.