本文關(guān)鍵詞：稀疏，由筆耕文化傳播整理發(fā)布。

聲明：本人屬于絕對的新手，剛剛接觸“稀疏表示”這個領(lǐng)域。之所以寫下以下的若干個連載，是鼓勵自己不要急功近利，而要步步為贏！所以下文肯定有所紕漏，敬請指出，我們共同進(jìn)步！

踏入“稀疏表達(dá)”（Sparse Representation）這個領(lǐng)域，純屬偶然中的必然。之前一直在研究壓縮感知（Compressed Sensing）中的重構(gòu)問題。照常理來講，首先會找一維的稀疏信號（如下圖）來驗證CS理論中的一些原理，性質(zhì)和算法，如測量矩陣為高斯隨機矩陣，貝努利矩陣，亞高斯矩陣時使用BP,MP,OMP等重構(gòu)算法的異同和效果。然后會找來二維稀疏信號來驗證一些問題。當(dāng)然，就像你所想的，這些都太簡單。是的，接下來你肯定會考慮對于二維的稠密信號呢，如一幅lena圖像？我們知道CS理論之所以能突破乃奎斯特采樣定律，使用更少的采樣信號來精確的還原原始信號，其中一個重要的先驗知識就是該信號的稀疏性，不管是本身稀疏，還是在變換域稀疏的。因此我們需要對二維的稠密信號稀疏化之后才能使用CS的理論完成重構(gòu)。問題來了，對于lena圖像這樣一個二維的信號，其怎樣稀疏表示，在哪個變換域上是稀疏的，稀疏后又是什么？于是竭盡全力的google…后來發(fā)現(xiàn)了馬毅的“Image Super-Resolution via Sparse Representation”（IEEE Transactions on Image Processing，Nov.2010）這篇文章，于是與稀疏表達(dá)的緣分開始啦！ [break] 談到稀疏表示就不能不提下面兩位的團(tuán)隊，Yi Ma ANDElad Michael，國內(nèi)很多高校（像TSinghua，USTC）的學(xué)生直奔兩位而去。（下圖是Elad M的團(tuán)隊，后來知道了CS界大牛Donoho是Elad M的老師，怪不得…）其實對于馬毅，之前稍有了解，因為韋穗老師，我們實驗室的主任從前兩年開始著手人臉識別這一領(lǐng)域并且取得了不錯的成績，人臉識別這個領(lǐng)域馬毅算是大牛了…因此每次開會遇到相關(guān)的問題，韋老師總會提到馬毅，于是通過各種渠道也了解了一些有關(guān)他科研和個人的信息。至于Elad.M，恕我直言，我在踏入這個領(lǐng)域其實真的完全不知道，只是最近文章看的比較多，發(fā)現(xiàn)看的文章中大部分的作者都有Elad，于是乎，好奇心驅(qū)使我了解了這位大牛以及他的團(tuán)隊成員…也深深的了解到了一個團(tuán)隊對一個領(lǐng)域的貢獻(xiàn)，從Elad.M那兒畢業(yè)的學(xué)生現(xiàn)在都成了這個領(lǐng)域中的佼佼者…不禁感嘆到：一個好的導(dǎo)師是多么的重要��！下面舉個簡單的例子，說說二維信號的稀疏性，也為后面將稀疏表示做個鋪墊。我們以一幅大小為256×256的Lena圖像為例，過完備字典（Dictionary，具體含義見后文，先理解為基吧，其實不完全等同）選擇離散余弦變換DCT，字典大小選擇64×256，，對圖像進(jìn)行分塊處理，由于僅僅為了說明稀疏性的概念，所以不進(jìn)行重疊處理，每塊大小8×8（pixel）…簡單的稀疏表示后的稀疏如下圖�？梢钥闯�，絕大多數(shù)的稀疏集中在0附近。當(dāng)然，這里僅僅是簡單的說明一下。后面我們有更好的選擇，比如說，字典的選擇，圖像塊的選擇等等…

 

| JustPlus

