【压缩感知合集9】压缩感知的OMP算法（算法步骤分析、举例分析、说明总结和缺陷）

0 前情提要

0.1 数学模型和总体框图如下

给定输入信号$\mathbf{X}\in {\mathbb{R}}^{N\times 1}$，最终想要得到压缩信号$\mathbf{A}\in {\mathbb{R}}^{M\times 1}$，$K<<N$

0.2 压缩过程图例分析如下

整个压缩过程也可以被称为感知过程
$$\mathbf{A}=\mathbf{\Phi }\mathbf{X}=\mathbf{\Phi }\mathbf{\Psi }\mathbf{Y}=\mathbf{\Theta }\mathbf{Y}$$

$\mathbf{\Theta }$即为感知过程的核心命名为感知矩阵

符号	含义	维度	属性
$\mathbf{X}$	输入信号；待压缩信号	${\mathbb{R}}^{N\times 1}$	未知，需要恢复
$\mathbf{\Phi }$	观测矩阵；测量矩阵	${\mathbb{R}}^{M\times N}$	已知（非自适应性）
$\mathbf{\Psi }$	变换矩阵；变换基矩阵；稀疏基矩阵；稀疏矩阵；正交基字典矩阵	${\mathbb{R}}^{N\times N}$	已知（非自适应性）
$\mathbf{Y}$	正交基变换后的稀疏表示	${\mathbb{R}}^{N\times 1}$	未知，需要恢复
$\mathbf{\Theta }$	感知矩阵，传感矩阵	${\mathbb{R}}^{M\times N}$	已知（非自适应性）
$\mathbf{A}$	观测压缩所得到压缩信号	${\mathbb{R}}^{M\times 1}$	已知

0.3 算法重构恢复过程如下

在得到已经压缩完的采样信号 $\mathbf{A}$ 后，根据确定的固定性观测矩阵 $\mathbf{\Phi }$ 和稀疏矩阵 $\mathbf{\Psi }$ 的先验信息进行恢复，数学表达如下
$$\check{\mathbf{X}}=f(\mathbf{A},\mathbf{\Theta })$$
给定被压缩的信号 $\mathbf{A}$ 和感知矩阵 $\mathbf{\Theta }$，求解输入原始信号 $\check{\mathbf{X}}$ 的过程称为重构。

重构问题相较于压缩问题来说是一个更加困难的一个任务。

由于$M<<N$， 已知条件所能构成的方程是欠定的，无法轻易求出数值解的，幸运的是，如果原始信号是稀疏的，那么这个问题可以被很好地解决。

解释一下为什么是欠定的：（ $\mathbf{X}$ 满足的约束如下）
$$\begin{gathered} \mathbf{A}=\mathbf{\Psi }\mathbf{X} \\ \left[\begin{array}{c}{a}_1\\ ⋮\\ a_M\end{array}\right]=\left[\begin{array}{cccc}{\psi }_{11}& {\psi }_{12}& ...& {\psi }_{1N}\\ {\psi }_{2N}& {\psi }_{2N}& ...& {\psi }_{2N}\\ ⋮& ⋮& ⋮& ⋮\\ {\psi }_{M1}& {\psi }_{M2}& ...& {\psi }_{MN}\end{array}\right]\left[\begin{array}{c}{x}_1\\ x_2\\ ⋮\\ x_N\end{array}\right] \end{gathered}$$
实际使用过程中需要求解出 $N$ 个未知数，但是只有 $M$ 个方程，未知数的个数远远大于方程的个数。

若$N=M$，则可轻松由$\mathbf{A}$解出$\mathbf{X}$和$\mathbf{Y}$

而$M<<N$，可根据稀疏表示下的信号$\mathbf{Y}$和矩阵所具有的RIP性质重构

一般可以抽象为如下求解任务
$$\begin{gathered} \min {\Vert {\mathbf{\Psi }}^T\mathbf{X}\Vert }_0 \\ s.t.\mathbf{\Theta }\mathbf{X}=\mathbf{\Phi }\mathbf{\Psi }\mathbf{X}=\mathbf{A} \end{gathered}$$

