• 壓縮感知系列博客：
• 壓縮感知入門①從零開始壓縮感知
• 壓縮感知入門②信號的稀疏表示和約束等距性
• 壓縮感知入門③基于ADMM的全變分正則化的壓縮感知重構(gòu)算法
• 壓縮感知入門④基于總體最小二乘的擾動壓縮感知重構(gòu)算法

1. Problem

信號壓縮是是目前信息處理領(lǐng)域非常成熟的技術(shù)，其主要原理是利用信號的稀疏性。一個稀疏信號的特征是，信號中有且僅有少量的位置是有值的，其它位置都是零。對于一個稀疏的信號，在存儲時只需要記錄有值的位置，從而實現(xiàn)對原始信號的壓縮。

對于原本不稀疏的信號，可以利用一種字典（正交變換基，例如傅里葉、小波）對信號進行線性表示，得到原始信號在這個稀疏域上的稀疏表示系數(shù)，只需要記錄這些少量的系數(shù)就能夠?qū)崿F(xiàn)對信號的壓縮存儲。

在信號重建時，首先對信號補零得到原始維度的信號，由于采用的變換字典是正交的，可以通過正交的變換得到原始信號。

壓縮感知與傳統(tǒng)的信號壓縮有著異曲同工之妙，而不同之處在于，壓縮感知的信號壓縮過程是將原始的模擬信號直接進行壓縮（采樣即壓縮）。而傳統(tǒng)信號壓縮通常是先將信號采樣后再進行壓縮（采樣后壓縮），這種壓縮方式的主要問題在于采用較高的采樣資源將信號采集后，又在信號壓縮過程將費盡心思采集到的信號丟棄了，從而造成資源的浪費¹。

以經(jīng)典的圖像采集為例，對于一副$m\times n$的圖像信號$I$進行采集，根據(jù)正交變換的思想，至少要對其進行$N=m\times n$次采樣才能夠得到這副圖像。

若該圖像是稀疏的，根據(jù)壓縮感知理論，可以至少進行$M$次采樣就能夠采集到該信號，其中$M<<N$，$M$的值有信號的稀疏性決定。

2. Formulation

一個壓縮感知采樣過程可以表示為

$$Ax=b+e\text{\hspace{1em}(1)}$$

其中，$A\in {\mathbb{R}}^{M\times N}$表示觀測矩陣，$x\in {\mathbb{R}}^{N\times 1}$表示$I\in {\mathbb{R}}^{m\times n}$的向量形式，$b\in {\mathbb{R}}^{M\times 1}$表示測量向量，$e\in {\mathbb{R}}^{M\times 1}$表示測量噪聲?？梢钥吹?，一個$N$維的信號在采樣過程中被壓縮成$M$維的信號，這里的$M<N$。

信號重構(gòu)過程可以表示為一個優(yōu)化問題，利用信號的梯度稀疏性質(zhì)，可以構(gòu)建目標(biāo)函數(shù)：

$$\underset{x}{\min }||Dx|{|}_{1}s.t.Ax=b+e\text{\hspace{1em}(2)}$$

其中，$x$為待求解信號，$D$為全變分算子²。去除約束有

$$\underset{x}{\min }\lambda ||Dx|{|}_1+\frac{1}{2}||Ax?b|{|}_2^2\text{\hspace{1em}(3)}$$

該問題是一個非凸、不光滑問題，無法直接采用梯度下降法求解。引入變量$d$，將問題(3)轉(zhuǎn)化為ADMM的一般形式³

$$\underset{x}{\min }\frac{1}{2}||Ax?b|{|}_2^2+\lambda ||d|{|}_{1}s.t.Dx?d=0\text{\hspace{1em}(4)}$$

利用增廣拉格朗日法引入凸松弛，同時去除約束條件，有