-----------------------------------------------------

壓縮感知（CS），或許你最近聽說的比較多，不錯，CS最近比較火，什么問題不管三七二十一就往上粘連，先試試能不能解決遇到的問題，能的話就把文章發(fā)出來忽悠大家，這就是中國學(xué)術(shù)浮躁的表現(xiàn)…我們沒有時間去思考的更多，因為你一思考，別人可能就把“你的東西”搶先發(fā)表了…不扯了，反正也干預(yù)不了…稀疏表示的現(xiàn)狀有點像CS，能做很多事，也不能做很多事…但是它確實是解決一些棘手問題的方法，至少能提供一種思路…目前用稀疏表示解決的問題主要集中在圖像去噪（Denoise），代表性paper：Image Denoise Via Sparse and Redundant Representations Over Learned Dictionaries(Elad M. and Aharon M. IEEE Trans. on Image Processing,Dec,2006);Image Sequence Denoising Via Sparse and Redundant Representations(Protter M. and Elad M.IEEE Trans. on Image Processing,Jan,2009)， 還有超分辨率（Super-Resolution OR Scale-Up），代表性paper：Image Super-Resolution via Sparse Representation（Jianchao Yang, John Wright, Thomas Huang, and Yi Ma，IEEE Transactions on Image Processing, Nov,2010）,A Shrinkage Learning Approach for Single Image Super-Resolution with Overcomplete Representations( A. Adler, Y. Hel-Or, and M. Elad,ECCV,Sep,2010)…. 另外還有inpait，deblur，F(xiàn)ace Recognition，compression等等..更多應(yīng)用參考Elad M的書，google能找到電子檔，這里不提
供下載地址

當(dāng)然Elad M.和Yi Ma的團(tuán)隊不僅僅在應(yīng)用上大做文章，在理論上也是不斷革新…就拿字典學(xué)習(xí)為例，一開始將固定字典（如過完備 DCT，Contourlet，Wavelet字典）用在去噪，超分辨率上，效果不明顯，可能某些情況下還不如空域或者頻域的去噪效果好（我拿Lena圖像和DCT實驗了，以PSNR為標(biāo)準(zhǔn)，比時域的去噪方法要好，但是比小波去噪相比，稍稍遜色），但是速度快是它的優(yōu)點。于是他們開始研究自適應(yīng)的字典學(xué)習(xí)算法，開始使用很多圖像進(jìn)行學(xué)習(xí)，后來采用單幅圖像進(jìn)行學(xué)習(xí)來提高運算速度（使用很多圖像進(jìn)行學(xué)習(xí)屬于半自適應(yīng)的學(xué)習(xí)，對于自然圖像的處理需要學(xué)習(xí)自然圖像，對遙感圖像的處理需要學(xué)習(xí)遙感圖像，但是對自然圖像或遙感圖像的去噪，超分辨率處理，都可以使用已經(jīng)訓(xùn)練好的相應(yīng)的字典）；同時學(xué)習(xí)的方法也不盡相同，開始使用MOD，后來就是一直比較流行的K-SVD，最近又出來了Online，總體而言O(shè)nline比較快。下面是我提到的幾種字典的例子，所有的字典都是64×256大小的，依次為DCT，globally（訓(xùn)練圖像是：標(biāo)準(zhǔn)圖像 lena，boat，house，barbara，perppers），K-SVD（單幅含躁lena，噪聲標(biāo)準(zhǔn)差為25），online（單幅含躁 lena，噪聲標(biāo)準(zhǔn)差為25），其中g(shù)lobally的訓(xùn)練方法是將訓(xùn)練圖像分成8×8的overlap patch，平均取，共取10000塊，K-SVD和online也是分成相同的重疊塊，取所有可能的塊。

       