1 正交匹配跟踪算法（OMP）

1.0 问题设置

考虑以下情况，给定 $\mathbf{X}=\left[\begin{array}{c}-1.2\\ 1\\ 0\end{array}\right]$ 和 $\mathbf{D}=\left[\begin{array}{ccc}-0.707& 0.8& 0\\ 0.707& 0.6& -1\end{array}\right]$ ，

我们可以计算 $\mathbf{A}=\mathbf{D}\mathbf{X}$ ，这是非常简单地，也就是
$$\mathbf{A}=\mathbf{D}\mathbf{X}=\left[\begin{array}{ccc}-0.707& 0.8& 0\\ 0.707& 0.6& -1\end{array}\right]\left[\begin{array}{c}-1.2\\ 1\\ 0\end{array}\right]=\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]$$
现在让我们介绍比较难的部分，给定 $\mathbf{A}=\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]$ 和 $\mathbf{D}=\left[\begin{array}{ccc}-0.707& 0.8& 0\\ 0.707& 0.6& -1\end{array}\right]$，如何找到最接近的 $\mathbf{X}$ ？

解释一下为什么是欠定的：
$$\begin{array}{rl}\mathbf{A}& =\mathbf{D}\mathbf{X}\\ \left[\begin{array}{c}1.65\\ -0.25\end{array}\right]& =\left[\begin{array}{ccc}-0.707& 0.8& 0\\ 0.707& 0.6& -1\end{array}\right]\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right]\end{array}$$
例子中需要求解出3个未知数，但是只有2个方程

实际使用过程中需要求解的未知数有 $N$ 个，但是只有 $M$ 个方程，未知数的个数远远大于方程的个数 $M<<N$ 。

1.1 步骤描述

输入：字典矩阵 $\mathbf{D}$ ，采样所得向量 $\mathbf{A}$ ，稀疏度 $K$

输出： $\mathbf{X}$ 的 $K$ 稀疏逼近 $\check{\mathbf{X}}$

初始化：残差 ${\mathbf{f}}_0=\mathbf{A}$ ，索引集 ${\Lambda }_0=\varnothing$ ， $t=1$

循环执行步骤1-5：

• 找出残差 ${\mathbf{f}}_t$ 和字典矩阵的列 ${\mathbf{d}}_{r_t}$ 内积中最大值所对应的脚标 $r_t$ ，即 $r_i=\arg {\max }_{r_t=1,\cdots ,N}|⟨\mathbf{f},{\mathbf{d}}_{r_t}⟩|$
• 更新索引集 ${\Lambda }_t={\Lambda }_{t-1}\cup r_t$ ，记录找到的字典矩阵中的重建原子集合 ${\mathbf{D}}_t=\left[{\mathbf{D}}_{t-1},{\mathbf{d}}_{r_t}\right]$
• 由最小二乘得到 $\check{\mathbf{X}}=\arg {\min }_{\mathbf{w}}\Vert \mathbf{A}-{\mathbf{D}}_t\mathbf{w}\Vert$
• 更新残差 ${\mathbf{f}}_t=\mathbf{A}-{\mathbf{A}}_t\check{\mathbf{X}}$ ， $t=t+1$
• 判断是否满足 $t>K$ ，若满足，（或者判断残差的大小）则迭代停止；若不满足，则继续循环

1.2 匹配追踪的概念

这一部分是和MP算法相同的，正好也在此说明匹配追踪的算法。

字典由一系列基信号或者原子组成的，也就是说
$$\mathbf{D}=\left[\begin{array}{ccc}-0.707& 0.8& 0\\ 0.707& 0.6& -1\end{array}\right]=\left[\begin{array}{c}{\mathbf{d}}_1,{\mathbf{d}}_2,{\mathbf{d}}_3\end{array}\right]$$
其中，
$${\mathbf{d}}_1=\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right],{\mathbf{d}}_2=\left[\begin{array}{c}0.8\\ 0.6\end{array}\right],{\mathbf{d}}_3=\left[\begin{array}{c}0\\ -1\end{array}\right]$$
被称为原子或者基信号。

