压缩感知——浅析匹配追踪算法（Matching Pursuit，MP）与正交匹配追踪算法（Orthogonal Matching Pursuit，OMP）

前言

考虑下列数学模型：

在压缩感知（Compressed sensing）术语中，从X和A中找到Y称为压缩（compression）。 相反，从Y和A中找到原始X称为重构（reconstruction）。
x被称为原始信号（original signal）或原始向量（original vector），A被称为压缩矩阵（Compression matrix）或感知矩阵(sensing matrix)，而Y被称为压缩信号(compressed signal)或压缩向量(compressed vector)。

第一个例子

给定压缩矩阵A和原始向量X：

很容易就可以得到压缩向量Y：

然而，反过来给定压缩向量Y和原始向量X：

很难得到找到原始的(或尽可能接近的) X。

基本概念

首先介绍一些基本概念。

我们将压缩矩阵A看作是列向量的集合，即：

其中：

这些列向量(b₁,b₂,b₃)被称为基(basis)或原子(atoms)。
令：

那么：

上述方程告诉我们，Y是用X中的系数(x₁, x₂, x₃)表示的基(b₁,b₂,b₃)的线性组合，且已知x₁=-1.2，x₂=1，x₃=0。换句话说，在上述方程中，基b₁对生成Y的贡献分别是-1.2，基b₂对生成Y的贡献分别是1，基b₃对生成Y的贡献分别是0。

OMP算法求解过程

OMP的工作原理与MP相似，但MP是找出哪一个基对Y贡献最大，然后哪一个基贡献第二，哪一个基贡献第三，依此类推。OMP计算过程需要进行N次迭代，其中N是A中基的个数。在上述的示例中，N为3。

1. 计算基的贡献度
上述例子中有三个基，分别是：

以及

所以Y中每个基的贡献度分别是：

显然，基b₁的贡献最大，其贡献度为-1.34（使用绝对值）。
在第一轮迭代中，选择基b₁，该基的贡献度为-1.34。
当然，可以在一个单独的步骤中计算每一个基的点积，即:

2. 计算残差r₁
首先，选择第一个基，即基b₁和贡献度λ₁=-1.34。残差为：

如图2所示，残差r₁与第一个基b₁是垂直的。

3. 重复迭代
第一次迭代后，我们得到所选基b₁。把这些基重新组合成一个新的基A_new。即：

重写贡献度值，即：

值-1.34在第一个元素因为它是第一个基b₁的贡献。
第一次迭代后残差为：

现在我们要从剩下的基b₂和b₃中选择对残差r贡献最大的基，在一步中：

因此，基b₂的贡献更大，即0.98。
我们把选定的基b₂添加到矩阵A_new中，

现在，这一步不同于MP。这里，我们需要重新计算基b₁和b₂共同对Y的贡献(而不是基b₂对r₁的贡献)

为了解决这个贡献问题，我们必须求出最小值平方问题，很容易表述为:

在数学术语中，这个公式可以写成:

根据最小二乘原理，上述公式的解为：

其中，A_new⁺是A_new的广义逆矩阵。

在上述例子中

在Matlab中计算广义逆矩阵很简单，只需使用pinv()命令即可。

计算λ之后，更新x_rec

在这里我们把λ放在x_rec在第一和第二行,因为它对应于第一和第二被选择的基b₁和b₂。
重新得到A_new和λ之后,现在我们可以更新新的残差r₂

到此，残差为0。我们可以停在这里，或者继续下一步。(如果在这里停止，我们可以节省一些计算工作) 。

4. 最后迭代

由于残差已经消失，因此不需要此步骤。（许多OMP实现软件需要信号稀疏度K的参数。这是该软件将迭代K次的线索。此后它将停止，无论剩下多少残差）
因此，我们的重构信号为:

它和原始信号是一样的.

以上便是OMP的整个计算流程。

MP算法求解过程

1. 计算贡献度
Y中每个基的贡献度分别是：

显然，基b₁的贡献最大，其贡献度为-1.34（使用绝对值）。

==================================================================================

第一轮迭代

选择基b₁，该基的贡献度为-1.34
重写贡献度值，即：

计算残差

接下来我们计算这个残差和基b₂和基b₃的点积。(不需要用基b₁计算，因为这个残差必定垂直于基b₁).

