OAMP的理解

1 前言

本次博文主要介绍了OAMP论文,同时加了一些粗浅的理解。一方面,前两部分涉及了一些AMP相关的知识以及我自己给出的解释,另一方面,博文所推公式可能与论文稍有偏差,笔者才疏学浅,又是第一次写博客,这两个部分都不敢保证能理解到位和描述清晰,难免解释可能会有所偏差,甚至有歪曲原文的不当之处,还希望各位读者能够包涵。下一篇博文可能会是对VAMP的解读,希望这可以在年内完成,之后博客可能会推出一些与消息传递算法较相关的论文介绍,以及自认为有价值的积累。

2 绪论

简短回顾一下Approximate Message Passing (AMP) 的问题模型

y

=

A

x

+

n

(1a)

pmb y=pmb A pmb x + pmb n tag{1a}

yyy=AAAxxx+nnn(1a)

x

j

P

X

(

x

)

,

x

(1b)

mathit x_{j} sim P_{X}(x), forall x tag{1b}

xjPX(x),x(1b)

其中

A

C

M

×

N

pmb A in mathbb C^{M times N}

AAACM×N是感知矩阵,

n

C

M

×

1

pmb n in mathbb C^{M times 1}

nnnCM×1是均值为

0

0

0,方差为

σ

2

sigma^{2}

σ2的高斯向量。AMP的一个重要性质就是其算法性能可以用state evolution来精确刻画,就是说给定真实的

x

pmb x

xxx,我们可以依据state evolution的结果

τ

t

tau^{t}

τt,预测AMP在第

t

t

t次迭代过程中估计得到的

x

^

t

{hat {pmb x}}^{t}

xxx^t与真实

x

pmb x

xxx的均方差(MSE)

E

[

x

x

^

t

2

2

]

mathbb E[{Vert pmb x-{hat {pmb x}}^{t} Vert}^2_{2} ]

E[xxxxxx^t22],表示为

τ

t

2

1

N

x

x

^

t

2

2

,

(

N

)

(2)

{ tau_{t} }^{2} rightarrow frac{1}{N} {Vert pmb x-{hat {pmb x}}^{t} Vert}^2_{2}, (N rightarrow infty) tag{2}

τt2N1xxxxxx^t22,(N)(2)

事实上,式(2)为趋近于而不是严格等于,因为

1

N

x

x

^

t

2

2

frac{1}{N} {Vert pmb x-{hat {pmb x}}^{t} Vert}^2_{2}

N1xxxxxx^t22是关于

x

^

t

{hat {pmb x}}^{t}

xxx^t的二阶Lipschitz函数,满足如下收敛关系

1

N

x

x

^

t

2

2

E

[

x

x

^

t

2

2

]

,

(

N

)

(3)

frac{1}{N} {Vert pmb x-{hat {pmb x}}^{t} Vert}^2_{2} rightarrow mathbb E[{Vert pmb x-{hat {pmb x}}^{t} Vert}^2_{2} ] , (N rightarrow infty) tag{3}

N1xxxxxx^t22E[xxxxxx^t22],(N)(3)

但遗憾的是,大部分情况下,只有当

A

pmb A

AAA为高斯矩阵或者次高斯矩阵时,state evolution才能与AMP估计的结果统一,如果感知矩阵的特征值分布与高斯矩阵的特征值分布相差较远,AMP的性能就不能保证,甚至可能会出现不收敛的情况。为了解决该问题,Junjie Ma和Li Ping提出了OAMP[1]

AMP还有一个重要的点是其线性迭代过程中含有"onsager"这一项,它的作用是为了消除迭代过程中感知矩阵

A

pmb A

AAA与估计结果

x

^

t

{hat {pmb x}}^{t}

xxx^t之间的相关性。虽然OAMP把AMP中的“onsager”这一项给去掉了,但是为了补偿"onsager"原来的作用,OAMP在非线性估计中加入了divergence-free的约束(可能divergence-free这个概念有点抽象,简单理解为导数为0就行)。

备注:其实OAMP并没有严格意义的数学推导,作者先是给了两个独立性的假设(假设1和假设2),而OAMP-state evolution就是由该假设条件推出来的。让OAMP迭代开始(

t

=

0

t=0

t=0时刻)先保证独立性,满足两个假设条件之一,如果之后的每一次迭代假设1和假设2都能互相推出对方,那么我们可以认为每一次迭代的程始终可以保证独立条件成立,也就保证了state evolution。但略有遗憾,假设1和假设2只能“部分”地互相推出对方,可以保证不相关(由正交推出),但是不能保证独立。幸运的是,仿真结果表明,虽然迭代过程中的独立条件不能始终保持,但是state evolution还是能够和OAMP估计结果统一,论文作者将此描述为:假设1和假设2是只是state evolution的充分条件。

3 AMP

在开始OAMP之前,先回顾一下AMP。

3.1 AMP算法

假设矩阵

A

=

[

a

1

,

a

2

,

,

a

N

]

pmb A=[pmb a_{1}, pmb a_{2}, text{…},pmb a_{N}]

AAA=[aaa1,aaa2,,aaaN]是列归一化的,即,

i

{

1

,

...

,

N

}

,

a

i

C

M

×

1

{forall i} in {1,text{...} ,N}, pmb a_{i} in mathbb C^{M times 1}

i{1,...,N},aaaiCM×1,满足

E

[

a

i

2

]

=

1

mathbb E[{Vert pmb a_{i} Vert}_2]=1

E[aaai2]=1。AMP的迭代过程如下

r

t

=

s

t

+

A

T

(

y

A

s

t

)

+

N

M

<

η

t

1

(

r

t

1

)

>