现在我们将基信号表达出来
$$\begin{array}{rl}\mathbf{A}& =\mathbf{D}\mathbf{X}=\left[\begin{array}{ccc}-0.707& 0.8& 0\\ 0.707& 0.6& -1\end{array}\right]\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right]=\left[\begin{array}{c}{x}_1{\mathbf{d}}_1,x_2{\mathbf{d}}_2,x_3{\mathbf{d}}_3\end{array}\right]\\ & =\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right]x_1+\left[\begin{array}{c}0.8\\ 0.6\end{array}\right]x_2+\left[\begin{array}{c}0\\ -1\end{array}\right]x_3\end{array}$$
上面的方程说明了 $\mathbf{A}$ 是 $\mathbf{D}$ 中原子的线性组合，而线性组合系数恰好是 $\mathbf{X}$ 。实际上，我们知道 $x_1=-1.2,x_2=1,x_3=0$ ，换句话说原子 ${\mathbf{d}}_1$ 对于 $\mathbf{A}$ 贡献了 $-1.2$ ， ${\mathbf{d}}_2$ 对于 $\mathbf{A}$ 贡献了 $1$ ，${\mathbf{d}}_3$ 对于 $\mathbf{A}$ 贡献了 $0$ 。

MP算法和OMP算法要找到对 $\mathbf{A}$ 贡献最大的那个原子，然后是贡献次之的原子，一直进行到贡献最小的原子。

第一步：求取最大相关

在矩阵A中有三个原子，
$${\mathbf{d}}_1=\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right],\mathbf{d}2=\left[\begin{array}{c}0.8\\ 0.6\end{array}\right],\mathbf{d}3=\left[\begin{array}{c}0\\ -1\end{array}\right]$$
压缩信号 $\mathbf{A}=\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]$

我们计算每个原子对于 $\mathbf{A}$ 的贡献（主要是点乘计算投影的大小）
$$\begin{array}{rl}⟨{\mathbf{d}}_1,\mathbf{A}⟩& =\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right]\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]=-0.707\times 1.65+0.707\times (-0.25)=-1.34\\ ⟨{\mathbf{d}}_2,\mathbf{A}⟩& =\left[\begin{array}{c}0.8\\ 0.6\end{array}\right]\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]=0.8\times 1.65+0.6\times (-0.25)=1.17\\ ⟨{\mathbf{d}}_3,\mathbf{A}⟩& =\left[\begin{array}{c}0\\ -1\end{array}\right]\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]=0\times 1.65+(-1)\times (-0.25)=0.25\end{array}$$

根据上述结果，可以看出 ${\mathbf{d}}_1$ 对于 $\mathbf{A}$ 的贡献最大，值为 $-1.34$ （这里只考虑大小，不用管负号，负号只是代表方向相反。我也想说下自己对于这块的理解，这里使用的是点积（内积），求的是 $\mathbf{A}$ 和 $\mathbf{D}$ 中每个原子的相关性）。当然也可以使用一步矩阵相乘直接计算出贡献
$$\mathbf{w}={\mathbf{D}}^T\mathbf{A}=\left[\begin{array}{cc}-0.707& 0.707\\ 0.8& 0.6\\ 0& -1\end{array}\right]\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]=\left[\begin{array}{c}-1.34\\ 1.17\\ 0.25\end{array}\right]$$
可以得到相同的结果，下图展示了 $\mathbf{D}$ 中的原子 ${\mathbf{d}}_1,{\mathbf{d}}_2,{\mathbf{d}}_3$ 和 $\mathbf{A}$ ，从图中可以看出来 ${\mathbf{d}}_1$ 的负向和 $\mathbf{A}$ 是最接近的，意味着贡献最大最相关。

第二步：求取残差

