摘要
快照压缩成像(SCI)是一种新型的压缩成像系统,它将多帧图像压缩成单个快照测量,具有低成本、低带宽和高速的传感率。然而,通过应用现有的SCI方法来处理高光谱图像,并不能充分利用底层结构,从而证明重建性能不令人满意。为了解决这一问题,本文旨在提出一种新的有效的方法,利用高光谱图像的两个内在先验,即深度图像去噪和全变差(TV)先验。具体来说,我们提出了一个优化目标来利用这两个先验。通过求解这一优化目标,我们的方法等价于合并一个加权FFDNet和一个2DTV或3DTV去噪器到即插即用框架。大量的数值实验证明,所提出的方法优于几种最先进的替代方案的性能。此外,我们提供了一个详细的收敛分析关于结果的即插即用算法在相对较弱的条件下,如不使用减小的步长。
1. Introduction
压缩传感[5,1]是一种流行的成像技术,可用于捕获视频和高光谱图像。最重要的压缩传感系统之一是所谓的快照压缩成像(SCI),准确地说,SCI使用二维传感器获取高维图像数据,并利用相应的算法重建所需的数据。与传统的压缩感知技术相比,SCI具有低内存、低功耗、低带宽、低成本等优点,因此可以用于高效地捕获高光谱图像。在现有的SCI系统中,编码孔径快照光谱成像(CASSI)[25]是一种具有代表性的高光谱SCI系统,它将不同波长的高光谱图像组合成一个单一的2D图像。
随着硬件的发展,各种针对SCI的重建算法已经被提出。GAP-TV[28]在广义交替投影(GAP)[13]框架下应用了全变分最小化。最近,DeSCI[14]在重建视频和高光谱图像数据方面展示了最先进的结果。如[31]所示,DeSCI可以看作是一种即插即用(PnP)算法,在重构过程中使用秩最小化作为重建的中间步骤。然而,DeSCI的低计算速度阻碍了它的应用。例如,DeSCI需要花费超过6个小时来重建一个大小为1021×703×24的高光谱图像,虽然GAP-TV是一种更快的算法,但它不能重建适合实际应用的高质量的图像。因此,[31]在PnP算法[3]中将FFDNet[34]等深度去噪网络纳入PnP算法。因为FFDNet可以在GPU上执行,所以与DeSCI相比,它运行得非常快。然而,[31]主要关注视频SCI重建,正如我们稍后将展示的,它在高光谱图像上的重建性能并不令人满意。
基本上,应用DeSCI进行高光谱图像重建需要执行GAP-TV来得到其初始值。数值实验表明,该初始值对DeSCI的性能至关重要。如果初始化稍差,DeSCI的性能可能会很差。因此,迫切需要开发新的有效方法来解决上述问题。
DeSCI使用GAP-TV重建结果作为初始化可以看作是一种结合不同先验的原始方式。考虑到不同的内在先验可以相互促进,我们提出了一种将深度图像去噪和TV先验相结合来增强高光谱图像重建的新方法。虽然使用去噪网络如FFDNet或TV正则化本身不能获得满意的结果,但是按照我忙得做法先验结合,可以同时利用去噪网络和TV,使两个先验相互促进,从而获得更多更好的重建结果。
本文所做出的贡献如下:
- 我们提出了一个通用的框架,可以结合深度去噪和全变分先验的SCI高光谱图像重建。在重建过程中,这两种先验可以相互促进,以提高重建图像的质量。
- 我们在模拟和真实数据集上进行了广泛的实验。数值结果验证了该算法的有效性和有效性。
- 我们证明了我们的算法在SCI重建中不减小步长的收敛性。更重要的是,我们是第一个证明加速PnP-GAP在SCI中是收敛的,这是在我们的模拟实验中使用的。
本文的其余部分组织如下。我们将在第2节中介绍SCI模型和相应的重建算法。我们提出的方法和结果算法的收敛性分析在第3节中展开介绍。第4节介绍了广泛的实验。
Related Work
如果我们想从SCI系统中重建图像,就需要使用不同的先验。目前已经提出了多种算法去重建SCI,如基于稀疏性的算法[29,33]、GMM[26,30]、GAPTV[28]和DeSCI[14],其中,DeSCI已经带来了最先进的结果,除了对SCI问题的先验的显式建模外,隐式先验也被用于重建。[31]将FFDNet[34]等各种去噪器集成到PnP算法中,用于视频SCI重建,并获得了良好的重新扫描结果。最近,一些深度学习方法也对SCI问题[11,12,17,18,21,16]取得了良好的效果。与这些方法不同,我们将深度去噪和全变分先验结合到PnP框架中,以促进高光谱图像的SCI重建。
[31]证明了PnP-GAP在使用减小步长的SCI重建中的收敛性。最近,[20]证明了PnP-ADMM的收敛性,而不使用减小步长的收敛性。受这两项工作的启发,我们证明了PnP-GAP和在不适用减小的步长情况下加速PnP-GAP的不动点收敛性。
2. 快照压缩光谱成像技术综述
2.1 SCI model
CASSI (coded aperture snapshot spectral imaging)是一种具有代表性的高光谱图像捕获SCI系统。图1显示了CASSI的传感过程。这个固定的掩模对光谱场景进行了空间编码。被编码的场景被棱镜以光谱的形式分散。最后,用灰度摄像机检测空间光谱编码场景。因此,相机上的快照可以编码场景的几十个光谱带。
Y
=
∑
b
=
1
B
C
b
⊙
X
b
+
Z
(1)
\mathbf{Y}=\sum_{b=1}^{B}\mathbf{C_b}\odot\mathbf{X_b}+\mathbf{Z} \tag{1}
Y=b=1∑BCb⊙Xb+Z(1)
我们重写公式(1)可以得到向量形式如下:
y
=
H
x
+
z
(2)
\mathbf{y}=\mathbf{H}\mathbf{x}+\mathbf{z } \tag{2}
y=Hx+z(2)
其中
y
=
V
e
c
(
X
)
∈
R
n
x
n
y
z
=
V
e
c
(
Z
)
∈
R
n
x
n
y
y=Vec(X)\in \mathbb{R^{n_xn_y}}\ z=Vec(Z)\in \mathbb{R^{n_xn_y}}
y=Vec(X)∈Rnxny z=Vec(Z)∈Rnxny ,高光谱图像向量
x
∈
R
n
x
n
y
B
x∈\mathbb{R^{n_xn_yB}}
x∈RnxnyB表示为
与传统的压缩感知[5]不同,高光谱图像SCI中的感知矩阵
H
∈
R
n
x
n
y
×
n
x
n
y
B
H∈R^{n_xn_y×n_xn_yB}
H∈Rnxny×nxnyB并不是一个密集的矩阵。H的特殊结构可以表示为
H
=
[
D
1
,
.
.
.
,
D
B
]
.
(4)
\mathbf{H}=[\mathbf{D_1},...,\mathbf{D_B}]. \tag{4}
H=[D1,...,DB].(4)
其中
D
b
=
d
i
a
g
(
V
e
c
(
C
b
)
)
∈
R
n
×
n
,
n
=
n
x
n
y
,
b
=
1
,
.
.
.
,
B
\mathbf{D}_b=diag(Vec(\mathbf{C}_b))\in \mathbb{R}^{n\times n},n=n_xn_y,b=1,...,B
Db=diag(Vec(Cb))∈Rn×n,n=nxny,b=1,...,B. 由于
x
∈
R
n
x
n
y
B
\mathbf{x}\in \mathbb{R}^{n_xn_yB}
x∈RnxnyB且
H
∈
R
n
x
n
y
×
n
x
n
y
B
\mathbf{H}\in \mathbb{R}^{n_xn_y\times n_xn_yB}
H∈Rnxny×nxnyB,SCI的采样率为1/B。[9,10]从理论上证明了当B>1时,x从y中恢复是可能的。
2.2. Plug-and-Play Algorithms for SCI reconstruction
SCI的数学模型可以表示为以下反问题:
其中g(x)为正则化,λ为正则化的调优参数。我们将介绍两种图像重建算法,即PnP-GAP和PnP-ADMM[31]。它们在无噪声条件下是等价的。在本节中,我们将简单地给出这些算法的具体步骤。如果想知道如何推导这些步骤,请参考[14]。
PnP-GAP首先在线性流形
y
=
H
x
\mathbf{y}=\mathbf{Hx}
y=Hx上投影数据,然后对投影数据进行去噪。在无噪声条件下,加速PnP-GAP的建议如下:
在无噪声条件下,加速PnP-GAP的性能远远优于PnPGAP。在[14]中讨论了PnP-GAP与加速PnP-GAP之间的关系。[31]进一步介绍了PnPADMM与PnP-GAP之间的关系。
3. Combining Deep Denoising and Total Variation Priors
如前所述,DeSCI对高光谱图像重建采用GAP-TV结果作为其初始值,这对其性能至关重要。然而,这种结合不同先验的方法是相对原始的。我们将在本节中介绍我们使用不同先验的方法。我们将首先介绍我们的方法背后的基本思想。假设在SCI重建算法的每一步中都存在一个最佳的先验。如果具有不同先验的优化算法具有相似的重构过程,则每次迭代的最佳先验在重构过程中应该相似。假设我们想要组合两个先验,我们的方法选择最接近这两个先验的先验,那么这样的先验很可能接近最佳先验。直观地说,不同的先验之间的紧密关系使他们有机会拥有彼此的优势。在此基础上,我们提出了一个基于PnP算法的通用框架,以结合FFDNet和TV。该方法的基本思想如图2所示。
在下面的3.1节中,我们将把PnP算法中的去噪步骤视为MAP估计,推导出去噪器的后验,并提出该方法的优化目标。
3.1. Denoiser into Posterior
去噪步骤
v
k
+
1
=
D
σ
(
x
k
+
1
)
\mathbf{v}^{k+1}=D_\sigma(\mathbf{x}^{k+1})
vk+1=Dσ(xk+1)等价于解决问题
v
k
+
1
=
a
r
g
m
i
n
v
λ
γ
g
(
v
)
+
1
2
∣
∣
x
(
k
+
1
)
−
v
∣
∣
2
2
(11)
\mathbf{v}^{k+1}=argmin_{\mathbf{v}}\frac{\lambda}{\gamma}g(\mathbf{v})+\frac{1}{2}||\mathbf{x^{(k+1)}-\mathbf{v}||_2^2} \tag{11}
vk+1=argminvγλg(v)+21∣∣x(k+1)−v∣∣22(11)
g(v)的近端算子可以作为FFDNet和TV去噪器,从概率的角度来看,公式(11)可以看作是最大的后验(MAP)估计
v
k
+
1
=
a
r
g
m
a
x
p
(
v
∣
x
k
+
1
,
σ
)
=
D
σ
(
x
k
+
1
)
(12)
\begin{aligned} \mathbf{v}^{k+1}&=argmaxp(\mathbf{v}|\mathbf{x}^{k+1},\sigma) \tag{12} \\ &=D_\sigma(\mathbf{x}^{k+1}) \end{aligned}
vk+1=argmaxp(v∣xk+1,σ)=Dσ(xk+1)(12)
现在我们假设去噪器超参数σ存在一个分布,并将其表示为
p
(
σ
∣
x
(
k
+
1
)
)
p(σ|x^{(k+1)})
p(σ∣x(k+1))。利用超参数的分布,我们整合σ来消除超参数
我们很容易看出,这个(14)是我们从去噪器Dσ中得到的后验。
由于我们的方法是被提出来组合两个去噪器的,所以我们将另一个去噪器表示为Dt。对上面的Dt做同样的事情,我们有
然后我们的目标是最小化
p
(
v
∣
x
(
k
+
1
)
)
p(v|x^{(k+1)})
p(v∣x(k+1))和
q
(
v
∣
x
(
k
+
1
)
)
q(v|x^{(k+1)})
q(v∣x(k+1))之间的距离,即,
为了尽量减少这种距离,我们列出了两个挑战及其相应的解决方案:
- p ( v ∣ x ( k + 1 ) ) p(v|x^{(k+1)}) p(v∣x(k+1))和 q ( v ∣ x ( k + 1 ) ) q(v|x^{(k+1)}) q(v∣x(k+1))是未知的分布,这是需要进行估计的;
- p ( v ∣ x ( k + 1 ) ) p(v|x^{(k+1)}) p(v∣x(k+1))和 q ( v ∣ x ( k + 1 ) ) q(v|x^{(k+1)}) q(v∣x(k+1))很难计算,应该将方程(14)和(15)中的积分离散化。
接下来,我们将解决这些挑战,并使问题(16)能够得到解决。
3.2. Minimizing Distance between Posteriors
我们现在通过最小化后验之间的距离来解决上述两个挑战。对应于两个后验
p
(
v
∣
x
)
p(v|x)
p(v∣x)和
q
(
v
∣
x
)
q(v|x)
q(v∣x)的去噪器为FFDNet 和TV。然后,我们对分布
p
(
v
∣
x
,
σ
)
p(v|x,σ)
p(v∣x,σ)和
q
(
v
∣
x
,
t
)
q(v|x,t)
q(v∣x,t)进行建模。
由于FFDNet的训练数据是一对干净图片和噪声图片,噪声分布是高斯分布,因此很自然地将FFDNet的后验
p
(
v
∣
x
,
σ
)
p(v|x,σ)
p(v∣x,σ)建模为高斯分布,因此,我们有
F
F
D
σ
(
x
)
FFD_σ(x)
FFDσ(x)是以
x
x
x为输入的FFDNet,
σ
I
σI
σI 是高斯分布的协方差矩阵。请注意,如果我们汇总所有光谱波段的去噪结果,则
F
F
D
σ
(
x
)
FFD_σ(x)
FFDσ(x)的大小为
n
x
×
n
y
×
B
n_x\times n_y\times B
nx×ny×B。
我们也可以利用高斯分布来模拟TV去噪器的后验分布。因此
T
V
t
(
x
)
∈
R
n
x
×
n
y
×
B
TV_t(x)\in \mathbb{R}^{n_x\times n_y\times B}
TVt(x)∈Rnx×ny×B是均值,
Σ
t
\Sigma_t
Σt是协方差矩阵。为简单起见,我们假设
Σ
t
Σ_t
Σt是一个常数矩阵。稍后我们会看到的,
Σ
t
Σ_t
Σt的实际值在我们的算法中并不重要,所以我们不需要花费更多的精力来建模
q
(
v
∣
x
,
t
)
q(v|x,t)
q(v∣x,t) 的特定方差。
我们现在可以应对第二个挑战。从本质上讲,计算
p
(
v
∣
x
)
p(v|x)
p(v∣x)和
q
(
v
∣
x
)
q(v|x)
q(v∣x)的困难来自于σ和t的连续性。因此,一个非常简单的方法是离散化σ和t。为了离散σ,我们考虑σ∈A,其中A是一个具有固定元素的集合,p(σ|x)是一个离散分布。用同样的方法,我们考虑t∈B,其中B也是一个具有有限元素的集合。因此,我们有两个离散分布
p
(
σ
∣
x
)
p(σ|x)
p(σ∣x)和
q
(
t
∣
x
)
q(t|x)
q(t∣x)。因此,(14)和(15)可以重写为
然后,关键是在问题(16)中指定距离函数
d
i
s
t
(
⋅
,
⋅
)
dist(·,·)
dist(⋅,⋅)。基本上,可以使用各种指标作为距离,如KL散度、L2距离和MMD[7]。为简单起见,我们选择使用
p
(
v
∣
x
)
p(v|x)
p(v∣x)和
q
(
v
∣
x
)
q(v|x)
q(v∣x)的第一矩之间的
L
2
L_2
L2距离,因为它的计算效率和良好的实验结果。换句话说,p(v|x)和q(v|x)的平均值期望接近,因此,我们提出了以下的优化问题:
通过对上述优化问题进行重新排列,并表示
w
σ
f
f
d
=
p
(
σ
∣
x
)
,
σ
∈
A
w_σ^{ffd}=p(σ|x),σ∈A
wσffd=p(σ∣x),σ∈A和
w
t
t
v
=
q
(
t
∣
x
)
,
t
∈
B
w_t^{tv}=q(t|x),t∈B
wttv=q(t∣x),t∈B,我们可以得到以下优化问题:
上述优化问题(19)可以重写为
其中
W
=
[
w
σ
1
f
f
d
,
.
.
.
,
w
σ
∣
A
∣
f
f
d
,
w
t
1
t
v
,
.
.
.
,
w
t
∣
B
∣
t
v
]
T
W=[w_{\sigma _1}^{ffd},...,w_{\sigma _{|A|}}^{ffd},w_{t_1}^{tv},...,w_{t_{|B|}}^{tv}]^T
W=[wσ1ffd,...,wσ∣A∣ffd,wt1tv,...,wt∣B∣tv]T.我们令
X
f
=
[
F
F
D
σ
1
(
x
)
,
.
.
.
,
F
F
D
σ
∣
A
∣
(
x
)
]
,
X
t
=
[
T
V
t
1
(
x
)
,
.
.
.
,
T
V
t
∣
B
∣
]
T
X_f=[FFD_{\sigma_1}(x),...,FFD_{\sigma_{|A|}}(x)],X_t=[TV_{t_1}(x),...,TV_{t_{|B|}}]^T
Xf=[FFDσ1(x),...,FFDσ∣A∣(x)],Xt=[TVt1(x),...,TVt∣B∣]T,然后将半正定矩阵P表示为
优化问题(20)是一个二次规划。由于|A|+|B|在我们的实验中相对较小,因此可以有效地解决(20)问题。
通过求解问题(20),我们得到了两个具有相互优势的相似的后验问题。然后我们简单地将这些后验平均合并,并从合并后验的平均值中得到一个加权去噪器。这种加权去噪器可以很容易地应用于PnP算法中,如下小节所示。
3.3. Algorithm
实际上,我们可以将我们的方法应用于任何即插即用算法,并以PnP-GAP为例。确切地说,我们的算法有以下四个步骤:
- 欧几里得投影。这个过程与公式(9)相同;
- 根据集合A和集合B中的超参数获得去噪图像,分别表示为 { v σ : σ ∈ A } \{v_σ:σ∈A\} {vσ:σ∈A}和 { v t : t ∈ B } \{v_t:t∈B\} {vt:t∈B};
- 解决问题(20),得到 w ^ f f d \hat{w}^{ffd} w^ffd和 w ^ t v \hat{w}^{tv} w^tv的最佳加权系数;
- 获得去噪结果
与经典的PnP-GAP相比,我们的算法需要执行|A|+|B|步骤来去噪并解决一个二次规划问题。该方法的伪代码如算法1所示。
3.4. Fixed-point Convergence
在[20]和[31]的激励下,我们将证明所提出的PnP-GAP收敛于一个不动点。与[31]不同,我们可以在不使用减小步长的情况下证明其收敛性。此外,我们还证明了在模拟实验中所采用的加速PnP-GAP的收敛性。为此,我们对去噪器的假设如下。由于该算法等价于PnP算法中的加权去噪器,因此我们证明了该加权去噪器也符合引理1中的假设1。然后,我们将所提出的PnP-GAP和加速的PnP-GAP分别转换为两个算符。定理1和定理2指出它们的算子在某些条件下是收缩的。我们首先引入在证明中使用的去噪器的假设。
Assumption 1我们假设在我们的方法中使用的所有去噪器
D
σ
:
R
d
→
R
d
D_σ:\mathbb{R}^d→\mathbb{R}^d
Dσ:Rd→Rd都满足
对于
∀
x
,
y
∈
R
d
\forall x,y\in \mathbb{R}^d
∀x,y∈Rd,对于一些
ϵ
>
0
\epsilon >0
ϵ>0.
基本上,我们可以选择小的σ,使
D
σ
D_σ
Dσ接近单位映射。因此,假设假设1是合理的。如果所有的去噪器都符合假设1,则可以证明加权去噪器也符合假设1。
Lemma 1.
S
=
{
D
σ
:
σ
∈
S
}
S=\{D_σ:σ∈S\}
S={Dσ:σ∈S}是一组满足假设1和
∣
S
∣
<
∞
|S|<∞
∣S∣<∞的去噪器,则S的加权去噪器:
D
w
(
x
)
=
∑
σ
∈
S
w
σ
D
σ
(
x
)
D_w(x)=\sum \limits_{\sigma \in S}w_\sigma D_\sigma (x)
Dw(x)=σ∈S∑wσDσ(x)也满足假设1.其中
∑
σ
∈
S
w
σ
=
1
,
w
σ
≥
1
,
∀
σ
∈
S
\sum \limits_{\sigma \in S}w_\sigma=1,w_\sigma \ge1,\forall \sigma \in S
σ∈S∑wσ=1,wσ≥1,∀σ∈S.
为了证明其收敛性,我们还进一步需要以下假设。
Assumption 2 假设
{
R
j
}
j
=
1
n
>
0
\{R_j\}^n_{j=1}>0
{Rj}j=1n>0 表示对于每个空间位置j,该位置的B帧调制掩模至少有一个非零项。我们进一步假设
R
m
a
x
>
R
m
i
n
R_{max}>R_{min}
Rmax>Rmin。
在假设2中,
R
=
H
H
T
=
d
i
a
g
(
R
1
,
⋅
⋅
⋅
,
R
n
)
R=HH^T=diag(R_1,···,R_n)
R=HHT=diag(R1,⋅⋅⋅,Rn),我们定义
3.4.1 Convergence of PnP-GAP
我们表示P为欧几里得投影,Dσ表示去噪器。如下定理说明了PnP-GAP的收敛性。
Theorem 1. 假设
H
H
H满足假设2,则以下运算符
G
=
D
σ
∘
P
G=D_\sigma \circ P
G=Dσ∘P是收缩的,如果Dσ满足假设1且
Remark 1. 在定理1中,G首先投影数据,然后对相当于PnP-GAP的投影数据进行去噪。G是一个收缩,意味着PnP-GAP收敛到一个不动点。
3.4.2 Convergence of accelerated PnP-GAP
为了证明加速PnP-GAP的收敛性,我们需要证明以下算子是一个收缩,即
T
=
1
2
I
+
1
2
(
2
P
−
I
)
(
2
D
σ
−
I
)
(23)
T=\frac{1}{2}I+\frac{1}{2}(2P-I)(2D_\sigma-I) \tag{23}
T=21I+21(2P−I)(2Dσ−I)(23)
其中P是线性流形上的欧几里得投影
y
=
H
x
,
D
σ
y=Hx,D_σ
y=Hx,Dσ是去噪器。
z
(
k
+
1
)
=
T
(
z
(
k
)
)
z^{(k+1)}=T(z^{(k)})
z(k+1)=T(z(k))是PnP-DRS(即插即用的道格拉斯-拉赫福德分裂),其收敛性等价于PnP-ADMM的收敛性。PnP-ADMM的收敛性相当于加速PnP-GAP的收敛性(这些等价物在补充材料中提出)。下面的定理说T是收缩的
**Theorem 2.**假设H满足假设2,设P是线性流形y=Hx上的欧几里得投影,则
T
=
1
2
I
+
1
2
(
2
P
−
I
)
(
2
D
σ
−
I
)
T=\frac{1}{2}I+\frac{1}{2}(2P-I)(2D_\sigma-I)
T=21I+21(2P−I)(2Dσ−I)是收缩的,如果
D
σ
D_\sigma
Dσ满足假设1且
**Remark 2.**从补充材料中,我们可以知道,不同的算法之间存在以下关系:
定理2证明了PnP-DRS的收敛性,即加速的PnP-GAP收敛于一个不动点。