取绝对值，我们得到基b₂是第二强的影响。
贡献度值，即：

第一轮迭代完成

==================================================================================

第二轮迭代，计算残差：

接下来我们计算这个残差r₂和基b₁基b₃的点积。

在第二轮迭代中选择基b₁。

贡献度值，即：

第二轮迭代完成

==================================================================================

第三轮迭代，计算残差：

接下来我们计算这个残差r₃和基b₂基b₃的点积。

在第三轮迭代中选择基b₂。

贡献度值，即：

第三轮迭代完成

==================================================================================

第四轮迭代，计算残差：

重复上述过程，直到残差达到指定阈值。

==================================================================================
以上便是MP算法的计算流程。

第二个例子

现在我们继续看另一个例子，给定：

很容易就可以得到：

现在给定A和Y用OMP计算X。

步骤如下：

1.A中有四个基

2. 因为基的长度不是1，我们需要把它们都归一化。

3. 计算这些归一化基的贡献度:

4.第二个标准化的基b ̂₂的贡献最高。因此，我们将第二个基b₂放入矩阵A_new中。

5.计算最小二乘问题的解λ₁

因为这个λ₁是基b₂的系数，那么:

6.计算残差

7.重复第3步。计算基b ̂₁，b ̂₃，b ̂₄中哪一个基对残差r₁的贡献最大。

在这里，在基b ̂₃的贡献最大的，是1.7902。
8.将基b₃而不是(b ̂₃）添加到A_new中。

4. 计算最小二乘问题的解λ₂

由于λ₂对应基b₂和基b₃，因此：

10.计算残差
11.重复第3步，直到最后一次迭代。
12.计基b₁和基b₄的贡献度

因此，基b₄的贡献度最大，为0.7308.
5. 将基b₄而不是(b ̂₄）添加到A_new中。

6. 计算最小二乘问题的解λ₃.
因为λ₂对应基b₂、基b₃和基b₄，因此，

15.计算残差

16.迭代到此停止，因为残差已经全部是0了。重构的原始向量x_rec为：

MP和OMP的缺点

OMP很快。与凸优化相比，速度非常快。这是贪婪算法的一个优点。但是OMP也遭受矩阵A高度相干性的困扰。矩阵中的相干性是什么？简单来说，矩阵A中的相干性表示矩阵A中两列之间的相似性(similarity)。
请看下面两个矩阵：

矩阵A₂具有非常高的相干性，因为第二类和第一列非常相似。矩阵的相干性较低，因为第二列与第一列和第三列不太相似。相干值用μ表示：

μ的值在0和1之间。如果μ的值很大，则通常OMP给出一个错误的重构结果。

试试下面的例子：
已知

求出Y。
然后，给定A和Y，使用MP和OMP求X。
下面直接给出结果：

可以看出，MP和OMP有类似特性，且都不能计算出正确结果

后记

正如我们在上面的过程中看到的那样，我们现在注意到以下事实。
（1）我们从OMP计算得出的最高贡献来自标准化基础。并非来自最初的非标准化基础。
（2）如果测量或重构矩阵A中的所有基础都已被归一化，则无需进行归一化步骤即可继续进行计算。
（3）每个基的贡献是由残基和归一化基的点积完成的。
（4）在MP的情况下，重建系数x_rec是根据基和残基的点积计算得出的，而在OMP上，则是根据A_new到Y的最小二乘解计算得出的系数x_rec。最小二乘解λ表示基对应位置的x_rec。A_new本身是从非标准化的基得来的，这个过程需要时间。因此，OMP比MP更慢。
（5）残差r是由Y、A_new和λ计算得来的
（6）迭代将持续与A中的行数（即M）一样多。或者，如果信号k的稀疏度已知，则迭代重复k次。如果k <M，则知道k将有助于我们更快地进行计算。但是，如果k未知，则迭代次数最多应为M。

参考文献

[1]: Usman, Koredianto, (2017), Introduction to Orthogonal Matching Pursuit, Telkom University, Online : http://korediantousman.staff.telkomuniversity.ac.id.
[2]: Usman, Koredianto, 2017, Introduction to Matching Pursuit, Telkom University, Online : http://korediantousman.staff.telkomuniversity.ac.id/tutorial.