基追踪及其实现

基追踪

我们将(l_1)范数替换(l_0)范数以后，稀疏表征模型可以表示为：

[min |alpha|_1 quad mathrm{s.t.} ; Phialpha = s ]

这是一个二次规划问题，如何将(l_1)范数优化问题转为线性规划问题呢？

参考Atomic Decomposition by Basis Pursuit中的方法，可以将(l_1)范数优化问题转化为一个常见的线性规划问题，然后我们可以用单纯形法或者内点法来求解.

(l_1)范数优化转换为线性规划问题

上面的(l_1)范数优化模型等价于下面的线性规划问题：

[min c^{mathrm{T}}x quad Ax = b, ; x geq 0 ]

转换规则（(p)为(Phi)的列数）：

[m Leftrightarrow 2p qquad A Leftrightarrow (Phi, Phi) qquad b Leftrightarrow s qquad c Leftrightarrow (1; 1) qquad x Leftrightarrow (u; v) qquad alpha Leftrightarrow u - v ]

MATLAB仿真

基追踪实现（使用MATLAB内置函数linprog进行线性优化）

function [ alpha ] = bpalg( s, Phi )
% 使用BP思想计算稀疏系数
% 参考文献：Chen, S.S., Donoho, D.L. and Saunders, M.A., 2001. Atomic decomposition by basis pursuit. SIAM review, 43(1), pp.129-159.
% s：信号向量
% Phi：字典矩阵
% alpha：稀疏系数向量
    
    p = size(Phi, 2);
    c = ones(2 * p, 1);
    A = [Phi, -Phi];
    lb = zeros(2 * p, 1); 
    x = linprog(c, [], [], A, s, lb);
    alpha = x(1 : p) - x(p + 1 : 2 * p);

end

图像重建实现

% 利用稀疏理论进行图像重建测试
clear; clc;

im = imread('https://upload.wikimedia.org/wikipedia/en/2/24/Lenna.png');  % 读入图像
figure(),
subplot(131),
imshow(im);
title('原始彩色图像');


imgrey = rgb2gray(im);              % 转为灰度图像
subplot(132),
imshow(imgrey);
title('原始灰度图像');

bsize = 8;
imcols = im2col(im2double(imgrey), [bsize, bsize], 'distinct');

codebook = dctmtx(bsize ^ 2);                   % DCT字典矩阵

% 稀疏求解，需要对imcols的行进行遍历
cols = size(imcols, 2);
sparse = zeros(size(imcols));
for i = 1 : cols
    sparse(:, i) = bpalg(imcols(:, i), codebook);
end

% 图像重建
imrecons = codebook * sparse;
imrecons = col2im(imrecons, [bsize bsize], size(imgrey), 'distinct');
subplot(133),
imshow(imrecons);
title('稀疏重建图像');

运行结果：

参考资料

压缩感知的常见稀疏基名称及离散傅里叶变换基

Atomic Decomposition by Basis Pursuit