根据上一节中的计算结果，我们选择贡献最大的原子 ${\mathbf{d}}_1=\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right]$ ，相关的系数是 $⟨{\mathbf{d}}_1,\mathbf{A}⟩=-1.34$ ，如果我们从 $\mathbf{A}$ 中减去，那么剩余的残差是
$$\mathbf{f}=\mathbf{A}-⟨{\mathbf{d}}_1,\mathbf{A}⟩{\mathbf{d}}_1=\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]-(-1.34)\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right]=\left[\begin{array}{c}0.7\\ 0.7\end{array}\right]$$
那么所得的残差由什么意义呢？从 $\mathbf{A}$ 中减去了所有与原子 ${\mathbf{d}}_1$ 有关的信息，为什么说是所有呢，因为刚才计算的系数同时也代表着投影的意义，$\mathbf{A}$ 减去这个长度的投影，所剩残差 $\mathbf{f}$ 与 ${\mathbf{d}}_1$ 正交，也表示残差不能再由 ${\mathbf{d}}_1$ 线性表示。

第三步： 重复迭代

在第一次迭代时，我们选择了原子 ${\mathbf{d}}_1$ ，将其作为一个基放入新的压缩矩阵 ${\mathbf{D}}_{\text{new}}$ 中，也就是 ${\mathbf{D}}_{\text{new}}=[{\mathbf{d}}_1]=\left[\begin{array}{c}-0.707\\ 0.707\end{array}\right]$ ，并且将贡献系数写入到重构信号 $\check{\mathbf{X}}$ 中，$\check{\mathbf{X}}=\left[\begin{array}{c}-1.34\\ 0\\ 0\end{array}\right]$ ， $-1.34$ 被放置到第一个元素的位置是因为这个贡献系数来自于 $\mathbf{D}$ 中的第一个基 ${\mathbf{d}}_1$ 。残差 $\mathbf{f}$ 计算为 $\mathbf{f}=\left[\begin{array}{c}0.7\\ 0.7\end{array}\right]$ 。现在我们要从剩余的原子 ${\mathbf{d}}_2$ 或者 ${\mathbf{d}}_3$ 中选择出对残差贡献最大的，
$$\mathbf{w}={\left[\begin{array}{cc}{\mathbf{d}}_2& {\mathbf{d}}_3\end{array}\right]}^T\mathbf{f}=\left[\begin{array}{cc}0.8& 0.6\\ 0& -1\end{array}\right]\left[\begin{array}{l}0.7\\ 0.7\end{array}\right]=\left[\begin{array}{c}0.98\\ -0.7\end{array}\right]$$
由于 ${\mathbf{f}}_2$ 的贡献比较大（$0.98>0.7$，忽略负号），所以选择 ${\mathbf{f}}_2$ 。

现在我们将已经选择了的基 ${\mathbf{d}}_1,{\mathbf{d}}_2$ 都放到新的压缩矩阵中 ${\mathbf{D}}_{\text{new}}$ ，
$${\mathbf{D}}_{\text{new}}=\left[\begin{array}{cc}{\mathbf{d}}_1& {\mathbf{d}}_2\end{array}\right]=\left[\begin{array}{cc}-0.707& 0.8\\ 0.707& 0.6\end{array}\right]$$

第四步：与MP算法不同的OMP算法细节

接下来这一步是与之前提出的另一种方法匹配跟踪算法不同的地方，

计算 ${\mathbf{D}}_{\text{new}}$ 对于 $\mathbf{A}$ 的贡献，得到系数（而不是MP中的做法，计算 ${\mathbf{d}}_{r_2}$ ，此时 ${\mathbf{d}}_{r_2}即为{\mathbf{d}}_2$ 对残差的贡献，得到系数）。为了得到OMP的新系数，OMP会去解一个最小二乘问题，如下
$$w_1\cdot \left[\begin{array}{c}-0.707\\ 0.707\end{array}\right]+w_2\cdot \left[\begin{array}{c}0.8\\ 0.6\end{array}\right]\text{ as close as possible to }y=\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]$$

