正交匹配追踪算法(OMP)简介与详解

秒懂

匹配追踪算法(MP)

为了了解OMP我们先从MP开始说起。

匹配追踪最早是时频分析的分析工具，目的是要将一已知讯号拆解成由许多被称作为原子讯号的加权总和，而且企图找到与原来讯号最接近的解。

在机器学习领域，我们可以将它抽象为$X\omega =Y$，其中$X$是我们已知的向量，$Y$是我们要去拟合或者称为追踪的向量， MP算法的目的就是找到公式中的$\omega$。MP算法采用逐步迭代的方式，每一步寻找到$X$中与$Y$内积最大的向量，将其从$Y$中剔除（相减）得到$Y^{\prime }$（被称为残差），这样一步步迭代下去就能找到较优的近似解。

举个例子，我想要使用$x_1=(-1,1)$和$x_2=(1,0)$两个向量来拟合向量$y=(1,1)$，如果让我们口算，我们可以很容易地计算出$1x_1+2x_2=y$，那么使用MP算法呢？

• 第一步迭代：
  计算$x_1$和$x_2$与$y$的内积$m_1=<x_1,y>=0$，$m_2=<x_2,y>=1$
  比较绝对值发现$x_2$对$y$的影响更大，因此$y^{\prime }=y-m_2{x}_2=(0,1)$，记录$\omega =(\_,1)$。

• 第二步迭代：
  计算$x_1$和$x_2$与$y^{\prime }$的内积$m_1=<x_1,y^{\prime }>=1$，$m_2=<x_2,y^{\prime }>=0$
  比较绝对值发现$x_1$对$y^{\prime }$的影响更大，因此$y^{\prime \prime }=y^{\prime }-m_1{x}_1=(1,0)$，记录$\omega =(1,1)$。

• 第三步迭代：
  计算$x_1$和$x_2$与$y^{\prime }$的内积$m_1=<x_1,y^{\prime \prime }>=-1$，$m_2=<x_2,y^{\prime \prime }>=1$
  比较绝对值发现$x_1$与$x_2$对$y^{\prime }$的影响一样大，因此随便选一个$y^{\prime \prime \prime }=y^{\prime \prime }-m_2{x}_2=(0,0)$，记录$\omega =(1,2)$，得到正确答案。

but！如果随机选选到了$x_1$怎么办，这就可能导致出现迭代时间变长、收敛速度变慢甚至死循环等问题。于是引出正交匹配追踪算法。

正交匹配追踪算法

正交匹配追踪算法是匹配追踪算法的升级版，每一步迭代时保证剔除后的$Y^{\prime }$与被剔除的向量正交，这样之后的迭代就不需要再计算已经剔除掉的向量，极大地减少了时间开销。
例如在上一节的问题中，使用$x_1=(-1,1)$和$x_2=(1,0)$两个向量来拟合向量$y=(1,1)$

• 第一步迭代：
  计算$x_1$和$x_2$与$y$的内积$m_1=<x_1,y>=0$，$m_2=<x_2,y>=1$
  比较绝对值发现$x_2$对$y$的影响更大，因此$y^{\prime }=y-k_1{x}_2=(1-k_1,1)s.t.y^{\prime }\perp x_2$，解得$k_1=2,y^{\prime }=(-1,1)$，记录$\omega =(\_,2)$。

• 第二步迭代：
  由于正交性，不需要再考虑$x_2$,因此$y^{\prime \prime }=y^{\prime }-k_2{x}_1=(-1-k_2,1-k_2)s.t.y^{\prime }\perp x_1\wedge y^{\prime }\perp x_2$，解得$k_2=1$记录$\omega =(1,1)$，得到正确答案。

正交匹配追踪算法的优势可见一斑。

详解

算法流程

简单介绍已经说清楚了，即通过$X\omega =Y$中已知的$X$和$Y$来求出$\omega$，我们将其拆开来写
$${\omega }_1{x}_1+{\omega }_2{x}_2+...+{\omega }_n{x}_n=y$$

1. 第一步计算每个向量$x$与$y$的内积，挑选出绝对值最大的设为$x_i$，令$y^{\prime }=y-{\omega }_i{x}_i$，${\omega }_i$满足$y^{\prime }\perp x_i$，将$y$替换成$y^{\prime }$（即令$y=y^{\prime }$）。
2. 随后，计算除$x_i$之外每个向量$x$与$y$的内积，挑选出绝对值最大的设为$x_j$，令$y^{\prime }=y-{\omega }_j{x}_j$，${\omega }_j$满足$y^{\prime }\perp x_j\wedge y^{\prime }\perp x_i$，将$y$替换成$y^{\prime }$。
3. 重复步骤2，值得注意的是，这里的政教关系需要用施密特正交化方法来计算。
4. 求出的$\omega$就是OMP算法的解。

代码实现

废话不多说，直接上sklearn

from sklearn import linear_model


def data_generate():
	# 这里随机生成了一些数据
    a = [[i + 1, i - 2, i + 4, i - 3] for i in range(5)]
    a.extend([[i - 1, i - 2, i, i] for i in range(5)])
    b = [sum(i) * 2 + 1 for i in a]
    # b = [[i + random.randint(-10, 10) for _ in range(1)] for i in range(5)]
    # b.extend([[i + random.randint(-10, 10) for _ in range(1)] for i in range(10, 15)])
    return a, b


if __name__ == "__main__":
    x, y = data_generate()
    reg = linear_model.OrthogonalMatchingPursuit()
    # 可以和下面这一行的标准线性回归做一个对比
    # reg = linear_model.LinearRegression()
    reg.fit(x, y)
    for i, j in zip(x, y):
        print(i, j)
    print('---------------')
    # 懒得再构造数据了，使用训练集预测一遍，效果。。。一般，主要是因为数据量太少了。
    # 当数据量增加时，OMP可以超过标准线性回归的效果
    print(reg.predict(x))
    print(reg.score(x, y))
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26

算法优劣

通过上文中的介绍可以发现
OMP算法相较于MP算法收敛快，效率高。
作为压缩感知领域的算法，可以应用在回归问题上，效果还不错。