姓名：方文 19021210911 ? ?
轉(zhuǎn)載自?https://blog.csdn.net/jbb0523/article/details/45099655
【嵌牛導(dǎo)讀】壓縮感知算法中匹配追蹤是一個(gè)基礎(chǔ)概念跨嘉，可以作為進(jìn)一步理解其他重構(gòu)算法的基礎(chǔ)
【嵌牛鼻子】稀疏表示 匹配追蹤
【嵌牛提問(wèn)】 通過(guò)信號(hào)的稀疏表示理解匹配追蹤
【嵌牛正文】
如果研究壓縮感知重構(gòu)算法，那么匹配追蹤絕對(duì)是一個(gè)繞不過(guò)去的概念纠屋，雖然最原始的匹配追蹤算法應(yīng)用于重構(gòu)的文獻(xiàn)并不多，但基于匹配追蹤改進(jìn)后的正交匹配追蹤算法在研究壓縮感知的人群中可以說(shuō)是家喻戶(hù)曉变抽、人盡皆知的亚皂。既然匹配追蹤是與重構(gòu)有關(guān)的，為什么本篇取名為《稀疏表示與匹配追蹤》呢廊移？一個(gè)重構(gòu)算法與稀疏表示有什么關(guān)系扭粱？一般文獻(xiàn)中講到匹配追蹤所引參考文獻(xiàn)均為文獻(xiàn)[1]舵鳞，而文獻(xiàn)[1]就是在講解基于匹配追蹤算法的信號(hào)稀疏表示問(wèn)題，這到底是怎么回事呢琢蛤？以下思路參考或直接引用了文獻(xiàn)[2]中的內(nèi)容蜓堕。
正交基稀疏表達(dá)一個(gè)信號(hào)也有很多缺點(diǎn)，因?yàn)橐唤M基表達(dá)信號(hào)的能力取決于信號(hào)的特征是否與基向量的特征相吻合博其；例如套才，光滑連續(xù)信號(hào)可以被傅里葉基稀疏的表達(dá)，但脈沖信號(hào)就不行慕淡；再如背伴，帶有孤立不連續(xù)點(diǎn)的平滑信號(hào)可以被小波基稀疏表達(dá)，但小波基在表達(dá)傅里葉頻譜中有窄帶高頻支撐的信號(hào)時(shí)卻是無(wú)效的。現(xiàn)實(shí)世界中的信號(hào)經(jīng)常包含有用單一基所不能表達(dá)的特征傻寂。對(duì)于這些信號(hào)息尺，你或許希望可以選擇來(lái)自不同基的向量（如用小波基和傅里葉基來(lái)聯(lián)合表達(dá)一個(gè)信號(hào)）。因?yàn)槟阆氡ＷC你可以表達(dá)一個(gè)信號(hào)空間的所有信號(hào)向量疾掰，所以由所有可選向量組成的字典應(yīng)該能夠張成這個(gè)信號(hào)空間搂誉。然而由于這組字典中的向量來(lái)自不同的基，它們可能不是線(xiàn)性獨(dú)立的静檬。由于這組字典中的向量不是線(xiàn)性獨(dú)立的炭懊，會(huì)造成用這組字典做信號(hào)表達(dá)時(shí)系數(shù)不唯一。然而如果創(chuàng)建一組冗余字典拂檩，你就可以把你的信號(hào)展開(kāi)在一組可以適應(yīng)各種時(shí)頻或時(shí)間-尺度特性的向量上侮腹。這樣構(gòu)造的字典可以極大地增加你稀疏表達(dá)各種特性信號(hào)的能力。
? ? ? ??例如三維空間中有基向量vx=(1,0,0)广恢，vy=(0,1,0)凯旋，vz=(0,0,1)，三維空間中的任意向量均可由三個(gè)基向量唯一表示钉迷，如v1=(1,2,3)= vx+2 vy +3 vz。然而若取三維空間的一組冗余字典va=(1,1,1)钠署，vb=(1,1,2)糠聪，vc=(2,1,1)，vd=(1,2,4)谐鼎，這時(shí)若用這四個(gè)向量表示v1則有無(wú)窮多種表示方法舰蟆，如v1=2va+ vb－vc或v1=va－vb + vd。當(dāng)然狸棍，這個(gè)用矩陣的方法理解更容易身害，第一種情況求系數(shù)是三個(gè)獨(dú)立的方程三個(gè)未知數(shù)，解唯一草戈；第二種情況是三個(gè)方程四個(gè)未知數(shù)塌鸯，必然有無(wú)窮多解。問(wèn)題來(lái)了唐片，既然在冗余字典里信號(hào)表達(dá)時(shí)的系數(shù)不唯一丙猬，那么是否存在一種最好的表達(dá)方式呢？
? ? ? ??定義表達(dá)你的信號(hào)空間的歸一化基本模塊作為字典费韭。這些歸一化向量叫做原子茧球。如果字典的原子張成了整個(gè)信號(hào)空間，那么字典就是完全的星持。如果有原子之間線(xiàn)性相關(guān)抢埋，那么字典就是冗余的。在大多數(shù)匹配追蹤的應(yīng)用中，字典都是完全且冗余的揪垄。
? ? ? ??假設(shè){φk}表示字典的原子穷吮。假設(shè)字典是完全且冗余的。使用字典的線(xiàn)性組合表達(dá)信號(hào)將是不唯一的福侈。
? ? ? ??一個(gè)重要的問(wèn)題是是否存在一種最好的表達(dá)方式酒来，一種直觀(guān)的最好方式是選擇φk使得近似信號(hào)和原始信號(hào)有最大的內(nèi)積，如最好的φk滿(mǎn)足
即對(duì)于正交原子肪凛，為投影到由φk張成的子空間上的幅值堰汉。匹配追蹤的中心問(wèn)題是你如何選擇信號(hào)在字典中最優(yōu)的M個(gè)展開(kāi)項(xiàng)。
? ? ? ??用Φ={φk}表示一個(gè)原子歸一化的字典伟墙，x表示信號(hào)翘鸭。匹配追蹤算法流程如下：
? ? ? ??(1)首先初始化殘差e0=x；
? ? ? ??(2)匹配追蹤的第一步是從字典中找到與e0的內(nèi)積絕對(duì)值最大的原子戳葵，表示為φ1就乓；
? ? ? ??(3)通過(guò)從e0減去其在φ1所張成空間上的正交投影得到殘差e1；
其中<e0, φ1>表示e0與φ1的內(nèi)積拱烁。由于已經(jīng)說(shuō)明Φ為原子歸一化的字典生蚁，即φ1為單位列向量，所以e0在φ1所張成空間上的正交投影可以表示為<e0, φ1>φ1 (由于為一個(gè)向量戏自，所以e0在φ1所張成空間上的正交投影即為e0在φ1方向上的正交投影分量)邦投，若φ1不是單位列向量，則e0在φ1方向上的正交投影分量為<e0, φ1>φ1/(φ1Tφ1)擅笔，其中上標(biāo)T表示轉(zhuǎn)置志衣。
? ? ? ??(4)對(duì)殘差迭代執(zhí)行(2)、(3)步猛们；
其中φm+1是從字典中找到與em的內(nèi)積絕對(duì)值最大的原子念脯。
? ? ? ??(5)直到達(dá)到某個(gè)指定的停止準(zhǔn)則后停止算法。
? ? ? ??第(4)步對(duì)應(yīng)到文獻(xiàn)[1]中為文獻(xiàn)中的式(12)：
? ? ? ??這里要從數(shù)學(xué)上說(shuō)明一點(diǎn)：由于內(nèi)積<em, φm+1>實(shí)際上為一個(gè)數(shù)(標(biāo)量)弯淘，所以
若令P=φm+1 (φTm+1φm+1)-1φTm+1 绿店，則第(4)步的迭代公式為
注意矩陣P每次迭代都是不同的。
? ? ? ??在匹配追蹤(MP)中耳胎，字典原子不是相互正交的向量惯吕。因此上面減去投影計(jì)算殘差的過(guò)程中會(huì)再次引入與前面使用的原子不正交的成分。
? ? ? ??為了更形象的看出這點(diǎn)怕午，下面的例子非常清晰地展現(xiàn)了每一個(gè)細(xì)節(jié)（這是一個(gè)類(lèi)比废登，不代表匹配追蹤過(guò)程）∮粝В考慮歐幾里得二維空間的歸一化框架字典
假設(shè)有如下信號(hào)
這個(gè)例子如圖所示堡距，字典原子以紅色表示甲锡，信號(hào)向量藍(lán)色表示
? ? ? ??在matlab中如下創(chuàng)建此字典和信號(hào)：
dictionary = [1 0; 1/2 sqrt(3)/2;-1/sqrt(2) -1/sqrt(2)]';
x = [1 1/2]';
? ? ? ??計(jì)算信號(hào)和字典原子的內(nèi)積：
scalarproducts = dictionary'*x
? ? ? ??輸出為：
scalarproducts =
? ? 1.0000
? ? 0.9330
? -1.0607
? ? ? ??內(nèi)積絕對(duì)值最大的是第3個(gè)原子[－1/sqrt(2); －1/sqrt(2)]。這也很容易理解羽戒，因?yàn)閮烧咧g的夾角約為π缤沦。從信號(hào)中減去信號(hào)在此原子上的正交投影得到殘差：
residual = x-scalarproducts(3)*dictionary(:,3)
? ? ? ??輸出為：
residual =
? ? 0.2500
? -0.2500
? ? ? ??可以驗(yàn)證此時(shí)residual與原子[－1/sqrt(2); －1/sqrt(2)]正交。再次計(jì)算信號(hào)和字典原子的內(nèi)積：
scalarproducts = dictionary'*x
? ? ? ?輸出為：
?scalarproducts =
? ? 0.2500
? -0.0915
? -0.0000
? ? ? ? 第三個(gè)值為0易稠，這也說(shuō)明residual與第3個(gè)原子[-1/sqrt(2);-1/sqrt(2)]正交「追希現(xiàn)在內(nèi)積絕對(duì)值最大的是第1個(gè)原子[1;0]（內(nèi)積值為0.25）。從信號(hào)中減去信號(hào)在此原子上的正交投影得到殘差：
residual = residual-scalarproducts(1)*dictionary(:,1)
? ? ? ??輸出為：
residual =
? ? ? ? 0
? -0.2500
? ? ? ??可以驗(yàn)證此時(shí)residual與第1個(gè)原子[1;0]正交驶社。再次計(jì)算信號(hào)和字典原子的內(nèi)積：
scalarproducts = dictionary'*x
? ? ? ??輸出為：
scalarproducts =
? ? ? ? 0
? -0.2165
? ? 0.1768
? ? ? ??第一個(gè)值為0企量，這也說(shuō)明residual與第1個(gè)原子[1;0]正交；但同時(shí)也應(yīng)當(dāng)注意到第3個(gè)值不再為零亡电，也就是說(shuō)此時(shí)residual與第3個(gè)原子[-1/sqrt(2);-1/sqrt(2)]不再正交届巩，這意味著這個(gè)原子有可能被再次選擇到。這個(gè)事實(shí)及與收斂相關(guān)的問(wèn)題在正交匹配追蹤中得以解決》萜梗現(xiàn)在內(nèi)積絕對(duì)值最大的是第2個(gè)原子[1/2;sqrt(3)/2]（內(nèi)積值為-0.2165）恕汇。從信號(hào)中減去信號(hào)在此原子上的正交投影得到殘差：
residual = residual-scalarproducts(2)*dictionary(:,2)
? ? ? ??輸出為：
residual =
? ? 0.1083
? -0.0625
? ? ? ??到這里可能你會(huì)發(fā)現(xiàn)一個(gè)問(wèn)題：三個(gè)原子都選擇了一遍，但現(xiàn)在的residual仍然不為零或辖；而二維空間中任意兩個(gè)線(xiàn)性無(wú)關(guān)的二維向量即可生成整個(gè)二維空間瘾英，在這個(gè)例子中，三個(gè)列向量（原子）兩兩都是線(xiàn)性無(wú)關(guān)的颂暇，即我們要稀疏表達(dá)的信號(hào)x實(shí)際上可以由三個(gè)原子中任意兩個(gè)的線(xiàn)性組合表示方咆。將三個(gè)原子中的任意兩個(gè)再加上信號(hào)x組合一個(gè)矩陣，用Matlab函數(shù)rref將矩陣化簡(jiǎn)為行階梯最簡(jiǎn)形(Reduced row echelon form)：
rref([dictionary(:,1) dictionary(:,2) x])
rref([dictionary(:,1) dictionary(:,3) x])
rref([dictionary(:,2) dictionary(:,3) x])
? ? ? ??三行代碼的輸出分別為：
ans =
? ? 1.0000? ? ? ? 0? ? 0.7113
? ? ? ? 0? ? 1.0000? ? 0.5774
ans =
? ? 1.0000? ? ? ? 0? ? 0.5000
? ? ? ? 0? ? 1.0000? -0.7071
ans =
? ? 1.0000? ? ? ? 0? -1.3660
? ? ? ? 0? ? 1.0000? -2.3801
? ? ? ??由以上輸出可以知道用其中任意兩個(gè)原子表示x時(shí)的組合系數(shù)蟀架。下面把上面所有程序給出：
dictionary = [1 0; 1/2 sqrt(3)/2;-1/sqrt(2) -1/sqrt(2)]';
x = [1 1/2]';
scalarproducts = dictionary'*x
residual = x-scalarproducts(3)*dictionary(:,3)
scalarproducts = dictionary'*residual
residual = residual-scalarproducts(1)*dictionary(:,1)
scalarproducts = dictionary'*residual
residual = residual-scalarproducts(2)*dictionary(:,2)
rref([dictionary(:,1) dictionary(:,2) x])
rref([dictionary(:,1) dictionary(:,3) x])
rref([dictionary(:,2) dictionary(:,3) x])
? ? ? ??雖然所有3個(gè)原子均選擇了一遍后殘差仍不為零，但如果繼續(xù)迭代榆骚，殘差將會(huì)越來(lái)越小片拍，程序如下（以下迭代5*N=5*3=15次）：
dictionary = [1 0; 1/2 sqrt(3)/2;-1/sqrt(2) -1/sqrt(2)]';
x = [1 1/2]';
[M,N] = size(dictionary);
residual = zeros(M,5*N+1);%殘差矩陣，保存每次迭代后的殘差
residual(:,1) = x;%初始化殘差為x
L = size(residual,2);%得到殘差矩陣的列
pos_mm = zeros(1,L);%用來(lái)保存每次選擇的列序號(hào)
resi_norm = zeros(1,L);%用來(lái)保存每次迭代后的殘差的2范數(shù)
resi_norm(1) = norm(x);%因?yàn)榍懊嬉殉跏蓟瘹埐顬閤
for mm = 1:length(residual)
? ? %求出dictionary每列與上次殘差的內(nèi)積
? ? scalarproducts = dictionary'*residual(:,mm);
? ? %找到內(nèi)積中最大的列及其內(nèi)積值
? ? [val,pos] = max(abs(scalarproducts));
? ? %更新殘差
? ? residual(:,mm+1) = residual(:,mm) - ...
? ? ? ? scalarproducts(pos)*dictionary(:,pos);
? ? %計(jì)算殘差的2范數(shù)（平方和再開(kāi)根號(hào)）
? ? resi_norm(mm+1) = norm(residual(:,mm+1));
? ? %保存選擇的列序號(hào)
? ? pos_mm(mm+1) = pos;
end
%繪出殘差的2范數(shù)曲線(xiàn)
plot(resi_norm);grid;
? ? ? ??將pos_mm輸出：pos_mm=[0 3 1 2 1 2 1 2 1 2 1 2 1 2 1 2 1]妓肢，從第4次往后捌省，交替選擇第1個(gè)原子和第2個(gè)原子。有關(guān)采用匹配追蹤算法在冗余字典中稀疏表達(dá)實(shí)際的信號(hào)可參見(jiàn)參考文獻(xiàn)[3]碉钠。
? ? ? ??在正交匹配追蹤OMP中纲缓，殘差是總與已經(jīng)選擇過(guò)的原子正交的。這意味著一個(gè)原子不會(huì)被選擇兩次喊废，結(jié)果會(huì)在有限的幾步收斂祝高。OMP的算法如下
? ? ? ??(1)用x表示你的信號(hào)，初始化殘差e0=x污筷；
? ? ? ??(2)選擇與e0內(nèi)積絕對(duì)值最大的原子工闺，表示為φ1；
? ? ? ??(3)將選擇的原子作為列組成矩陣Φt，定義Φt列空間的正交投影算子為
通過(guò)從e0減去其在Φt所張成空間上的正交投影得到殘差e1陆蟆；
? ? ? ??(4)對(duì)殘差迭代執(zhí)行(2)雷厂、(3)步；
其中I為單位陣叠殷。需要注意的是在迭代過(guò)程中Φt為所有被選擇過(guò)的原子組成的矩陣改鲫，因此每次都是不同的，所以由它生成的正交投影算子矩陣P每次都是不同的林束。有關(guān)正交投影的概念可以參見(jiàn)《壓縮感知中的數(shù)學(xué)知識(shí)：投影矩陣（projectionmatrix）》像棘。
? ? ? ??(5)直到達(dá)到某個(gè)指定的停止準(zhǔn)則后停止算法。
? ? ? ??可以看出OMP與MP的不同根本在于殘差更新過(guò)程：OMP減去的Pem是em在所有被選擇過(guò)的原子組成的矩陣Φt所張成空間上的正交投影诊县，而MP減去的Pem是em在本次被選擇的原子φm所張成空間上的正交投影讲弄。基于此依痊，OMP可以保證已經(jīng)選擇過(guò)的原子不會(huì)再被選擇避除。
? ? ? ??有關(guān)正交匹配追蹤的進(jìn)一步分析見(jiàn)下一篇《施密特(Schimidt)正交化與正交匹配追蹤》，有關(guān)究竟為什么把稀疏表示與匹配追蹤放在一起參見(jiàn)下下篇《正交匹配追蹤(OMP)在稀疏分解與壓縮感知重構(gòu)中的異同》胸嘁。
? ? ? ??推薦一篇有關(guān)匹配追蹤與正交匹配追蹤的詳細(xì)分析博文瓶摆，參見(jiàn)文獻(xiàn)[4]，寫(xiě)的很不錯(cuò)性宏。
【1】S Mallat, Z Zhang. Matchingpursuit with time-frequency dictionaries[J]. IEEE Transactions on SignalProcessing, 1993, 41(12): 3397-3415.
【2】了凡春秋. Matlab匹配追蹤(Matching Pursuit) 之一群井，https://chunqiu.blog.ustc.edu.cn/
【3】了凡春秋. Matlab匹配追蹤(Matching Pursuit) 之二，https://chunqiu.blog.ustc.edu.cn/
【4】逍遙劍客. MP算法和OMP算法及其思想毫胜，http://blog.csdn.net/scucj/article/details/7467955