写成数学公式如下：
$$\min {\Vert {\mathbf{D}}_{\text{new}}\mathbf{w}-\mathbf{A}\Vert }_2$$
在这个例子中，
$$\min {\Vert \left[\begin{array}{cc}-0.707& 0.8\\ 0.707& 0.6\end{array}\right]\cdot \left[\begin{array}{c}{w}_1\\ w_2\end{array}\right]-\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]\Vert }_2$$

得到 $w_1,w_2$ 。我们知道
$$\min {\Vert {\mathbf{D}}_{\text{new}}\mathbf{w}-\mathbf{A}\Vert }_2$$

的解是 $\mathbf{w}={\mathbf{D}}_{\text{new}}^{+}\mathbf{A}$ ，其中是 ${\mathbf{D}}_{\text{new}}^{+}$ 是 ${\mathbf{D}}_{\text{new}}$ 的伪逆，也就是说 ${\mathbf{D}}_{new}^{+}={\left({\mathbf{D}}_{new}^T{\mathbf{D}}_{new}\right)}^{-1}{\mathbf{D}}_{new}^T$ ，在我们这个例子中
$${\mathbf{D}}_{new}^{+}={\left[\begin{array}{cc}-0.707& 0.8\\ 0.707& 0.6\end{array}\right]}^{+}=\left[\begin{array}{cc}-0.6062& 0.8082\\ 0.7143& 0.7143\end{array}\right]$$
可以在MATLAB中使用pinv()来计算伪逆。计算完伪逆之后，我们得到了
$$\mathbf{w}={\mathbf{D}}_{\text{new}}^{+}\mathbf{A}=\left[\begin{array}{cc}-0.6062& 0.8082\\ 0.7143& 0.7143\end{array}\right]\left[\begin{array}{c}1.65\\ -0.25\end{array}\right]\approx \left[\begin{array}{c}-1.2\\ 1\end{array}\right]$$
在得到更新后的 ${\mathbf{D}}_{\text{new}}$ 和 $\mathbf{w}$ ，我们更新残差
$$\mathbf{f}=\mathbf{A}-{\mathbf{D}}_{\text{new}}\cdot \mathbf{w}\approx \left[\begin{array}{c}1.65\\ -0.25\end{array}\right]-\left[\begin{array}{cc}-0.707& 0.8\\ 0.707& 0.6\end{array}\right]\left[\begin{array}{c}-1.2\\ 1\end{array}\right]\approx \left[\begin{array}{l}0\\ 0\end{array}\right]$$
现在我们得到的残差约等于 $0$ ，所以停止迭代。

我们将得到的 $w_1,w_2$ 分别放入 $\check{\mathbf{X}}$ 中第一和第二个位置，因为它们分别对应了我们选择的第一和第二个原子。接着更新重构信号
$$\check{\mathbf{X}}=\left[\begin{array}{c}-1.2\\ 1\\ 0\end{array}\right]$$

第五步 最后一次迭代

因为残差已经为 $0$ ，所以这一步并不是必须的。许多OMP算法需要设置一个关于信号稀疏性的参数 $K$ ，这告诉算法它需要迭代 $K$ 次，即使残差为 $0$ ，也需要迭代 $K$ 次。

1.3 第二个例子说明

给定 $\mathbf{X}=\left[\begin{array}{c}0\\ 3\\ 1\\ 2\end{array}\right]$ 和 $\mathbf{D}=\left[\begin{array}{cccc}-0.8& 0.3& 1& 0.4\\ -0.2& 0.4& -0.3& -0.4\\ 0.2& 1& -0.1& 0.8\end{array}\right]$ ，有压缩信号 $\mathbf{A}=\mathbf{D}\mathbf{X}=\left[\begin{array}{l}2.7\\ 0.1\\ 4.5\end{array}\right]$ ，现在给定 $\mathbf{A}$ 和 $\mathbf{D}$ ，使用OMP算法求解 $\mathbf{X}$。

