OMP 算法

問題描述

給定一個過完備字典矩陣，其中它的每列表示一種原型信號的原子。給定一個信號y，它可以被表示成這些原子的稀疏線性組合。信號 y 可以被表達為 y = Dx ，或者。 字典矩陣中所謂過完備性，指的是原子的個數遠遠大於信號y的長度(其長度很顯然是n)，即n<<k。

主要問題

數學模型

輸入輸出

• 輸入：測量值y，字典D；

• 輸出：稀疏稀疏向量x；

算法求解

算法流程

代碼實現

function [ theta ] = CS_OMP( y,A,t )
%CS_OMP Summary of this function goes here
%Version: 1.0 written by jbb0523 @2015-04-18
% Detailed explanation goes here
% y = Phi * x
% x = Psi * theta
% y = Phi*Psi * theta
% 令 A = Phi*Psi, 則y=A*theta
% 現在已知y和A，求theta
    [y_rows,y_columns] = size(y);
    if y_rows<y_columns
        y = y';%y should be a column vector
    end
    [M,N] = size(A);%傳感矩陣A為M*N矩陣
    theta = zeros(N,1);%用來存儲恢復的theta(列向量)
    At = zeros(M,t);%用來迭代過程中存儲A被選擇的列
    Pos_theta = zeros(1,t);%用來迭代過程中存儲A被選擇的列序號
    r_n = y;%初始化殘差(residual)為y
    for ii=1:t%迭代t次，t為輸入參數
        product = A'*r_n;%傳感矩陣A各列與殘差的內積
        [val,pos] = max(abs(product));%找到最大內積絕對值，即與殘差最相關的列
        At(:,ii) = A(:,pos);%存儲這一列
        Pos_theta(ii) = pos;%存儲這一列的序號
        A(:,pos) = zeros(M,1);%清零A的這一列，其實此行可以不要，因為它與殘差正交
        %y=At(:,1:ii)*theta，以下求theta的最小二乘解(Least Square)
        theta_ls = (At(:,1:ii)'*At(:,1:ii))^(-1)*At(:,1:ii)'*y;%最小二乘解
        %At(:,1:ii)*theta_ls是y在At(:,1:ii)列空間上的正交投影
        r_n = y - At(:,1:ii)*theta_ls;%更新殘差 
    end
    theta(Pos_theta)=theta_ls;%恢復出的theta
end

總結

參考文獻

1. MP算法和OMP算法及其思想:http://blog.csdn.net/scucj/article/details/7467955

2. Orthogonal Matching Pursuit(OMP)正交匹配追蹤算法學習筆記:http://blog.sciencenet.cn/blog-810210-653094.html