(

r

t

1

s

t

1

)

o

n

s

a

g

e

r

 

t

e

r

m

(4a)

pmb r^{t}=pmb s^{t} + pmb A^{T}(pmb y - pmb A pmb s^{t}) +\ underbrace {frac {N}{M} <{eta_{t-1}}^{'}( pmb r^{t-1})>( pmb r^{t-1}-pmb s^{t-1} )}_{onsager text{ } term} tag{4a}

rrrt=ssst+AAAT(yyyAAAssst)+onsager term

MN<ηt1(rrrt1)>(rrrt1ssst1)(4a)

s

t

+

1

=

η

t

(

r

t

)

(4b)

pmb {s^{t+1}}={eta}_{t}(pmb r^t) tag{4b}

st+1st+1st+1=ηt(rrrt)(4b)

其中

η

t

{eta}_{t}

ηt是关于

r

t

pmb r^t

rrrt的一个Lipschitz连续函数(component-wise),

s

t

pmb s^t

ssst是最后的估计,

r

t

s

t

pmb r^{t}-pmb s^{t}

rrrtssst其实就是一般AMP表达的“残差”(residual error)。式(4a)的最后一项就是前言所描述的"onsager"项,这一项是根据松弛信念传播算法(relaxed-Belief-Propagation, relaxed-BP)在极限条件下

M

,

N

M,N rightarrow infty

M,N时依据大数定律、中心极限定理和消息近似补偿(为了补偿经过近似的

O

(

1

N

)

O(frac {1}{sqrt N})

O(N

1)的消息),以及泰勒展开逐步推导出来的。

"onsager"项隐含的意义:"onsager"项的存在确保了AMP-state evolution的正确性,但是在推导AMP之前,relaxed-BP里边的方差项其实跟state-evolution是对应的(只是对应,并不相等),区别在于当

M

,

N

<

M,N lt infty

M,N<时,relaxed-BP的方差项(指消息传递过程中,消息分布的方差逐渐积累)一般记为

Σ

m

n

(

t

)

Sigma^n_{m}(t)

Σmn(t),这里不去管具体的

n

,

m

n,m

n,m是什么含义了,

t

t

t是指迭代到第

t

t

t步,只是它隐式地强调了

A

pmb A

AAA和估计

s

t

pmb s^{t}

ssst具有相关性。而

M

,

N

M,N rightarrow infty

M,N时借助大数定理和中心极限定理进一步推导可以使得

lim

M

,

N

Σ

m

n

(

t

)

Σ

(

t

)

lim_{M,Nto infty}Sigma^n_{m}(t) rightarrow Sigma(t)

limM,NΣmn(t)Σ(t),即去掉了相关性。

3.2 AMP-state evolution与等效信号模型

为了方便之后OAMP的讨论,下面的叙述几乎都是基于OAMP原论文的,表达式跟一般的AMP表达式有些差异,是为了类比,在后面方便证明正交性,本质上并无差异。

定义两种误差,第一种误差是非线性估计结果

s

t

pmb s^t

ssst与真实值的差,第二种误差是线性估计结果

r

t

pmb r^t

rrrt与真实值的差,分别定义为

q

t

=

s

t

x

(5a)

pmb q^t = pmb s^t - pmb x tag{5a}

qqqt=ssstxxx(5a)

h

t

=

r

t

x

(5b)

pmb h^t = pmb r^t - pmb x tag{5b}

hhht=rrrtxxx(5b)

结合(5), (4)可以被展开为

h

t

=

(

I

A

T

A

)

q

t

+

A

T

n

+

N

M

<

η

t

1

(

x

+

h

t

1

)

>

(

h

t

1

q

t

1

)

(6a)

pmb h^t = (pmb I - pmb A^T pmb A)pmb q^t + pmb A^Tpmb n + \frac {N}{M}<eta^{'}_{t-1}(pmb x + pmb h^{t-1})>(pmb h^{t-1}-pmb q^{t-1}) tag{6a}

hhht=(IIIAAATAAA)qqqt+AAATnnn+MN<ηt1(xxx+hhht1)>(hhht1qqqt1)(6a)

q

t

+

1

=

η

t

(

x

+

h

t

)

x

(6b)

pmb q^{t+1}=eta_{t}(pmb x + pmb h^t)-pmb x tag{6b}

qqqt+1=ηt(xxx+hhht)xxx(6b)

式(6)并不是要给迭代算法,因为真实值

x

pmb x

xxx作为一个已知的量参与其中。

与AMP对应的evolution由下式直接给出,这里不展开描述,因为AMP-evolution的证明过程极其复杂,但是如果只是理解,有一个简单的思路,就是如上面所述,从relaxed-BP入手,着眼于其方差演变,可以在一定程度上帮助理解。

τ

t

2

=

N

M

v

t

2

+

σ

2

(7a)

tau^2_{t}=frac {N}{M} v^2_{t}+sigma^2 tag{7a}

τt2=MNvt2+σ2(7a)

v

t

+

1

2

=

E

{

[

η

t

(

X

+

τ

t

Z

)

X

]

2

}

(7b)

v^2_{t+1}=mathbb E{[eta_{t}(X+tau_{t}Z)-X]^2} tag{7b}

vt+12=E{[ηt(X+τtZ)X]2}(7b)

注意式(7b)中的

X

,

Z

X,Z

X,Z表示随机变量,并且

Z

N

(

0

,

1

)

Z sim mathcal N(0,1)

ZN(0,1)

X

X

X独立。这里的理解对AMP以及AMP-state evolution是至关重要的,因为其独立性,我们可以将每次迭代估计得到的结果模型等效为

X

^

t

=

X

+

τ

t

2

Z

,

Z

N

(

0

,

1

)

(8)

hat X^t=X+tau^2_{t}Z, Z sim mathcal N(0,1) tag{8}

X^t=X+τt2Z,ZN(0,1)(8)

这也意味着

X

^

t

N

(

X

,

τ

t

2

)

(9)

hat X^t sim mathcal N(X,tau^2_{t}) tag{9}

X^tN(X,τt2)(9)

τ

t

2

tau^2_{t}

τt2指的是估计后的方差,结合式(7a)可以看出,方差项由高斯噪声的方差

σ

2

sigma^2

σ2和非线性估计结果与真实值之间的MSE构成。这反映了两点,一方面,高斯噪声与其他各个变量都保持独立,另一方面,AMP从始至终虽然都旨在接近真实值,因为

τ

t

tau^t

τt一直在变小,但是依然没有克服高斯噪声。

等效模型还有一个重要的作用就是,在部分

η

eta

η函数的推导过程中,使用等效模型的概念可以一定程度上简化推导过程,比如推导MMSE函数,或者退化为LMMSE等。以及,

τ

t

2

tau^2_{t}

τt2可能作为其中的一个参数,然而实际的仿真中我们并不知道真实

τ

t

2

tau^2_{t}

τt2,这个时候就需要将其近似为

τ

^

t

2

hat tau^2_{t}

τ^t2,具体表示为

τ

^

t

2

=

1

N

r

t

s

t

2

2

(10)

hat tau^2_{t}=frac {1}{N} {Vert pmb r^t - pmb s^t Vert}^2_{2} tag {10}

τ^t2=N1rrrtssst22(10)

4 OAMP

4.1 OAMP产生的动机

在阐述动机之前,先描述论文给出的一个定义

定义1: Divergence-free
对函数

η

:

R

R

eta: mathbb R mapsto mathbb R

η:RR,如果满足

E

[

η

(

R

)

]

=

0

mathbb E [eta^{'}(R)]=0

E[η(R)]=0

那么就认为

η

eta

η是devergence-free

根据定义1,一个divergence-free的函数可以被构造成

η

(

r

)

=

C

(

η

^

(

r

)

E

[

η

^

(

R

)

]

r

)

(11)

eta(r)=C cdot (hat eta( r)-mathbb E[hat eta(R)] cdot r) tag{11}

η(r)=C(η^(r)E[η^(R)]r)(11)

其中

η

^

hat eta

η^是任意一个一阶可导的函数,C是常数。

如果把AMP迭代公式(4)中的"onsager"项给去掉,

η

eta

η函数按照式(11)的方式给出,其中

η

^

hat eta

η^是软阈值函数(压缩感知常用来恢复稀疏信号的函数),在这样的设置下,作者发现即使感知矩阵

A

pmb A

AAA不是高斯矩阵,而是一个离散余弦变换矩阵,state evolution的结果却意外地与去"onsager"的AMP迭代结果一致。这个发现也就引出了OAMP。

4.2 去相关的线性估计

先给出两个定义

定义2: 酉不变矩阵(Unitarily-Invariant Matrix)
如果矩阵

U

,

V

,

Σ

pmb U,pmb V,pmb Sigma

UUU,VVV,ΣΣΣ三者之间相互独立,并且

U

,

V

pmb U,pmb V

UUU,VVV满足Haar分布(随机各向同性)是正交阵,

Σ

pmb Sigma

ΣΣΣ是对角阵,那么认为

A

=

U

Σ

V

pmb A= pmb U pmb Sigma pmb V

AAA=UUUΣΣΣVVV是酉不变的。

定义3: 去相关矩阵
如果感知矩阵

A

=

U

Σ

V

pmb A= pmb U pmb Sigma pmb V

AAA=UUUΣΣΣVVV是酉不变的,矩阵

W

pmb W

WWW如果满足

t

r

(

I

W

A

)

=

0

tr(pmb I - pmb W pmb A)=0

tr(IIIWWWAAA)=0,就说

W

pmb W

WWW是关于

A

pmb A

AAA的一个去相关矩阵。指定准去相关矩阵

W

^

=

U

G

V

T

(12)

hat {pmb W}=pmb U pmb G pmb V^T tag{12}

WWW^=UUUGGGVVVT(12)

那么满足

t

r

(

I

W

A

)

tr(pmb I - pmb W pmb A)

tr(IIIWWWAAA)的去相关矩阵

W

pmb W

WWW可以被构建为

W

=

N

t

r

(

W

^

A

)

W

^

(13)

pmb W= frac {N}{tr(hat {pmb W} pmb A)} hat {pmb W} tag{13}

WWW=tr(WWW^AAA)NWWW^(13)

其实定义3在论文中并没有直接给出,是我抽出来的,而且跟论文稍有出入,但是问题不大。我刚开始读到这里的去相关概念的时候非常不理解,如果矩阵A的映射是单射或者双射,可能稍微懂个大概,满射可能把

W

pmb W

WWW乘在右边,但是概念依然很模糊,后来请教了一下数院的学长,大概意思就说作者就是这么叫了而已。所以这里也就理解个大概吧,没有再细究。

下面给出三个常用的准去相关矩阵
(1)匹配滤波器(Matched Filter, MF)

W

^

M

F

=

A

T

(14a)

hat {pmb W}^{MF}=pmb A^T tag{14a}

WWW^MF=AAAT(14a)

(2)伪逆(Pseudo-inverse, PINV)

W

^

P

I

N

V

=

{

A

T

(

A

A

T

)

1

;

M

<

N

(

A

T

A

)

1

A

T

;

M

>

N

(14b)

hat {pmb W}^{PINV}= left{ begin{array}{lr} pmb A^T (pmb A pmb A^T)^{-1}; M<N \ (pmb A^T pmb A)^{-1} pmb A^T; M>N end{array} right tag{14b}.

WWW^PINV={AAAT(AAAAAAT)1;M<N(AAATAAA)1AAAT;M>N(14b)

(3)LMMSE

W

^

L

M

M

S

E

=

A

T

(

A

A

T

+

σ

2

v

2

I

)

1

(14c)

hat {pmb W}^{LMMSE}=pmb A^T(pmb A pmb A^T+frac {sigma^2}{v^2} pmb I)^{-1} tag{14c}

WWW^LMMSE=AAAT(AAAAAAT+v2σ2III)1(14c)

其实式(14c)中

v

2

v^2

v2的含义有些微妙,将在后面合适的地方做更深的阐述。

LMMSE的简述
假设模型为

y

=

A

x

+

n

pmb y = pmb A pmb x + pmb n

yyy=AAAxxx+nnn,假设随机向量

X

,

Y

X,Y

X,Y的均值都为0(实际当中不满足的话可以先减去均值),

n

C

N

(

0

,

σ

2

I

)

pmb n sim mathcal C mathcal N(0,sigma^2 pmb I)

nnnCN(0,σ2III),那么LMMSE的估计为

x

^

=

Σ

x

y

Σ

y

1

y

hat {pmb x} = Sigma_{xy} Sigma^{-1}_{y} pmb y

xxx^=ΣxyΣy1yyy

其中

Σ

x

y

=

E

[

x

y

H

]

=

E

[

x

(

A

x

+

n

)

H

]

=

E

[

x

x

H

A

+

x

n

H

]

=

Σ

x

A

Sigma_{xy} = mathbb E[pmb x pmb y^{H}]=mathbb E[pmb x pmb (pmb A pmb x+pmb n)^{H}]=mathbb E[pmb x pmb x^Hpmb A+pmb xpmb n^H]=Sigma_{x} pmb A

Σxy=E[xxxyyyH]=E[xxx(((AAAxxx+nnn)H]=E[xxxxxxHAAA+xxxnnnH]=ΣxAAA

Σ

y

=

E

[

y

y

H

]

=

E

[

(

A

x

+

n

)

(

A

x

+

n

)

H

]

=

A

Σ

x

A

H

+

σ

2

I

begin{aligned} Sigma_{y} &= mathbb E[pmb y pmb y^{H}] \ &=mathbb E[(pmb A pmb x+pmb n) pmb (pmb A pmb x+pmb n)^{H}] \ &=pmb A Sigma_{x} pmb A^H+sigma^2 pmb I end{aligned}

Σy=E[yyyyyyH]=E[(AAAxxx+nnn)(((AAAxxx+nnn)H]=AAAΣxAAAH+σ2III

那么LMMSE的估计可以转化为

x

^

=

Σ

x

A

(

A

Σ

x

A

H

+

σ

2

I

)

1

y

hat {pmb x} = Sigma_{x} pmb A (pmb A Sigma_{x} pmb A^H+sigma^2 pmb I)^{-1}y

xxx^=ΣxAAA(AAAΣxAAAH+σ2III)1y

该结果还可以继续延申,根据

(

E

+

B

C

D

)

1

=

E

1

E

1

B

(

C

1

+

D

E

1

B

)

1

D

E

1

(pmb E + pmb B pmb C pmb D)^{-1}=pmb E^{-1}- pmb E^{-1} pmb B (pmb C^{-1}+ D pmb E^{-1} pmb B)^{-1} pmb D pmb E^{-1}

(EEE+BBBCCCDDD)1=EEE1EEE1BBB(CCC1+DEEE1BBB)1DDDEEE1

将上述公式代入

Σ

y

Sigma_{y}

Σy项即可。

4.3 OAMP算法

OAMP的迭代公式

r

t

=

s

t

+

W

t

(

y

A

s

t

)

(15a)

pmb r^{t}=pmb s^{t} + pmb W_{t}(pmb y - pmb A pmb s^{t}) tag{15a}

rrrt=ssst+WWWt(yyyAAAssst)(15a)

s

t

+

1

=

η

t

(

r

t

)

(15b)

pmb s^{t+1} = eta_{t}(pmb r^t) tag{15b}

ssst+1=ηt(rrrt)(15b)

其中

W

t

pmb W_{t}

WWWt是去相关矩阵,

η

t

eta_{t}

ηt是满足divergence-free的约束,即式(11)。将该式与AMP迭代式(4)做比较,可以发现线性估计中的矩阵

A

T

pmb A^T

AAAT变得更一般化,不再局限于匹配滤波,而且末尾缺少了"onsager"项,把“onsager”的作用加在了divergence-free约束里,这也跟OAMP动机部分所阐述的内容一致。

4.4 估计误差迭代与OAMP-state evolution

依然保持式(5a, 5b)所述的误差符号

h

t

,

q

t

pmb h^t, pmb q^t

hhht,qqqt,可以类比AMP中的式(6)写出OAMP的误差迭代,如下

h

t

=

B

t

q

t

+

W

t

n

(16a)

pmb h^t = pmb B_{t} pmb q^t + pmb W_{t} pmb n tag{16a}

hhht=BBBtqqqt+WWWtnnn(16a)

q

t

+

1

=

η

t

(

x

+

h

t

)

x

(16b)

pmb q^{t+1}=eta_{t}(pmb x+pmb h^t) - pmb x tag{16b}

qqqt+1=ηt(xxx+hhht)xxx(16b)

其中

B

t

=

I

W

t

A

pmb B_{t} = pmb I-pmb W_{t} pmb A

BBBt=IIIWWWtAAA,然后如AMP一样,指定

τ

t

2

=

1

N

E

[

h

t

2

2

]

(17a)

tau^2_{t}=frac {1}{N} mathbb E[{Vert pmb h^t Vert}^2_{2}] tag{17a}

τt2=N1E[hhht22](17a)

v

t

+

1

2

=

1

N

E

[

q

t

+

1

2

2

]

(17b)

v^2_{t+1}=frac {1}{N} mathbb E[{Vert pmb q^{t+1} Vert}^2_{2}] tag{17b}

vt+12=N1E[qqqt+122](17b)

式(17)就是所谓的state evolution,可以对式(17a)做进一步推导

τ

t

2

=

1

N

E

[

h

t

2

2

]

=

1

N

E

[

B

t

q

t

+

W

t

n

2

2

]

=

1

N

{

E

[

B

t

q

t

2

2

]

+

E

[

W

t

n

2

2

]

}

=

1

N

{

E

[

t

r

(

(

q

t

)

T

B

t

T

B

t

q

t

)

)

]

+

E

[

t

r

(

n

T

W

t

T

W

t

n

)

]

}

=

1

N

{

E

[

t

r

(

B

t

T

B

t

)

]

E

[

q

t

2

2

]

+

E

[

t

r

(

W

t

T

W

t

)

]

E

[

n

2

2

]

}

=

1

N

{

E

[

t

r

(

B

t

T

B

t

)

]

N

v

t

2

+

E

[

t

r

(

W

t

T

W

t

)

]

M

σ

2

}

=

E

[

t

r

(

B

t

T

B

t

)

]

v

t

2

+

M

N

E

[

t

r

(

W

t

T

W

t

)

]

σ

2

begin{aligned} tau^2_{t}&= frac {1}{N} mathbb E[{Vert pmb h^t Vert}^2_{2}] \ &=frac {1}{N} mathbb E[{Vert pmb B_{t} pmb q^t + pmb W_{t} pmb n Vert}^2_{2}] \ &=frac {1}{N} { mathbb E [{Vert pmb B_t pmb q^t Vert}^2_{2}] + mathbb E [{Vert pmb W_t pmb n Vert}^2_{2}] } \ &=frac {1}{N} { mathbb E [tr ( (pmb q^t)^T pmb B^T_t pmb B_t pmb q^t) )] + mathbb E [tr(pmb n^T pmb W^T_t pmb W_t pmb n)] } \ &=frac {1}{N} { mathbb E [tr(pmb B^T_t pmb B_t )] cdot E[{Vert pmb q^{t} Vert}^2_{2}] + mathbb E [tr(pmb W^T_t pmb W_t )] cdot E[{Vert pmb n Vert}^2_{2}] } \ &=frac {1}{N} { mathbb E [tr(pmb B^T_t pmb B_t )] cdot N v^2_t + mathbb E [tr(pmb W^T_t pmb W_t )] cdot M sigma^2 } \ &=mathbb E [tr(pmb B^T_t pmb B_t )] cdot v^2_t + frac {M}{N} mathbb E [tr(pmb W^T_t pmb W_t )] cdot sigma^2 end{aligned}

τt2=N1E[hhht22]=N1E[BBBtqqqt+WWWtnnn22]=N1{E[BBBtqqqt22]+E[WWWtnnn22]}=N1{E[tr((qqqt)TBBBtTBBBtqqqt))]+E[tr(nnnTWWWtTWWWtnnn)]}=N1{E[tr(BBBtTBBBt)]E[qqqt22]+E[tr(WWWtTWWWt)]E[nnn22]}=N1{E[tr(BBBtTBBBt)]Nvt2+E[tr(WWWtTWWWt)]Mσ2}=E[tr(BBBtTBBBt)]vt2+NME[tr(WWWtTWWWt)]σ2

注意:上面推导出来的结果与OAMP论文式(23)给出的在系数上有差异,感觉应该上面是正确的,不管如何,思路应该没什么问题。那么,据此,就可以轻易地写出OAMP-state evolution

τ

t

2

=

E

[

t

r

(

B

t

T

B

t

)

]

v

t

2

+

M

N

E

[

t

r

(

W

t

T

W

t

)

]

σ

2

(18a)

tau^2_{t}=mathbb E [tr(pmb B^T_t pmb B_t )] cdot v^2_t + frac {M}{N} mathbb E [tr(pmb W^T_t pmb W_t )] cdot sigma^2 tag{18a}

τt2=E[tr(BBBtTBBBt)]vt2+NME[tr(WWWtTWWWt)]σ2(18a)

v

t

+

1

2

=

E

[

η

t

(

X

+

τ

t

Z

)

X

2

]

(18b)

v^2_{t+1}=mathbb E[{vert eta_t(X+tau_t Z) - X vert}^2] tag{18b}

vt+12=E[ηt(X+τtZ)X2](18b)

其中

X

P

X

(

x

)

X sim P_X(x)

XPX(x),与

Z

N

(

0

,

1

)

Z sim mathcal N(0,1)

ZN(0,1)独立。

4.5 关于OAMP的合理性以及两个重要假设

关于OAMP的合理性,前言部分已经做了简短的铺垫,这里再详细展开。论文作者提出了两个假设,虽然OAMP的证明并不严格,但是基于两个假设展开的讨论还是比较合理的。

假设1:式(16a)中的

h

t

N

(

0

,

τ

t

2

)

pmb h^t sim mathcal N(0,tau^2_t)

hhhtN(0,τt2),并且独立于真实值

x

pmb x

xxx
假设2:式(16b)中的

q

t

+

1

pmb q^{t+1}

qqqt+1里的元素是独立同分布的(i.i.d),并且独立于

A

,

n

pmb A, pmb n

AAA,nnn

一般的条件会有

x

pmb x

xxx是i.i.d.的,独立于

A

,

n

pmb A, pmb n

AAA,nnn,在OAMP中,当迭代次数

t

=

1

t=-1

t=1时,

q

0

=

x

pmb q^0=- pmb x

qqq0=xxx,因此假设2在

t

=

1

t=-1

t=1时是成立的。虽然OAMP的线性迭代式(15a)比AMP的线性迭代式(4a)少了"onsager"这一项,但是只有我们能证明假设2在迭代过程中一直成立,那么式(15a)便是合理的,因为"onsager"的作用就是去除迭代估计结果与

A

pmb A

AAA的相关性。接下来会提出两个推论来更直观地理解OAMP。

4.5.1 从假设2看假设1

推论1:如果假设2是成立的,矩阵

A

pmb A

AAA是酉不变的,

W

t

pmb W_t

WWWt是去相关矩阵,那么就有

h

t

pmb h^t

hhht的元素与

x

pmb x

xxx的元素不相关,以及,

h

t

pmb h^t

hhht的元素彼此之间不相关,而且它们拥有共同的方差,均值为0。

证明:从式(16b)

q

t

+

1

=

η

t

(

x

+

h

t

)

x

pmb q^{t+1}=eta_{t}(pmb x+pmb h^t) - pmb x

qqqt+1=ηt(xxx+hhht)xxx可以看出,

q

t

pmb q^t

qqqt

x

pmb x

xxx具有相关性,这可能会进一步导致

h

t

pmb h^t

hhht

x

pmb x

xxx的相关,因为式(16a)中

h

t

pmb h^t

hhht

q

t

pmb q^t

qqqt生成。但去相关矩阵

W

t

pmb W_t

WWWt的引入可以抑制此相关性,具体描述如下。
因为

A

=

V

Σ

U

T

pmb A=pmb V Sigma pmb U^T

AAA=VVVΣUUUT

W

t

=

U

G

t

V

T

pmb W_t=pmb U pmb G_t pmb V^T

WWWt=UUUGGGtVVVT

B

=

I

W

t

A

=

U

(

I

G

t

Σ

)

U

T

pmb B = pmb I- pmb W_t pmb A = pmb U(pmb I-pmb G_t Sigma) pmb U^T

BBB=IIIWWWtAAA=UUU(IIIGGGtΣ)UUUT,那么

E

U

[

(

B

t

)

i

,

j

]

=

m

=

1

N

E

[

U

i

,

m

U

j

,

m

]

(

1

g

m

λ

m

)

mathbb E_U[(pmb B_t)_{i,j}] = sum_{m=1}^N mathbb E[U_{i,m} U_{j,m}] cdot (1-g_m lambda_m)

EU[(BBBt)i,j]=m=1NE[Ui,mUj,m](1gmλm)

其中

λ

m

,

g

m

lambda_m, g_m

λm,gm分别是矩阵

A

,

W

t

pmb A, pmb W_t

AAA,WWWt的奇异值,并且当

m

>

M

m>M

m>M时,

λ

m

=

g

m

=

0

lambda_m=g_m=0

λm=gm=0
对于一个Haar分布的矩阵

U

pmb U

UUU,有

E

[

U

i

,

m

U

j

,

m

]

=

{

0

;

      

i

j

N

1

;

 

i

=

j

mathbb E[U_{i,m} U_{j,m}]= left{ begin{array}{lr} 0; i neq j \ N^{-1}; i=j end{array} right.

E[Ui,mUj,m]={0;      i=jN1; i=j

那么就有

E

U

[

(

B

t

)

i

,

j

]

=

{

0

;

              

i

j

1

N

t

r

(

B

t

)

;

 

i

=

j

mathbb E_U[(pmb B_t)_{i,j}]= left{ begin{array}{lr} 0; i neq j \ frac {1}{N} tr ( pmb B_t ); i=j end{array} right.

EU[(BBBt)i,j]={0;              i=jN1tr(BBBt); i=j

因为

W

t

pmb W_t

WWWt是去相关矩阵,根据定义3,有

t

r

(

B

t

)

=

t

r

(

I

W

t

A

)

=

0

tr(pmb B_t)=tr(pmb I - pmb W_t pmb A)=0

tr(BBBt)=tr(IIIWWWtAAA)=0,所以

E

[

B

t

]

=

0

mathbb E[pmb B_t] = 0

E[BBBt]=0

假设2给出了条件,

q

t

pmb q^t

qqqt独立于

A

pmb A

AAA,显然也独立于

B

t

pmb B_t

BBBt,那么

E

[

h

t

]

=

E

[

B

t

q

t

]

+

E

[

W

t

n

]

=

E

[

B

t

]

E

[

q

t

]

+

E

[

W

t

]

E

[

n

]

=

0

N

(19)

begin{aligned} mathbb E[pmb h^t] &= mathbb E[pmb B_t pmb q^t] + mathbb E[pmb W_t pmb n] \ &=mathbb E[pmb B_t] mathbb E[pmb q^t] + mathbb E[pmb W_t] mathbb E[pmb n] \ &=pmb 0_N end{aligned} tag{19}

E[hhht]=E[BBBtqqqt]+E[WWWtnnn]=E[BBBt]E[qqqt]+E[WWWt]E[nnn]=000N(19)

h

t

=

B

t

q

t

+

W

t

n

pmb h^t = pmb B_{t} pmb q^t + pmb W_{t} pmb n

hhht=BBBtqqqt+WWWtnnn,要证

x

pmb x

xxx

h

t

pmb h^t

hhht不相关,只需要证

x

pmb x

xxx

B

t

q

t

pmb B_{t} pmb q^t

BBBtqqqt不相关(已经有

E

[

B

t

q

t

]

=

0

N

mathbb E[pmb B_t pmb q^t]=pmb 0_N

E[BBBtqqqt]=000N,所以满足正交性即可)

E

[

h

t

x

T

]

=

E

[

B

t

q

t

x

T

]

=

E

[

B

t

]

E

[

q

t

x

T

]

=

0

N

×

N

(20)

mathbb E[pmb h^t pmb x^T]= mathbb E[pmb B_t pmb q^t pmb x^T]=mathbb E[pmb B_t] mathbb E[ pmb q^t pmb x^T]=pmb 0_{N times N} tag{20}

E[hhhtxxxT]=E[BBBtqqqtxxxT]=E[BBBt]E[qqqtxxxT]=000N×N(20)

所以,根据

E

[

h

t

]

=

0

N

mathbb E[pmb h^t]=pmb 0_N

E[hhht]=000N

E

[

h

t

x

T

]

=

0

N

×

N

mathbb E[pmb h^t pmb x^T]=pmb 0_{N times N}

E[hhhtxxxT]=000N×N(正交性),必然有

x

pmb x

xxx

h

t

pmb h^t

hhht不相关。要证

h

t

pmb h^t

hhht的元素彼此之间不相关,而且它们拥有共同的方差,均值为0,只需证明

h

t

pmb h^t

hhht为对角阵的系数,这里省略。

事实上,因为

E

[

B

t

]

=

E

[

n

]

mathbb E[pmb B_t]=mathbb E[pmb n]

E[BBBt]=E[nnn]都为0,式(19)隐含了

h

t

,

q

t

pmb h^t ,pmb q^t

hhht,qqqt的正交性

E

[

h

t

(

q

t

)

T

]

=

0

N

×

N

(21)

mathbb E[pmb h^t (pmb q^t)^T]=pmb 0_{N times N} tag{21}

E[hhht(qqqt)T]=000N×N(21)

从推论1的证明过程中可以看出,去相关矩阵

W

t

pmb W_t

WWWt在里边起到了重要的去相关作用(间接地借助了

B

=

I

W

t

A

pmb B = pmb I- pmb W_t pmb A

BBB=IIIWWWtAAA)。除此之外,还可以看出OAMP对矩阵

A

pmb A

AAA的特征值没有任何束缚,所以潜在的应用范围会更广一些相对于AMP。

4.5.2 从假设1看假设2

在这一部分,我们尝试基于假设1,来推出假设2,从式(16)可以看出,如果

q

t

+

1

pmb q^{t+1}

qqqt+1

h

t

pmb h^t

hhht独立,那么

q

t

+

1

pmb q^{t+1}

qqqt+1

A

,

n

pmb A, pmb n

AAA,nnn也独立,因为(16b)中

q

t

+

1

pmb q^{t+1}

qqqt+1只跟

h

t

pmb h^t

hhht

x

pmb x

xxx有关系,那么就存在这样一条马尔可夫链(注意上标

t

t

t)

A

,

n

h

t

q

t

+

1

pmb A, pmb n rightarrow pmb h^t rightarrow pmb q^{t+1}

AAA,nnnhhhtqqqt+1

也就是说,如果我们能够证明

q

t

+

1

pmb q^{t+1}

qqqt+1

h

t

pmb h^t

hhht独立,那么假设2就自然而然成立。但遗憾的是我们并不能证明独立性,只能证明正交性,再推广到不相关(推论2将阐述)。

回到定义1里边所阐述的divergence-free函数,利用Lipchitz函数期望的收敛性,我们可以构建一个近似divergence-free函数,如下

η

t

(

r

t

)

=

C

{

η

^

t

(

r

t

)

(

1

N

j

=

1

N

η

^

t

(

r

j

t

)

)

r

t

}

(22)

eta_t(pmb r^t)=C cdot { hat eta_t(pmb r^t)-mathbb (frac{1}{N} sum_{j=1}^N hat eta^{'}_t( r^t_j)) cdot pmb r^t } tag{22}

ηt(rrrt)=C{η^t(rrrt)(N1j=1Nη^t(rjt))rrrt}(22)

备注:式(22)的近似divergence-free函数与OAMP的迭代公式(15)结合,便是真正的实际所采用的OAMP算法。

推论2:如果

η

eta

η是divergence-free函数,那么

E

[

τ

t

Z

η

(

X

+

τ

t

Z

)

]

=

0

(23)

mathbb E[tau_t Z cdot eta(X+tau_t Z )]=0 tag{23}

E[τtZη(X+τtZ)]=0(23)

在证明上式之前,先引入一个重要引理
Stein引理:对

ψ

:

R

R

forall psi: R mapsto R

ψ:RR,使得下式中的期望存在,有

E

[

Z

ψ

(

Z

)

]

=

E

[

ψ

(

Z

)

]

(24)

mathbb E[Z cdot psi(Z)]=mathbb E[psi^{'}(Z)] tag{24}

E[Zψ(Z)]=E[ψ(Z)](24)

ϕ

(

Z

)

=

η

t

(

X

+

τ

t

Z

)

phi(Z)=eta_t(X+tau_t Z)

ϕ(Z)=ηt(X+τtZ),根据Stein引理,

E

[

τ

t

Z

η

t

(

X

+

τ

t

Z

)

]

=

τ

t

E

X

{

E

Z

X

[

Z

η

t

(

X

+

τ

t

Z

)

]

}

=

τ

t

2

E

X

{

E

Z

X

[

η

t

(

X

+

τ

t

Z

)

]

}

=

τ

t

2

E

[

η

t

(

X

+

τ

t

Z

)

]

begin{aligned} mathbb E[tau_t Z cdot eta_t(X+tau_t Z )] &= tau_t cdot mathbb E_X { mathbb E_{Z|X} [Z cdot eta_t(X+tau_t Z )] } \ &=tau^2_t cdot mathbb E_X { mathbb E_{Z|X} [ eta^{'}_t(X+tau_t Z )] } \ &=tau^2_t cdot mathbb E [ eta^{'}_t(X+tau_t Z )] end{aligned}

E[τtZηt(X+τtZ)]=τtEX{EZX[Zηt(X+τtZ)]}=τt2EX{EZX[ηt(X+τtZ)]}=τt2E[ηt(X+τtZ)]

因为

η

t

eta_t

ηt是divergence-free函数,所以

E

[

η

t

(

X

+

τ

t

Z

)

]

=

0

mathbb E [ eta^{'}_t(X+tau_t Z )]=0

E[ηt(X+τtZ)]=0,也就证明了(23)。而(23)等价于

E

[

(

R

t

X

)

(

η

t

(

R

t

)

X

)

]

=

0

(25)

mathbb E[(R^t-X) cdot (eta_t(R^t)-X)]=0 tag{25}

E[(RtX)(ηt(Rt)X)]=0(25)

其中

R

t

=

X

+

τ

t

Z

R^t= X+ tau_t Z

Rt=X+τtZ

式(25)直接推广可以得到

E

[

h

t

(

q

t

+

1

)

T

]

=

0

N

×

N

(26)

mathbb E[pmb h^t (pmb q^{t+1})^T]=pmb 0_{N times N} tag{26}

E[hhht(qqqt+1)T]=000N×N(26)

强调为什么叫OAMP:式(21)说明了线性估计(16a)中,“input-error”

q

t

pmb q^t

qqqt和"output-error"

h

t

pmb h^t

hhht是正交的(由推论1得出);式(26)说明了非线性估计(16b)中,"before-error"

h

t

pmb h^t

hhht和"after-error"

q

t

+

1

pmb q^{t+1}

qqqt+1是正交的(由推论2得出)。这也就是OAMP名字里正交的由来。

4.6 MSE估计和state evolution仿真

在开始这部分之前,先回顾一下式(14c)LMMSE中的参数

v

2

v^2

v2,OAMP迭代公式

r

t

=

s

t

+

W

t

(

y

A

s

t

)

=

s

t

+

W

t

(

A

(

x

s

t

)

+

n

)

pmb r^{t}=pmb s^{t} + pmb W_{t}(pmb y - pmb A pmb s^{t})=pmb s^{t} + pmb W_{t}(pmb A ( pmb x - pmb s^{t} )+ pmb n)

rrrt=ssst+WWWt(yyyAAAssst)=ssst+WWWt(AAA(xxxssst)+nnn),所以LMMSE中的参数

v

2

v^2

v2表示

E

[

(

x

s

t

)

(

x

s

t

)

T

]

=

v

2

I

mathbb E[( pmb x - pmb s^{t} ) ( pmb x - pmb s^{t} )^T]=v^2 pmb I

E[(xxxssst)(xxxssst)T]=v2III

两个MSE:

v

t

2

=

1

N

E

[

q

t

2

2

]

v^2_t=frac{1}{N}{mathbb E[{Vert pmb q^t Vert}^2_2]}

vt2=N1E[qqqt22]

τ

t

2

=

1

N

E

[

h

t

2

2

]

tau^2_t=frac{1}{N}{mathbb E[{Vert pmb h^t Vert}^2_2]}

τt2=N1E[hhht22],它们可以看作是去相关矩阵

W

t

pmb W_t

WWWt和divergence-free函数

η

t

eta_t

ηt的两个参数(这也就是为什么

W

t

pmb W_t

WWWt有下标

t

t

t的原因)

4.6.1 非线性均方误差的估计

非线性均方误差

v

t

2

v^2_t

vt2的估计表达式

v

^

t

2

=

1

N

y

A

s

t

2

2

M

σ

2

t

r

(

A

T

A

)

(27)

hat {v}^2_t=frac{1}{N} frac {{Vert {pmb y - pmb A pmb s^t} Vert}^2_2 - M cdot sigma^2} {tr(pmb A^T pmb A)} tag{27}

v^t2=N1tr(AAATAAA)yyyAAAssst22Mσ2(27)

(27)的理解:

E

[

y

A

s

t

2

2

]

=

E

[

A

x

A

s

t

+

n

2

2

]

=

E

[

A

(

x

s

t

)

2

2

]

+

E

[

n

2

2

]

=

E

[

(

x

s

t

)

2

2

]

E

[

t

r

(

A

T

A

)

]

+

M

σ

2

=

E

[

q

t

2

2

]

E

[

t

r

(

A

T

A

)

]

+

M

σ

2

begin{aligned} mathbb E[{Vert {pmb y - pmb A pmb s^t} Vert}^2_2] &= mathbb E[{Vert {pmb A pmb x - pmb A pmb s^t + pmb n} Vert}^2_2] \ &=mathbb E[{Vert {pmb A( pmb x - pmb s^t )} Vert}^2_2] + mathbb E[{Vert {pmb n} Vert}^2_2] \ &= mathbb E[{Vert {( pmb x - pmb s^t )} Vert}^2_2] cdot mathbb E[tr(pmb A^T pmb A)] + M cdot sigma^2 \ &=mathbb E[{Vert pmb q^t Vert}^2_2] cdot mathbb E[tr(pmb A^T pmb A)] + M cdot sigma^2 \ end{aligned}

E[yyyAAAssst22]=E[AAAxxx