有4个基（原子）：
$${\mathbf{d}}_1=\left[\begin{array}{c}-0.8\\ -0.2\\ 0.2\end{array}\right]{\mathbf{d}}_2=\left[\begin{array}{c}0.3\\ 0.4\\ 1\end{array}\right]{\mathbf{d}}_3=\left[\begin{array}{c}1\\ -0.3\\ -0.1\end{array}\right]{\mathbf{d}}_4=\left[\begin{array}{c}0.4\\ -0.4\\ 0.8\end{array}\right]$$
由于基向量的长度不是 $1$ ，所以我们首先进行标准化，实现字典矩阵的标准化 $\mathbf{D}\to \hat{\mathbf{D}}=\left[\widehat{{\mathbf{d}}_1}\widehat{{\mathbf{d}}_2}\widehat{{\mathbf{d}}_3}\widehat{{\mathbf{d}}_4}\right]$ 。

$$\widehat{{\mathbf{d}}_1}={\mathbf{d}}_1/\Vert {\mathbf{d}}_1\Vert =\left[\begin{array}{c}-0.8\\ -0.2\\ 0.2\end{array}\right]/\sqrt{(-0.8{)}^2+(-0.4{)}^2+(0.2{)}^2}=\left[\begin{array}{c}-0.9428\\ -0.2357\\ 0.2357\end{array}\right]$$

$$\widehat{{\mathbf{d}}_2}={\mathbf{d}}_2/\Vert {\mathbf{d}}_2\Vert =\left[\begin{array}{c}0.2680\\ 0.3578\\ 0.8940\end{array}\right]$$

$$\widehat{{\mathbf{d}}_3}={\mathbf{d}}_3/\Vert {\mathbf{d}}_3\Vert =\left[\begin{array}{c}0.9535\\ -0.2860\\ 0.0953\end{array}\right]$$

$$\widehat{{\mathbf{d}}_4}={\mathbf{d}}_4/\Vert {\mathbf{d}}_4\Vert =\left[\begin{array}{c}0.4082\\ -0.4082\\ -0.8165\end{array}\right]$$

标准化的字典基向量 ${\mathbf{d}}_{r_i}$ 对 $\mathbf{A}$ 的贡献
$$\begin{array}{rl}{\hat{\mathbf{D}}}^{\mathrm{T}}\cdot \mathbf{A}& ={\left[\begin{array}{cccc}-0.9428& 0.2680& 0.9535& 0.4082\\ -0.2357& 0.3578& -0.2860& -0.4082\\ 0.2357& 0.9840& -0.0953& -0.8165\end{array}\right]}^{\mathrm{T}}\cdot \left[\begin{array}{c}2.7\\ 0.1\\ 4.5\end{array}\right]\\ & =\left[\begin{array}{c}-1.5085\\ 4.7852\\ 2.1167\\ 4.7357\end{array}\right]\end{array}$$
第二个基向量 $\widehat{{\mathbf{d}}_2}$ 贡献值最大，所以将 ${\mathbf{d}}_2$ 加入到 ${\mathbf{D}}_{\text{new}}$中，${\mathbf{D}}_{\text{new}}=\left[{\mathbf{d}}_2\right]=\left[\begin{array}{c}0.3\\ 0.4\\ 1\end{array}\right]$

利用最小二乘法计算 $\mathbf{w}$
$$\min {\Vert {\mathbf{D}}_{\text{new}}\mathbf{w}-\mathbf{A}\Vert }_2=\min {\Vert {\mathbf{d}}_2{w}_1-\mathbf{A}\Vert }_2$$
也就是求解 $w_1={\mathbf{D}}_{\text{new}}^{+}\mathbf{A}=4.28$

因为 $w_1$ 对应着第二个基向量 ${\mathbf{d}}_2$ ，所以 $\check{\mathbf{X}}=\left[\begin{array}{c}0\\ 4.28\\ 0\\ 0\end{array}\right]$