總之，Elad M.和Yi Ma為稀疏表示這個領(lǐng)域作出了很大的貢獻(xiàn)…向大牛們致敬��！最后稍微說一下國內(nèi)的研究現(xiàn)狀，國內(nèi)的很多研究還沒浮出水面，不知道是不是我想的的這樣（我比較疑惑的是，為什么國外研究了十幾年，國內(nèi)還沒大動靜？），至少從google學(xué)術(shù)以及IEEE的文章搜索上來看是這樣的…不過還是有幾位教授在這方面作出了很大的貢獻(xiàn)的…

-------------------------------------------------

我們來考慮信號的稀疏表示問題，假如我們有了過完備字典D，如何求出信號x在這個過完備字典上的稀疏表示？先來回顧一下在壓縮感知中常常會遇到的問題，信號x在經(jīng)過測量矩陣A后得到測量值y，即y=A*x，其中測量矩陣A_mxn(m遠(yuǎn)小于n)，那么怎樣從y中精確的恢復(fù)出x呢？

        

由于m遠(yuǎn)小于n，用m個方程求解n個未知數(shù)，因此y=A*x是個欠定方程，有無窮多個解。就像我們解優(yōu)化問題一樣，如果我們加上適當(dāng)?shù)南薅l件，或者叫正則項，問題的解會變得明朗一些！這里我們加上的正則項是norm(x,0)，即使重構(gòu)出的信號 x盡可能的稀疏（零范數(shù)：值為非0的元素個數(shù)），后來Donoho和Elad這對師徒證明了如果A滿足某些條件，那么argmin norm(x,0) s.t.y=A*x 這個優(yōu)化問題即有唯一解！唯一性確定了，仍然不能求解出該問題，后來就嘗試使用l1和l2范數(shù)來替代l0范數(shù)，華裔科學(xué)家陶哲軒和candes合作證明了在A滿足UUP原則這樣一個條件下，l0范數(shù)可以使用l1范數(shù)替代，所以優(yōu)化問題變成argmin norm(x,1) s.t.y=A*x這樣一個凸優(yōu)化問題，可以通過線性優(yōu)化的問題來解決�。▍⒖嘉墨I(xiàn)：Stable signal recovery from incomplete and inaccurate measurements（E. J. Candès, J. Romberg and T. Tao））至此，稀疏表示的理論已經(jīng)初步成型。至于之后的優(yōu)化問題，都是一些變形，像lasso模型，TV模型等….這里推薦一本stanford的凸優(yōu)化教材convex optimization，我準(zhǔn)備抽個時間好好看一看，搞稀疏表達(dá)這一塊的，優(yōu)化問題少不了…最近一直在感嘆：數(shù)學(xué)不好的人哪傷不起�。�！有木有��！[break]

        

我們再回到我們一開始提到的問題上來，一開始說字典D已知，求出y在過完備字典D上的稀疏表示x，這個在稀疏表示里面被稱作為Sparse Coding…問題的模型是 x=argmin norm(y-D*x,2)^2 s.t.norm(x,1)<=k; 我們下面來介紹一下解決這個問題的常用方法OMP(Orthogonal Matching Pursuit) .我們主要目標(biāo)是找出x中最主要的K個分量（即x滿足K稀疏），不妨從第1個系數(shù)找起，假設(shè)x中僅有一個非零元x(m)，那么 y0=D(:,m)*x(m)即是在只有一個主元的情況下最接近y的情況，norm(y-y0,2)/norm(y,2)<=sigma,換句話說在只有一個非零元的情況下，D的第m列與y最“匹配”，要確定m的值，只要從D的所有列與y的內(nèi)積中找到最大值所對應(yīng)的D的列數(shù)即可，然后通過最小二乘法即可確定此時的稀疏系數(shù)�？紤]非零元大于1的情況，其實是類似的，只要將余量r=y-y0與D的所有列做內(nèi)積，找到最大值所對應(yīng)D的列即可�！旅媸谴a和示例結(jié)果[其中圖1是lena圖像某個8x8的圖像塊在其自身訓(xùn)練得到的字典上的稀疏表示，稀疏值k=30，相對誤差norm(y- y0,2)/norm(y,2)=1.8724，圖2是相同的塊在相同的字典下的稀疏表述，稀疏值k=50，相對誤差為0.0051]