$$L(x,d,\mu )=\frac{1}{2}||Ax?b|{|}_2^2+\lambda ||d|{|}_1+{\mu }^T(Dx?d)+\frac{\delta }{2}||Dx?d|{|}_2^2\text{\hspace{1em}(5)}$$

其中$\mu$為拉格朗日乘子，$\delta >0$為拉格朗日懲罰項。為了使表達更簡潔，可做如下替換：

$$L(x,d,\mu )=\frac{1}{2}||Ax?b|{|}_2^2+\lambda ||d|{|}_1+\frac{\delta }{2}||Dx?d+p|{|}_2^2?\frac{\delta }{2}||p|{|}_2^2\text{\hspace{1em}(6)}$$

其中$p=\mu /\delta$。利用ADMM，問題(6)的求解可通過交替求解以下三個問題進行實現(xiàn)：

$$x_{n+1}=arg\underset{x}{\min }\frac{1}{2}||Ax?b|{|}_2^2+\frac{\delta }{2}||Dx?d_n+p_n|{|}_2^2\text{\hspace{1em}(7)}$$

$$d_{n+1}=arg\underset{u}{\min }\lambda ||d|{|}_1+\frac{\delta }{2}||D{x}_n?d+p_n|{|}_2^2\text{\hspace{1em}(8)}$$

$$p_{n+1}=p_n+(D{x}_{n+1}?d_{n+1})\text{\hspace{1em}(9)}$$

3. Simulation

測試圖像采用的是house，分別測試在20%、40%、60%、80%和100%采樣率時，壓縮感知重構(gòu)算法的圖像恢復(fù)結(jié)果。

從仿真結(jié)果可以看到，在20%采樣率時，信號的基本輪廓信息就被成功采集了，當(dāng)采樣率達到60%以上時，繼續(xù)增加采樣率并沒有使得圖像更加的清晰，也就是說，針對這副圖像，若采用傳統(tǒng)的正交采集的方式，有將近一半的采樣資源是被浪費的。

4. Algorithm

ratio=0.5;%采樣率
x=double(imread('house.bmp'));
[m,n]=size(x);
N=m*n;
M=floor(ratio*N);
A=rand(M,N);%觀測矩陣
e=50*randn(M,1);%測量噪聲
% 信號采集過程，利用線性投影對信號x進行采集，同時考慮了測量噪聲
b=A*reshape(x,N,1)+e;
% 信號重構(gòu)過程，利用僅有的M個測量值恢復(fù)維度為N的信號
x_r=ADMM_TV_reconstruct(A,b,300,500,100);
figure;
subplot(121);
imshow(uint8(x));
title('原始圖像');
subplot(122);
imshow(uint8(reshape(x_r,m,n)));
title(sprintf('重構(gòu)圖像（%d%%采樣率）',ratio*100));

function xp=ADMM_TV_reconstruct(A,b,delta,lambda,iteratMax)
    [~,N]=size(A);
    [Dh,Dv]=TVOperatorGen(sqrt(N));
    D=sparse([(Dh)',(Dv)']');
    d=D*ones(N,1);
    p=ones(2*N,1)/delta;
    invDD=inv(A'*A+delta*(D'*D));
    for ii=1:iteratMax
        x=invDD*(A'*b+delta*D'*(d-p));
        d=wthresh(D*x+p,'s',lambda/delta);
        p=p+D*x-d;
    end
    xp=x;
end

function [Dh,Dv]=TVOperatorGen(n)
    Dh=-eye(n^2)+diag(ones(1,n^2-1),1);
    Dh(n:n:n^2,:)=0;
    Dv=-eye(n^2)+diag(ones(1,n^2-n),n);
    Dv(n*(n-1)+1:n^2,:)=0;
end

Github link: https://github.com/dwgan/ADMM_TV_reconstruct

參考文獻

到了這里，關(guān)于壓縮感知入門③基于ADMM的全變分正則化的壓縮感知重構(gòu)算法的文章就介紹完了。如果您還想了解更多內(nèi)容，請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！