接下来计算残差
$$\mathbf{f}=\mathbf{A}-{\mathbf{D}}_{\text{new }}\cdot \mathbf{w}=\left[\begin{array}{l}2.7\\ 0.1\\ 4.5\end{array}\right]-\left[\begin{array}{c}0.3\\ 0.4\\ 1\end{array}\right]\cdot 4.28=\left[\begin{array}{c}1.416\\ -1.612\\ 0.22\end{array}\right]$$
接下来重复计算 $\widehat{{\mathbf{d}}_1},\widehat{{\mathbf{d}}_3},\widehat{{\mathbf{d}}_4}$ 对 $\mathbf{f}$ 的贡献

$${\left[\widehat{{\mathbf{d}}_1}\widehat{{\mathbf{d}}_3}\widehat{{\mathbf{d}}_4}\right]}^{\mathrm{T}}\cdot \mathbf{f}={\left[\begin{array}{ccc}-0.9428& 0.9535& 0.4082\\ -0.2357& -0.2860& -0.4082\\ 0.2357& -0.0953& -0.8165\end{array}\right]}^{\mathrm{T}}\cdot \left[\begin{array}{c}1.416\\ -1.612\\ 0.22\end{array}\right]=\left[\begin{array}{c}-0.9032\\ 1.7902\\ 1.4158\end{array}\right]$$
选择第二个贡献最大的基 $\widehat{{\mathbf{d}}_3}$ ，其贡献值为 $1.7902$

将选择的 ${\mathbf{d}}_3$ 加入到 ${\mathbf{D}}_{\text{new}}$ 中，${\mathbf{D}}_{\text{new}}=\left[{\mathbf{d}}_2{\mathbf{d}}_3\right]$

用最小二乘法计算 $\mathbf{w}={\mathbf{D}}_{\text{new}}^{+}\cdot \mathbf{A}=\left[\begin{array}{l}4.1702\\ 1.7149\end{array}\right]$

这个 $\mathbf{w}$ 对应着 ${\mathbf{d}}_2$ 和 ${\mathbf{d}}_3$ ，因此 $\check{\mathbf{X}}=\left[\begin{array}{c}0\\ 4.1702\\ 1.7149\\ 0\end{array}\right]$

接着计算残差

$$\mathbf{f}=\mathbf{A}-{\mathbf{D}}_{\text{new}}\cdot \mathbf{w}=\left[\begin{array}{l}2.7\\ 0.1\\ 4.5\end{array}\right]-\left[\begin{array}{cc}0.3& 1\\ 0.4& -0.3\\ 1& -0.1\end{array}\right]\cdot \left[\begin{array}{c}4.172\\ 1.7149\end{array}\right]=\left[\begin{array}{c}-0.266\\ -1.0536\\ 0.5012\end{array}\right]$$
重复计算 $\widehat{{\mathbf{d}}_1},\widehat{{\mathbf{d}}_4}$ 对 $\mathbf{f}$ 的贡献
$$\left[\widehat{{\mathbf{d}}_1}\widehat{{\mathbf{d}}_4}\right]\cdot \mathbf{f}={\left[\begin{array}{cc}-0.9428& 0.4082\\ -0.2357& -0.4082\\ 0.2357& -0.8165\end{array}\right]}^{\mathrm{T}}\cdot \left[\begin{array}{c}-0.266\\ -1.0536\\ 0.5012\end{array}\right]=\left[\begin{array}{l}0.6172\\ 0.7308\end{array}\right]$$
选择第三个贡献最大的基 $\widehat{{\mathbf{d}}_4}$ ，其贡献值为 $0.7308$

将 ${\mathbf{d}}_4$ 加入到 ${\mathbf{D}}_{\text{new}}$ 中，${\mathbf{D}}_{\text{new}}=\left[{\mathbf{d}}_2{\mathbf{d}}_3{\mathbf{d}}_4\right]$

利用最小二乘计算权重
$$\mathbf{w}={\mathbf{D}}_{\text{new}}^{+}\cdot \mathbf{A}=\left[\begin{array}{l}3\\ 2\\ 1\end{array}\right]$$
这个 $\mathbf{w}$ 对应着 ${\mathbf{d}}_2,{\mathbf{d}}_3$ 和 ${\mathbf{d}}_4$ ，因此 $\check{\mathbf{X}}=\left[\begin{array}{c}0\\ 3\\ 1\\ 2\end{array}\right]$