function A=OMP(D,X,L) % 輸入?yún)?shù): % D - 過完備字典，注意：必須字典的各列必須經(jīng)過了規(guī)范化 % X - 信號 %L - 系數(shù)中非零元個數(shù)的最大值（可選，默認(rèn)為D的列數(shù)，速度可能慢） % 輸出參數(shù): % A - 稀疏系數(shù) if nargin==2 L=size(D,2); end P=size(X,2); K=size(D,2); for k=1:1:P, a=[]; x=X(:,k); residual=x; indx=zeros(L,1); for j=1:1:L, proj=D'*residual; [maxVal,pos]=max(abs(proj)); pos=pos(1); indx(j)=pos; a=pinv(D(:,indx(1:j)))*x; residual=x-D(:,indx(1:j))*a; if sum(residual.^2) < 1e-6 break; end end; temp=zeros(K,1); temp(indx(1:j))=a; A(:,k)=sparse(temp); end; return;

 

----------------------------------------------------

回顧一下前面所說的OMP算法，前提條件是字典D已知，求一個信號在這個字典上的稀疏表示…那假如我們不想使用過完備 DCT，wavelet呢？因為它們沒有自適應(yīng)的能力，不能隨著信號的變化作出相應(yīng)的變化，一經(jīng)選定，對所有的信號一視同仁，當(dāng)然不是我們想要的。我們想通過訓(xùn)練學(xué)習(xí)的方法獲取字典D來根據(jù)輸入信號的不同作出自適應(yīng)的變化。那現(xiàn)在我們的未知量有兩個，過完備字典D，稀疏系數(shù)x，已知量是輸入信號y，當(dāng)然先驗知識是輸入信號在字典D上可以稀疏表示…我們再次列出sparse-land模型：[D,x]=argmin norm(y-D*x,2)^2 s.t.norm(x,1)<=k。如何同時獲取字典D和稀疏系數(shù)x呢？方法是將該模型分解：第一步將D固定，求出x的值，這就是你常聽到的稀疏分解（Sparse Coding），也就是上一節(jié)提到的字典D固定，求信號y在D上稀疏表示的問題；第二步是使用上一步得到的x來更新字典D，即字典更新（Dictionary Update）。如此反復(fù)迭代幾次即可得到優(yōu)化的D和x。

Sparse Coding：x=argmin norm(y-D*x,2) s.t.norm(x,1)<=k
Dictinary Update:D=argmin norm(y-D*x,2)^2

[break]我們主要通過實例介紹三種方法：MOD，K-SVD,Online… 首先是MOD(Method of Optimal Direction)。Sparse Coding其采用的方法是OMP貪婪算法，Dictionary Update采用的是最小二乘法，即D=argmin norm(y-D*x,2)^2 解的形式是D=Y*x'*inv(x*x’)。因此MOD算法的流程如下：

初始化： 字典D_mxn可以初始化為隨機分布的mxn的矩陣，也可以從輸入信號中隨機的選取n個列向量，下面的實驗我們選取后者。注意OMP要求字典的各列必須規(guī)范化，因此這一步我們要將字典規(guī)范化。根據(jù)輸入信號確定原子atoms的個數(shù)，即字典的列數(shù)。還有迭代次數(shù)。
主循環(huán)： Sparse Coding使用OMP算法； Dictionary Update采用最小二乘法。注意這一步得到的字典D可能會有列向量的二范數(shù)接近于0，此時為了下一次迭代應(yīng)該忽略該列原子，重新選取一個服從隨機分布的原子。

01

02

03

04

05

06

07

08

09

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

function [A,x]= MOD(y,codebook_size,errGoal)

%==============================

%input parameter

%   y - input signal

%   codebook_size - count of atoms

%output parameter

%   A - dictionary

%   x - coefficent

%==============================

if(size(y,2)

  本文關(guān)鍵詞：稀疏，由筆耕文化傳播整理發(fā)布。

本文編號：85885