施密特(Schimidt)正交化与正交匹配追踪

题目：施密特(Schimidt)正交化与正交匹配追踪

        文献[1]中给出了施密特(Schimidt)正交化的过程：

上面的的[x,y]表示向量内积，[x,y]=x^Ty=y^Tx=[x,y]。施密特正交化公式中的b_r实际上可写为：

分子之所以可以这么变化是由于[x,y]实际上为一个数，因此[x,y]x=x[x,y]= xx^Ty。

        上一篇《稀疏表示与匹配追踪》中详细的解释了匹配追踪（Matching Pursuit，MP）的流程，在最后给出了正交匹配追踪（Orthogonal Matching Pursuit，OMP）的流程，并指出OMP与MP的不同根本在于残差更新过程：OMP减去的Pe_m是在所有被选择过的原子组成的矩阵φ_t所张成空间上的正交投影，而MP减去的Pe_m是在本次被选择的原子φ_p所张成空间上的正交投影。正交投影及正交投影变换矩阵的概念可以参见《压缩感知中的数学知识：投影矩阵（projectionmatrix）》。

        首先给出一个结论：

        设OMP共从冗余字典中选择了r个原子，分别是a₁，a₂，……，a_r，根据正交匹配追踪的流程可以知道待分解信号x最后剩余的残差e_romp为

         (式1)

        该残差也可以表示为

         (式2)

其中矩阵A为选择的r个原子组成的矩阵，e_(r-1)omp为选择(r-1)个原子时的残差。

        将选择的r个原子a₁，a₂，……，a_r进行施密特正交化：

则残差e_romp还可以写为

         (式3)

(式1)一般出现在稀疏分解算法中，(式2)一般出现在重构算法中，(式3)是自己琢磨出来的（受到沙威的文档中提到的施密特正交化的启发，但沙威只限于向量情况下，详情可参见[3]，此处相当于一个推广）。

        而实际上(式1)、(式2)和(式3)三种正交匹配追踪算法的残差计算方法是等价的，纯数学证明比较复杂，这里给一个MATLAB程序，可以验证三者的等价性：

%% 验证OMP残差求解过程与Schmidt正交化的关系
M = 4;N = 10;
Phi = randn(M,N);
for nn = 1:N
    Phi(:,nn) = Phi(:,nn)/norm(Phi(:,nn));
end
b = randn(M,1);
e0 = b;%初始化残差为待稀疏信号b
%选第1列
c1 = Phi'*e0;%求出矩阵Phi每列与b的内积
[val1,pos1]=max(abs(c1));%找到内积中最大的列及其内积值
phit = [Phi(:,pos1)];%由所有选出的列组合的矩阵
Pphi = phit*(phit'*phit)^(-1)*phit';%正交投影变换矩阵
e1ompe = e0 - Pphi*e0;%OMP用上一次残差减去残差在phit列空间的正交投影
e1ompb = b - Pphi*b;%OMP用待稀疏信号b减去b在phit列空间的正交投影
x = Phi(:,pos1);%Schimidt正交化第一个向量
Px = x*(x'*x)^(-1)*x';
%实际上是b - Px*b
e1ompsmt = e0 - Px*b;
e1 = e1ompe;
norm(e1ompe-e1ompb)+norm(e1ompsmt-e1ompb)
%选第2列
c2 = Phi'*e1;%求出矩阵Phi每列与e1的内积
[val2,pos2]=max(abs(c2));%找到内积中最大的列及其内积值
phit = [Phi(:,pos1) Phi(:,pos2)];%由所有选出的列组合的矩阵
Pphi = phit*(phit'*phit)^(-1)*phit';%正交投影变换矩阵
e2ompe = e1 - Pphi*e1;%OMP用上一次残差减去残差在phit列空间的正交投影
e2ompb = b - Pphi*b;%OMP用待稀疏信号b减去b在phit列空间的正交投影
y = Phi(:,pos2) - Px*Phi(:,pos2);%Schimidt正交化第二个向量
Py = y*(y'*y)^(-1)*y';
%实际上是b - Px*b - Py*b
e2ompsmt = e1 - Py*b;%上一次残差减去b在第2列正交化所得z上的投影
e2 = e2ompe;
norm(e2ompe-e2ompb)+norm(e2ompsmt-e2ompb)
%选第3列
c3 = Phi'*e2;%求出矩阵Phi每列与e2的内积
[val3,pos3]=max(abs(c3));%找到内积中最大的列及其内积值
phit = [Phi(:,pos1) Phi(:,pos2) Phi(:,pos3)];%由所有选出的列组合的矩阵
Pphi = phit*(phit'*phit)^(-1)*phit';%正交投影变换矩阵
e3ompe = e2 - Pphi*e2;%OMP用上一次残差减去残差在phit列空间的正交投影
e3ompb = b - Pphi*b;%OMP用待稀疏信号b减去b在phit列空间的正交投影
z = Phi(:,pos3) - Px*Phi(:,pos3) - Py*Phi(:,pos3);%Schimidt正交化第三个向量
Pz = z*(z'*z)^(-1)*z';
%实际上是b - Px*b - Py*b - Pz*b
e3ompsmt = e2 - Pz*b;%上一次残差减去b在第3列正交化所得z上的投影
e3 = e3ompe;
norm(e3ompe-e3ompb)+norm(e3ompsmt-e3ompb)

程序最后一行的输出结果为零，说明三个残差是相等的。

        这里可以看出OMP与Schimidt正交化的关系：

        OMP分解过程，实际上是将所选原子依次进行Schimidt正交化，然后将待分解信号减去在正交化后的原子上各自的分量即可得残差。其实(式3)求残差的过程也是在进行施密特正交化。

        有个关键问题还是要说的，分解后在所选择各原子上的系数是多少呢？答案其实也很简单，各个系数是(A^TA)^-1A^Tx，即最小二乘解，这个解是一个列向量，每一个元素分别是组成矩阵A的各原子的线性组合系数，这个在下一篇《正交匹配追踪(OMP)在稀疏分解与压缩感知重构中的异同》也会明确再次说明。

        同理，若设MP共从冗余字典中选择了r个原子，分别是a₁，a₂，……，a_r，根据匹配追踪的流程可以知道待分解信号x每次迭代后剩余的残差e_rmp为

比较式(3)的第2个等号表示的e_romp与此处的e_rmp也可以体会出OMP与MP的区别吧。

【参考文献】

【1】同济大学数学系. 线性代数（第五版）[M].高等教育出版社，2007：114.

【2】Joel A. Tropp and AnnaC. Gilbert . Signal Recovery From Random Measurements Via Orthogonal MatchingPursuit[J]. IEEE Transactions on Information Theory, VOL. 53, NO. 12, DECEMBER2007.

【3】沙威. “压缩传感”引论.http://www.eee.hku.hk/~wsha/Freecode/Files/Compressive_Sensing.pdf