接着计算残差

$$\mathbf{f}=\mathbf{A}-{\mathbf{D}}_{\text{new}}\cdot \mathbf{w}=\left[\begin{array}{c}2.7\\ 0.1\\ 4.5\end{array}\right]-\left[\begin{array}{ccc}0.3& 1& 0.4\\ 0.4& -0.3& -0.4\\ 1& -0.1& 0.8\end{array}\right]\cdot \left[\begin{array}{l}3\\ 1\\ 2\end{array}\right]=\left[\begin{array}{l}0\\ 0\\ 0\end{array}\right]$$
迭代到此为止，因为此时残差已经为 $0$ ，重建的$\mathbf{X}$ 为 $\check{\mathbf{X}}=\left[\begin{array}{c}0\\ 3\\ 1\\ 2\end{array}\right]$ ，和原来的信号相同。

2 注意事项

• 在计算最大贡献的时候，字典矩阵 $\mathbf{D}$ 中的原子应当是归一化后的，而不是使用归一化前的原子。
• 在计算贡献时，计算的是残差 $\mathbf{f}$ 与矩阵 $\mathbf{D}$ 中每个归一化的原子的点积。
• 在匹配跟踪算法中，重构信号 $\check{\mathbf{X}}$ 是通过计算字典基信号 ${\mathbf{d}}_{r_i}$ 和残差 $\mathbf{f}$ 的点积得到的，对于正交匹配跟踪算法，重构信号 $\check{\mathbf{X}}$ 是通过计算最小二乘解 $\mathbf{w}={\mathbf{D}}_{new}^{+}\mathbf{A}$ 得到的，这个过程需要花费一些时间，所以OMP算法比MP算法要慢一些。
• 残差 $\mathbf{f}$ 是通过原始的 $\mathbf{A}$ 和 ${\mathbf{D}}_{\text{new}}\cdot \mathbf{w}$ 计算得到的。
• 迭代次数最多是 $\mathbf{D}$ 的列数，或者是信号 $\mathbf{X}$ 的稀疏 $K$ 已知，那么就迭代 $K$ 次。

3 OMP缺陷

当 $\mathbf{D}$ 中存在中两个原子有相关性时，OMP算法可能会得到一个错误的重构信号。

$\mathbf{D}1=\left[\begin{array}{lll}0.6& 0.8& 1\\ 0.8& 0.6& 0\end{array}\right]$ and $\mathbf{D}2=\left[\begin{array}{lll}0.6& 0.61& 1\\ 0.8& 0.79& 0\end{array}\right]$

矩阵 $\mathbf{D}2$ 具有非常高的一致性，因为第二列与第一列非常相似。矩阵 $\mathbf{D}1$ 的一致性较差，因为第二列与第一列和第三列不太相似。由 $\mu$ 表示的一致性值定义为
$$\mu =\underset{i,j;i\ne j}{\max }|⟨\mathbf{D}(:,i),(\mathbf{D}:,j)⟩|.$$
$\mu$ 的值介于 $0$ 和 $1$ 之间。如果 $\mu$ 非常高，则OMP通常会给出错误的重建结果。

4 说明总结

在正交匹配追踪中，残差总是与已经选择过的原子所形成的的线性表示向量正交的。这意味着一个原子不会被选择两次，结果会在有限的几步收敛。

LAST、参考文献

匹配追踪算法（MP）简介 - 子孑 - 博客园

正交匹配追踪_给永远比拿愉快-CSDN博客

压缩感知重构算法之正交匹配追踪(OMP)_彬彬有礼的专栏-CSDN博客_omp

MP算法和OMP算法及其思想_逍遥剑客的专栏-CSDN博客_omp算法

压缩感知重构算法之压缩采样匹配追踪(CoSaMP)_彬彬有礼的专栏-CSDN博客

稀疏分解中的MP与OMP算法 - 云+社区 - 腾讯云