基于非凸复合函数的稀疏信号恢复算法

作者：周洁容李海洋凌军陈浩彭济根来源：《自动化学报》日期：2022-08-27人气：838

根据奈奎斯特(Nyquist)采样定律, 想要实现信号的无失真输出, 采样频率必须在信号带宽的两倍以上. 然而, 在大多数情形下, 按此定律采样所获得的信息是冗余的, 这不仅造成采样的浪费, 而且处理较大带宽的信号时会给硬件系统带来巨大压力. 压缩感知(Compressed sensing, CS)[1]的出现使该问题的解决成为可能. 压缩感知是一种新型的信号采样及重构理论, 利用少量的测量值就可以实现稀疏或可压缩信号的精确重构. 而稀疏信号的重构在众多科学研究和工程应用中十分重要, 如物理学中的量子态层析成像[2]、天体物理学成像[3]、磁共振成像[4]、信号处理[5]、雷达成像[6]等, CS理论的引入加快了上述应用的研究和发展.

数学上, 稀疏信号重构是指从欠定线性系统 $b = A x$ 中恢复原始信号 $x,$ 其中 $b \in R^{n}$ 是测量向量, $A \in R^{n \times m} (n < m)$ 为感知矩阵, $x \in R^{m}$ 是未知的稀疏向量. 最稀疏信号的重构模型是:

其中, ${‖ x ‖}_{0}$ 表示向量 $x$ 中非零元素的个数[7].

已经证明, $(P_{0})$ 模型的直接求解是NP (Non-deterministic polynomial)难的[8], 其计算量会随着稀疏向量维数的增加而增大, 模型的抗噪能力也很差. 为此人们提出了多种方法对 $(P_{0})$ 模型进行求解, 主要方法可归类于启发式方法、凸松弛方法、非凸松弛方法三种类型. 典型的启发式算法有正交匹配追踪[9]、阈值算法[10]、子空间追踪算法[11]、分级正交匹配追踪[12]等, 这些算法重构理论简单、速度较快, 但往往需要更多观测, 算法的重构精度较低、收敛速度较慢, 并且在有噪声情况下信号的恢复精度较低, 因此其应用范围有限. 典型的凸松弛算法有梯度投影算法[13]、BP (Basis pursuit)算法[14]、BPDN (Basis pursuit denoising)算法[15]、IRLS (Iterative reweighed least squares)算法[16]、Bregma迭代法[17]等, 其中最经典的是BP算法, 文献[18]证明了在测量矩阵满足有限等距性质[19]的条件下BP模型与 $(P_{0})$ 模型等价, 但在多种情形下BP模型中的 $L_{1}$ 范数不能充分反映信号的稀疏性特征, 往往难以获得稀疏解. 为此, 人们提出了用非凸泛函替代 $L_{0}$ 范数的方法, 称之为非凸松弛方法.

相较于 $L_{1}$ 范数, 许多非凸泛函能更好地近似 $L_{0}$ 范数, 从而更好地反映信号的稀疏性特征. 典型的非凸松弛算法有NSL0 (Newt-on smooth $L_{0}$ norm)算法[20]、SL0 (Smoothed $L_{0}$ )算法[21]、CTNRAL0 (Composite trigono-metric function null-space reweighted appro-ximate $L_{0}$ norm)算法[22]、SCSA (Successive concave sparsity approximation)算法[23]等. 易见, 松弛方法的核心是寻找 $L_{0}$ 范数的逼近函数, 通过极小化该逼近函数寻得最稀疏解. 显然这种近似模型对 $(P_{0})$ 模型的逼近性能取决于所选取的近似函数对 $L_{0}$ 范数的逼近程度. 因此, 如何构造具有更优逼近性能的近似函数已成为稀疏信号重构问题研究中的重要问题.

最近, 文献[23−24]分别构造了两种近似函数 $g_{σ} (x) = 1 - e^{- | x | / σ} 、$ $h_{σ} (x) = | x | / (| x | + σ)$ 以实现对 $L_{0}$ 范数的逼近, 取得了较好的重构效果. 一个有趣的想法是, 如果将上述每个泛函对信号的作用看成是一次前向(Forward)处理过程, 那么我们是否可以通过泛函的复合实现对信号的深度作用, 从而提高信号重构的性能? 基于这种想法, 本文将上述两个泛函进行复合, 构建一种新的非凸松弛模型, 给出该模型的理论分析, 对该模型提出一种新的近似算法, 并通过数值实验验证算法的有效性以及相较于SL0算法、IRLS算法、BP算法、SCSA算法等经典算法的优越性.

1. 一种用于逼近 $L_{0}$ 范数的非凸复合指数泛函

设 $δ$ 为Kronecker函数, 即

(1)

则对任意的 $x \in R^{m},$ ${‖ x ‖}_{0}$ 可表示为 $\sum_{i = 1}^{m} [1 - δ (x_{i})] .$ 由于函数 $δ (x_{i})$ 是不连续的, 人们总希望找到连续函数 $f_{σ} : R \to R,$ 通过控制参数 $σ (0 < σ < 1)$ 来逼近 $1 - δ (x_{i}),$ 即

(2)

为方便计, 本文称这样的函数为Delta逼近函数(Delta approximating, DA). 若给定这样的逼近函数 $f_{σ}$ , 并对任意 $x \in R^{m}$ 定义 $F_{σ} (x) = \sum_{i = 1}^{m} f_{σ} (x_{i})$ , 则容易看出

这样, 下面的优化问题可看成是 $(P_{0})$ 问题的近似模型

(3)

文献[21, 23−24]分别采用以下三种DA函数来实现对 $L_{0}$ 范数的逼近

(4)

(5)

(6)

并通过对相应近似模型的求解实现了信号的精确重构.

显而易见, 近似式(3)对 $(P_{0})$ 模型的逼近性能取决于所选取的DA函数对 $L_{0}$ 范数的逼近程度. 因此, 如何构造具有更优逼近性能的DA函数已成为稀疏信号重构问题研究中的重要课题. 注意到, 上述DA函数 $h_{σ} (x)$ 和 $g_{σ} (x)$ 都能很好地实现对 $L_{0}$ 范数的逼近. 一个自然的问题是, 是否可以将这两个函数进行复合以实现对 $L_{0}$ 范数的深度逼近? 针对这个问题, 本文对上述两个函数 $h_{σ} (x)$ 和 $g_{σ} (x)$ 的复合进行考察

(7)

其中, “ $\circ$ ” 表示函数复合运算, 显然这是一个DA函数. 不仅如此, 可以从几何图像和理论分析上说明, 在固定参数 $σ$ 后该DA函数相对于 $h_{σ} (x)$ 、 $g_{σ} (x)$ 、 $p_{σ} (x)$ 具有更好的对 $L_{0}$ 范数的逼近性能.

1)几何图像分析

图1展示了 $p_{σ} (x)$ 、 $h_{σ} (x)$ 、 $g_{σ} (x)$ 以及复合函数 $f_{σ} (x)$ 这4种DA函数的几何图像, 其中参数 $σ = 0.1.$ 可以看出, 在 $x = 0$ 附近, $p_{σ} (x)$ 函数曲线最平坦, 而本文所提出的复合函数 $f_{σ} (x)$ 的图像曲线相对其他三种函数最为陡峭, 这表明函数 $f_{σ} (x)$ 具有对 $L_{0}$ 范数最好的逼近效果.

图 1 4种函数在

σ = 0.1

时的一元函数分布

Fig. 1 The unary distribution of the four functions at

σ = 0.1

图2展示了 $p_{σ} (x)$ 、 $h_{σ} (x)$ 、 $g_{σ} (x)$ 和 $f_{σ} (x)$ 4种DA函数在二维情形的等高几何图像, 即 $p_{σ} (x) = h_{σ} (x) = g_{σ} (x) = f_{σ} (x) = 1$ 的图像, 其中参数 $σ = 0.1.$ 可以看出, 在4种函数与直线 $x_{1} = x_{2}$ 的交点中, 函数 $f_{σ} (x)$ 与 $x_{1} = x_{2}$ 所表示的直线得到的交点的分量绝对值最小, 这表明复合指数函数 $f_{σ} (x)$ 能更好地逼近 $L_{0}$ 范数.

图 2

p_{σ} (x)

、

h_{σ} (x)

、

g_{σ} (x)

和函数

f_{σ} (x)

在

σ = 0.1

时的二元函数分布

Fig. 2 The bivariate distribution of

p_{σ} (x),

h_{σ} (x),

g_{σ} (x)

and the function

f_{σ} (x)

σ = 0.1

2)理论分析

上述从几何直观上展示了复合指数函数相对于其他DA函数的更好逼近性能. 不仅如此, 也可以从理论上证明, 对区间 $(0, 0.9)$ 内任意非零参数 $σ,$ 复合指数函数相对于 $p_{σ} (x)$ 、 $h_{σ} (x)$ 、 $g_{σ} (x)$ 都具有更好的逼近性能. 不失一般性, 仅在第一象限内进行讨论, 即假定 $x \geq 0.$ 对区间内任意 $σ (σ \neq 0),$ 因为

(8)

(9)

(10)

(11)

首先, 固定函数值比较对应变量的大小.

假设 $h_{σ} (x) = g_{σ} (x) = r$ , 其中 $0 < r < 1$ , $0 < σ < 0.9$ , 则由式(8)以及式(9)可得对应变量为

因为

令 $R (r) = \frac{r}{1 - r} + \ln (1 - r),$ 有

所以有 $R (r)$ 为增函数, 即对任意 $r$ 有 $R (r) > R (0) = 0,$ 得 $x_{h} > x_{g} .$ 这说明相较于 $h_{σ} (x),$ 函数 $g_{σ} (x)$ 更贴近坐标轴.

下面固定变量 $x,$ 比较函数值 $p_{σ} (x)$ 、 $g_{σ} (x)$ 、 $f_{σ} (x)$ 的大小. 因为

所以, 为讨论 $g_{σ} (x)$ 、 $p_{σ} (x)$ 的大小关系, 只需要比较指数的大小. 易得 $x$ 趋于 $0$ (至少在 $x \leq 2 σ$ 时), 对任意非零 $σ$ 恒有

即 $g_{σ} (x) \geq p_{σ} (x) .$

同理, 对

进行讨论, 由于

所以, 在 $x = 0$ 附近(至少在 $| x | \leq 0.1$ 范围内)有 $f_{σ} (x) \geq g_{σ} (x)$ . 因此有

说明 $f_{σ} (x)$ 的函数曲线更陡峭.

利用函数的对称性, 其他象限也可类似讨论. 所以由上述函数的几何图像和理论分析也能得出 $f_{σ} (x)$ 相对于 $h_{σ} (x)$ 、 $g_{σ} (x)$ 、 $p_{σ} (x)$ 具有更强的稀疏性.

2. 基于复合指数函数的稀疏模型

以上分析表明, 由式(7)定义的复合函数 $f_{σ} (x)$ 具有对 $L_{0}$ 范数更优的逼近性能. 为此, 考虑由该函数所诱导的优化模型

(12)

其中, $F_{σ} (x) = \sum_{i = 1}^{m} f_{σ} (x_{i}),$ 参数 $σ$ 控制着函数 $F_{σ} (x)$ 对 $L_{0}$ 范数的逼近程度.

下面对优化模型 $(C E)$ 与 $(P_{0})$ 之间的关系进行研究, 为此首先证明以下引理.

引理 1. 设模型 $(C E)$ 的最优解为 $x^{*},$ 其稀疏度为 $k,$ 则矩阵 $A$ 中对应解 $x^{*}$ 支撑集的子矩阵 $A_{k}$ 是列满秩的, 其中 $x^{*}$ 的支撑集是指集合 ${i | x_{i}^{*} \neq 0, i = 1, \dots, m} .$

证明. 利用反证法, 假设矩阵 $A_{k}$ 不是列满秩的, 即矩阵的列向量线性相关, 则存在非零向量 $v,$ 其支撑集包含在向量 $x^{*}$ 的支撑集内, 并且使得 $A_{k} v = 0.$ 显然, 向量 $x^{*} \pm v$ 满足约束条件 $A x = b .$ 不失一般性, 假设向量 $v$ 中最大的分量绝对值不超过向量 $x^{*}$ 中最小的分量绝对值.

根据函数 $f_{σ} (x)$ 的严格凹性, 有

于是, 有

即

这表明 $F_{σ} (x^{*} + v) < F_{σ} (x^{*})$ 或 $F_{σ} (x^{*} - v) < F_{σ} (x^{*}) .$ 这与 $x^{*}$ 为模型 $(C E)$ 的最优解相矛盾, 所以假设不成立, 即矩阵 $A_{k}$ 是列满秩的. □

上述引理1表明, 模型 $(C E)$ 的最优解 $x^{*}$ 满足 ${‖ x^{*} ‖}_{0} = k \leq n,$ 即模型的最优解包含在有限集 $Γ = {x | A x = b, A_{s u p p o r t (x)} 是列满秩矩阵}$ 中.

定理 1. 对任意感知矩阵 $A$ 和测量向量 $b,$ 存在一个依赖于 $A$ 和 $b$ 的常数 $σ (A, b),$ 当 $0 < σ < σ (A, b)$ 时, 模型 $(C E)$ 与模型 $(P_{0})$ 等价.

证明. 要证明模型 $(C E)$ 与模型 $(P_{0})$ 等价, 则只需证明模型 $(C E)$ 的最优解为模型 $(P_{0})$ 的最优解即可.

根据引理1, 可以得到以下等价关系

说明对任意 $σ,$ $(C E)$ 模型的最优解 $\tilde{x} \in Γ,$ 其中 $Γ$ 为有限集.

下面证明存在一个常值 $σ (A, b)$ , 当 $0 < σ < σ (A, b)$ 时, 存在唯一一个点 $\tilde{x} \in Γ$ 均是 $(C E)$ 最小化模型的解, 下面利用反证法进行求证.

取一个固定序列 ${σ_{i}^{(1)}} \to 0,$ 由于对每个 $σ,$ 集合 $Γ$ 中均存在一个点为模型对应的解, 再根据集合 $Γ$ 为有限集, 所以必然存在一个点能反复解决对应的 $(C E)$ 模型. 不妨假设该序列为 ${σ_{i}^{(1)}}$ 本身, 对应的解为 ${\tilde{x}}_{1} .$

又假设存在另一个序列 ${σ_{i}^{(2)}} \to 0$ 和一个异于 ${\tilde{x}}_{1}$ 的集合 $Γ$ 中的点 ${\tilde{x}}_{2},$ 使得当 $σ \in {σ_{i}^{(2)}}$ 时, ${\tilde{x}}_{2}$ 均为 $(C E)$ 模型的解.

结合序列 ${σ_{i}^{(1)}}$ 和 ${σ_{i}^{(2)}},$ 构造如下新序列 ${σ_{i}^{(3)}}$

(13)

显然地, 新序列 ${σ_{i}^{(3)}} \to 0,$ 因此必然存在集合 $Γ$ 中的一个点为对应模型 $(C E)$ 的解. 由于序列 ${σ_{i}^{(3)}}$ 由 ${σ_{i}^{(1)}}$ 和 ${σ_{i}^{(2)}}$ 构成, 则这个点一定为 ${\tilde{x}}_{1}$ 和 ${\tilde{x}}_{2}$ 中的一个. 但是无论取哪一个点, 均与 ${\tilde{x}}_{1}$ 为 $(C E)$ 模型对应序列 ${σ_{i}^{(1)}}$ 的解和 ${\tilde{x}}_{2}$ 为 $(C E)$ 模型对应序列 ${σ_{i}^{(2)}}$ 的解中的其中一种假设相矛盾.

因此存在一个常值 $σ (A, b),$ 当 $0 < σ < σ (A, b)$ 时, 存在唯一一个点 $\tilde{x} \in Γ,$ 均是 $(C E)$ 最小化模型的解. 即 $F_{σ} (\tilde{x}) = \underset{A x = b}{m i n} F_{σ} (x) \leq {‖ x ‖}_{0},$ 其中 $x$ 为可行域中的任一点. 令 $σ$ 趋于零, 有 ${lim}_{σ \to 0} F_{σ} (\tilde{x}) = {‖ \tilde{x} ‖}_{0} \leq {‖ x ‖}_{0},$ 即 ${‖ \tilde{x} ‖}_{0} \leq {‖ x ‖}_{0},$ 则有 $\tilde{x}$ 为模型 $(P_{0})$ 的最优解, 即 $\tilde{x}$ 同为模型 $(C E)$ 、模型 $(P_{0})$ 的最优解, 所以定理1成立. □

3. 模型转换

为克服函数 $F_{σ} (x)$ 在原点不可微带来的不便, 可以将稀疏优化模型 $(C E)$ 转换为以下优化问题. 令未知量 $x = u - v,$ 其中 $u, v \in R^{m},$ $u$ 拥有 $x$ 中所有正元, 其余元素为零; $v$ 拥有 $x$ 中所有负元的绝对值, 其余元素也为零. 用 $z = [u^{T}, v^{T}]^{T} \in R^{2 m}$ 表示拼接向量. 经过替换, 易得

此时, 约束条件 $A x = b$ 转换为 $[A, - A] z = b .$ 那么无噪模型 $(C E)$ 转换为

(14)

为了便于求解式(14), 利用文献[25]中介绍的 $M M$ 优化方法. $M M$ 优化方法指当目标函数较难实现优化时, 通常可以选择更容易优化的替代目标函数, 当替代函数满足一定的条件时, 其最优解能够无限逼近原目标函数的最优解. 通过利用凹函数 $F_{σ} (z)$ 的一阶判别条件以及 $M M$ 优化方法对目标函数 $F_{σ} (z)$ 进行放缩, 可得

记

其中, $z, \tilde{z}$ 是可行域中的点, $λ (λ > 1 / σ^{2})$ 为正常值, $H_{σ} (z, \tilde{z})$ 为 $F_{σ} (z)$ 的一个上界函数. 忽略常值后, 式(14)的解可由下式迭代求得

(15)

其中, $z_{k}^{σ} (k = 1, 2, 3, \dots)$ 为可行域中的点, 取值与 $σ$ 相关.

接下来讨论由式(15)产生的序列 ${z_{k}^{σ}}$ 的收敛性.

命题1. 式(15)产生的序列 ${z_{k}^{σ}}$ 收敛到式(14)的局部极小值.

证明. 1)首先证明式(15)迭代产生的函数序列 ${F_{σ} (z_{k}^{σ})}$ 的有界性和收敛性.

根据凹函数的一阶判别条件有

(16)

再由式(15)的定义有

(17)

在式(16)中令 $z = z_{k + 1},$ 结合式(17)得

即

(18)

其次, 因为

有

(19)

所以由式(18)、(19)可知, 序列 ${F_{σ} (z_{k}^{σ})}$ 是单调有界收敛的.

2) 再证明式(15)产生的序列 ${z_{k}^{σ}}$ 的收敛性.

根据式(16)、(17)得

所以有

当 $N \to \infty,$ 利用式(19)得

是有界的, 所以序列 ${z_{k}^{σ}}$ 收敛.

3) 最后证明式(15)产生的序列 ${z_{k}^{σ}}$ 收敛到式(14)的局部极小值. 假设式(15)迭代产生的序列 ${z_{k}^{σ}}$ 收敛到 ${\tilde{z}}^{σ},$ 即

(20)

下面证明 ${\tilde{z}}^{σ}$ 是

(21)

的最优解.

记

接下来讨论是否有 $z_{0}^{σ} = {\tilde{z}}^{σ}$ 成立. 由于

又注意到 $ℵ = {z_{k}^{σ}}_{k = 1}^{\infty} \cup {{\tilde{z}}^{σ}} \cup {z_{0}^{σ}}$ 有界, 所以存在 $M > 0,$ 使得对任意 $z \in ℵ,$ 有 ${‖ z ‖}_{2} < M .$ 所以有

(22)

再利用假设式(20)有

(23)

即对任意 $ε_{1} > 0,$ 存在 $K_{1} > 0,$ 对任意 $k > K_{1},$ 有

(24)

又根据初等函数 $\nabla F_{σ} (z)$ 在定义区间内关于变量 $z$ 连续, 有

(25)

由此可知, 对任意 $ε_{2} > 0,$ 存在 $K_{2} > 0,$ 对任意 $k > K_{2},$ 有

(26)

故结合式(23)、式(25)可知, 对任意 $z \in ℵ,$ 对于式(22)可得

(27)

即对任意 $z \in ℵ,$ 有任意 $ε_{3} > 0,$ 存在 $K_{3} > 0,$ 对任意 $k > K_{3},$ 有

(28)

根据 $z_{0}^{σ}$ 和 $G_{0} (z)$ 的定义以及不等式(28)有

所以有

则由此可得

(29)

令 $ε = m i n (ε_{1}, ε_{2}, ε_{3}), K = m a x (K_{1}, K_{2}, K_{3}),$ 当 $k > K$ 时, 对不等式(29)有

(30)

由式(30)可知

对左端加减同一项, 有

(31)

对式(31)移项, 有

(32)

再对式(32)不等号左端进行放缩, 得

(33)

最后结合式(32)、(33)得

即

从而由 $ε$ 的任意性, 有 ${‖ z_{0}^{σ} - {\tilde{z}}^{σ} ‖}_{1} = 0,$ 即有 $z_{0}^{σ} = {\tilde{z}}^{σ},$ 由此可知 ${\tilde{z}}^{σ}$ 是式(21)的最优解. 令

优化模型(21)简记为

(34)

约束优化问题式(34)的解为 ${\tilde{z}}^{σ},$ 说明在极小值点 ${\tilde{z}}^{σ}$ 处必不存在可行域内的下降方向, 即不存在非零向量 $p \in R^{2 m},$ 使得

其中, $i, j \in I ({\tilde{z}}^{σ})$ , $I ({\tilde{z}}^{σ})$ 为点 ${\tilde{z}}^{σ}$ 的起作用下标集.

再由

易得

所以对于稀疏优化模型(14)在点 ${\tilde{z}}^{σ}$ 处也不存在可行域内的下降方向, 即 ${\tilde{z}}^{σ}$ 为模型(14)的局部极小值. □

4. 算法设计

为给出模型(14)的有效求解算法, 需要对算法初始值的选择进行分析, 首先给出复合函数的如下性质.

引理 2. ${lim}_{σ \to \infty} σ^{2} F_{σ} (x) = {‖ x ‖}_{1} .$

证明. 已知 $f_{σ} (x) = 1 - e^{- \frac{| x |}{σ (| x | + σ)}},$ 将 $f_{σ} (x)$ 有关 $| x |$ 的麦克劳林公式记为

其中, $g_{σ} (x)$ 为麦克劳林公式的余项, $f_{σ}^{(n)} (0)$ 表示 $f_{σ} (x)$ 在 $0$ 点处的 $n$ 阶导数.

所以有

(35)

由于

其中, 关于 $e^{- | x | / σ (| x | + σ)} / {(| x | + σ)}^{3}$ 、 $e^{- | x | / σ (| x | + σ)} / {(| x | + σ)}^{4}$ 的任意阶导数满足有关 $(| x | + σ)^{- 1}$ 的次方数不低于 $3$ 次, 则 $f_{σ}^{(n)} (0) (n \geq 2)$ 关于参数 $σ^{- 1}$ 的次方数不低于 $3$ 次. 由此可得

(36)

利用式(35)、(36), 有

□

引理2表明, 在 $σ \to \infty$ 时, 最小化 $F_{σ} (x)$ 近似于最小化 ${‖ x ‖}_{1} .$ ${‖ x ‖}_{1}$ 的凸性使得算法陷入局部极值的可能性减少, 从而提高算法重构精度. 鉴于此, 本文在设计算法时取初始值为以下问题

的解.

在此基础上, 提出如下算法.

步骤 1. 输入: $A, b, ε_{1}, ε_{2}, ε_{3} .$

步骤 2. 初始化:

1)选取一组下降序列 $σ, σ_{k + 1} = c \cdot σ_{k}, c \in (0, 0.5) .$

2)内、外层循环的终止阈值 $ε_{1}, ε_{2}, ε_{3} > 0.$

3)设初始值

步骤 3. 算法迭代:

1) $i = 0,$ $σ = σ_{0} .$

2) $w h i l e d_{1} > ε_{1} d o$

　a) ${{\hat{z}}_{0} = z_{k}, l = 0, k = k + 1$

　b) $w h i l e d_{2} > ε_{2} d o$

　c) ${l = l + 1$

　d) ${\hat{z}}_{l} = a r g \underset{z}{m i n} {G_{l - 1} (z) | [A, - A] z = b, z \geq 0}$

　e) $d_{2} = \frac{{‖ {\hat{z}}_{l} - {\hat{z}}_{l - 1} ‖}_{2}}{{‖ {\hat{z}}_{l - 1} ‖}_{2}}}$

　f) $z_{k} = {\hat{z}}_{l}$

　g) $d_{1} = \frac{{‖ z_{k} - z_{k - 1} ‖}_{2}}{{‖ z_{k - 1} ‖}_{2}}$

　h) $σ = c \cdot σ}$

　i) $x_{k} = z_{k} (1 : m) - z_{k} (m + 1 : e n d)$

　j) $U = d i a g {\frac{1}{({x_{k - 1})}_{1} + ε_{3}}, \dots, \frac{1}{{(x_{k - 1})}_{m} + ε_{3}}}$

　k) $x_{k} = a r g \underset{x}{m i n} {{‖ U x ‖}_{1} | A x = b}$

步骤 4. 输出稀疏向量: $x = x_{k} .$

其中, 步骤 3 算法迭代环节由内、外两部分循环构成. 外循环为步骤 a) ~ b) 以及 f) ~ h), 循环利用参数 $σ$ 实现复合函数对 $L_{0}$ 范数的逐次逼近. 内循环即步骤c) ~ e), 结合外点罚函数法[26]和共轭梯度法迭代求解模型(14), 具体为利用外点罚函数法引入正则参数实现模型(14)的无约束转换, 再通过共轭梯度法求解无约束优化问题. 最后在步骤j) ~ k), 针对共轭梯度法求得的非稀疏解, 利用加权L1范数最小化[27]进行稀疏化处理. 其中, $d_{1}$ 和 $d_{2}$ 分别表示外循环和内循环连续迭代的解之间的相对误差, 并用于判断循环是否停止. 值得关注的是, 整个算法的核心环节即步骤d), 是有关本文模型(14)的凸优化求解过程, 具体步骤如下.

模型利用外点罚函数法引入正数 $M,$ 本文取 $M = 1,$ 其放大系数取为 $5,$ 将式(15)转化为以下无约束问题

(37)

其中

(38)

式(38)中, $I$ 是大小为 $2 m \times 1$ 的单位向量.

实验中, 函数 $L_{z_{k}^{σ}} (z, M)$ 有关 $z$ 的一阶、二阶次梯度[28]分别取为

(39)

(40)

其中, $V (z) = (m i n (0, z_{1}), \dots, m i n (0, z_{2 m}))^{T} .$

再利用共轭梯度法迭代求解式(37)的 $z_{k + 1}^{σ},$ 其迭代更新格式为

(4.7)

其中, 步长因子

(4.8)

下降方向

(4.9)

其中, 令 $z_{k}^{σ} = z_{0}$

(4.10)

综上, 本文提出的基于非凸复合函数的稀疏信号恢复算法(Non-convex composite sparse, NCCS)参考了文献[23]中的SCSA算法, 主要包括三个部分. 一是算法的初始值选择, 取值为 $L_{1}$ 范数最小化问题的解, 参考文献[29], 其计算量为 $O (m^{3}) .$ 二是无约束问题的求解, 通过外点罚函数法和共轭梯度法迭代求解最优值. 由式(39) ~ 式(44) 可以看出, 第二部分每次迭代的主要操作是内循环的两部分矩阵乘法, 即

计算复杂度分别为 $O (4 m + 2 m n)$ 和 $O (4 m^{2} + 4 m^{2} n),$ 则第二部分的计算量为 $O (l (4 m + 2 m n + 4 m^{2} + 4 m^{2} n)),$ 其中 $l$ 表示内循环的迭代次数. 三是迭代结果的稀疏化, 利用加权 $L_{1}$ 范数最小化[27]对凸优化结果进行稀疏化处理, 其计算量为 $O (m^{3}) .$

5. 实验仿真与结果分析

为验证本文提出的算法在重构性能上的优越性, 本节设计了几组有关SL0算法[21]、BP算法[15]、IRLS算法[16]、SCSA算法[23]和本文介绍的NCCS这5种算法在重构性能上的对照实验. 在仿真实验中, 感知矩阵 $A$ 和稀疏原信号 $x$ 的具体取值为: 矩阵 $A$ 的大小取为 $250 \times 500,$ 其中矩阵元素服从零均值、单位方差的高斯分布, 矩阵的列具有单位 $L_{2}$ 范数; 稀疏原信号 $x$ 的维数取为 $500,$ 非零分项服从正态分布. 实验将讨论信号稀疏度在区间 $[20, 110]$ 范围内各算法的重构情况. 对于等式约束 $A x = b,$ 在已知感知矩阵 $A$ 和测量向量 $b$ 的情况下, 分别利用上述5种算法对稀疏原信号进行 $100$ 次仿真实验, 讨论算法的平均重构性能. 为展示本文所提出的模型相对于最新文献[23]所提算法的先进性, 本文对5种算法采用了与文献[23]相同的参数选择, 即

1)对SL0算法: 令 $m u = 2, σ_{min} = 10^{- 4}, c = 0.8,$ $L = 8.$

2)对IRLS算法: 令 $p = 0.5.$

3)对BP算法: 令 $l_{1} = 1.$

4)对SCSA算法: 令 $ε_{1} = 10^{- 3}, ε_{2} = 10^{- 2}, c = 0.1.$

5)本文NCCS算法: 令 $ε_{1} = 10^{- 3}, ε_{2} = 10^{- 2}, ε_{3} = 1 0^{- 1}, c = 0 .1, σ_{0} = min {2 {‖ x_{0} ‖}_{0}, α},$ 其中 $α$ 是待定数. 仿真实验一将给出模拟实验, 以获得最佳实验数值.

实验中的重构性能包括:

1)重构信噪比(Signal noise ratio):

2)重构误差(Mean square error):

3)归一化均方差(Normalized mean square error):

4)支撑集恢复成功率(Recovery success rate of the support set):

其中, $\tilde{x}$ 表示稀疏原信号, $\hat{x}$ 表示算法的恢复信号; $T_{1}$ 为稀疏原信号的支撑集, $T_{2}$ 为恢复信号在对应支撑集 $T_{1}$ 上的支撑恢复集. 具体取值分别为

仿真实验一. 待定数 $α$ 的最佳选值.

首先对NCCS算法中待定数 $α$ 进行模拟实验, 观察不同大小的数值 $α$ 对算法运行时间的影响, 以获得更优的实验参数. 为节省计算时长, 考虑 $α$ 大小为 $0 .1$ 到 $0 .9,$ 间隔为 $0 .1,$ 分别进行 $100$ 次仿真实验. 从图3可以看出, 各 $α$ 值对算法运行时间的影响差距不大, 总体趋势是算法运行时间随 $α$ 的增加而减少; 在 $α = 0 .8$ 时, 算法在稀疏度 $k \leq 100$ 时皆能保证运行时间最短. 信号的稀疏度 $k$ 是指信号向量中非零元素的个数, 为便于编程, 在后续实验中待定数 $α$ 全部选定为 $0 .8 .$

图 3 待定数

α

对NCCS算法运行时间的影响

Fig. 3 The influence of undetermined number

α

on the running time of NCCS algorithm

仿真实验二. NCCS算法对稀疏原信号重构的实验.

图4是稀疏原信号以及由NCCS算法获得的恢复信号的仿真结果. 图4的稀疏度取值为 $k = 65.$ 由图4可见, 恢复信号和稀疏原信号基本吻合, 实验结果显示算法的重构误差为 $2 .6933 \times 10^{- 17} .$ 由此可以看出, 本文提出的算法恢复出的稀疏信号很接近稀疏原信号.

图 4 NCCS算法的一维信号重构仿真图, 信号大小为

500 \times 1,

稀疏度为

65

Fig. 4 One-dimensional signal reconstruction simulation diagram of NCCS algorithm, the signal size is 500 × 1, the sparsity is 65

仿真实验三. 各算法的重构性能和稀疏度的变化关系.

图5为5组实验算法的重构误差 $(M S E)$ 和稀疏度 $k$ 的变化关系. 其中IRLS算法的重构误差最高; 由于最速下降法在迭代过程中存在“锯齿效应”[30], 实验中SL0算法的重构误差也很大, 且SL0 与BP、SCSA这3种算法都存在重构误差随着稀疏度的增加而增大的趋势; 而NCCS算法的重构误差相对最小, 整体变化较稳定.

图 5 SL0、IRLS、BP、SCSA、NCCS 5种算法的重构误差和稀疏度的变化关系

Fig. 5 The relationship between the reconstruction error and sparsity of the five algorithms of SL0, IRLS, BP, SCSA, and NCCS

图6表示随着稀疏度在区间 $[20, 110]$ 的等量增加, 5种算法的重构信噪比 $(S N R)$ 的变化情况. 其中SL0、BP和IRLS 3种算法的稀疏度越高, 重构信噪比越低; NCCS算法与SCSA算法的重构信噪比在稀疏度区间 $[20, 110]$ 上表现较稳定. 在实验的各组稀疏度下, NCCS算法的重构信噪比远高于SL0、BP和IRLS这3种算法, 略高于SCSA算法. 所以使用NCCS算法对信号的恢复精确度有一定的提高.

图 6 SL0、IRLS、BP、SCSA、NCCS 5种算法的重构信噪比和稀疏度的变化关系

Fig. 6 The relationship between the reconstructed signal-to-noise ratio and sparsity of the five algorithms of SL0, IRLS, BP, SCSA, and NCCS

图7展示5种算法的重构运行时间和稀疏度 $k$ 的变化关系. 实验结果显示各算法的运行时间随稀疏度的增大而增加, 稀疏度区间在 $[20, 90]$ 时, 各算法的实验平均运行时间能保持在 $2$ s以内. 在稀疏度 $k > 90$ 后, SCSA算法重构运行时间迅速增加; 当稀疏度 $k > 120$ 后, 该算法的平均重构运行时间不低于 $20$ s. 本文提出的NCCS算法由于参数 $σ$ 值的多次衰减导致算法进行多次迭代, 因此也造成了运行时间相较于剩下的3种算法要高, 但对比最新文献[23]的SCSA算法, NCCS算法有更好的实验效果.

图 7 SL0、IRLS、BP、SCSA、NCCS 5种算法的运行时间和稀疏度的变化关系

Fig. 7 The relationship between the running time and sparsity of the five algorithms of SL0, IRLS, BP, SCSA, and NCCS

图8是各算法的支撑集恢复成功率 $(R S S)$ 和稀疏度 $k$ 的变化关系. 实验结果显示SL0、BP和IRLS这3种算法的稀疏向量支撑集恢复成功率随稀疏度的增加而下降, SCSA算法在稀疏度 $k = 100$ 、 $k = 110$ 时, 其支撑集恢复成功率分别为 $R S S = 0 .9980$ 、 $R S S = 0 .9886,$ 考虑为实验次数的不足而出现的偏差. 相比之下, 在NCCS算法的重构仿真实验中, $R S S$ 值保持稳定, 信号的支撑集在整个实验区间内实现完全恢复.

图 8 SL0、IRLS、BP、SCSA、NCCS 5种算法的支撑集恢复成功率和稀疏度的变化关系

Fig. 8 The relationship between the recovery success rate of the support set and sparsity of the five algorithms of SL0, IRLS, BP, SCSA, and NCCS

图9以及表1分别为实验中5组算法恢复向量的归一化均方差 $(N M S E)$ 和稀疏度 $k$ 的变化关系图像曲线和具体数值记录. 由图9和表1可见, NCCS算法和SCSA算法的实验结果远小于SL0算法、BP算法和IRLS算法的实验值.

图 9 SL0、IRLS、BP、SCSA、NCCS 5种算法的归一化均方差和稀疏度的变化关系

Fig. 9 The relationship between the normalized mean square error and sparsity of the five algorithms of SL0, IRLS, BP, SCSA, and NCCS

表 1 5种算法的归一化均方差的数值记录

Table 1 Numerical records of the normalized mean square error of five algorithms

算法	归一化均方差 $(N M S E)$
算法	k = 20	k = 30	k = 40	k = 50	k = 60	k = 70	k = 80	k = 90	k = 100	k = 110
SL0	4.11×10−10	2.86×10−10	2.06×10−10	1.55×10−10	1.16×10−10	9.38×10−11	1.46×10−4	1.28×10−3	9.22×10−3	3.13×10−2
IRLS	2.02×10−4	4.08×10−4	1.15×10−3	4.23×10−3	1.07×10−2	2.79×10−2	5.04×10−2	7.20×10−2	1.14×10−1	1.38×10−1
BP	1.65×10−22	4.49×10−22	5.68×10−22	4.95×10−22	4.19×10−22	4.38×10−15	1.18×10−5	9.52×10−4	1.03×10−2	3.57×10−2
SCSA	2.05×10−23	1.94×10−22	1.94×10−22	2.45×10−22	2.50×10−22	2.96×10−22	3.18×10−22	1.27×10−16	1.27×10−13	3.57×10−12
NCCS	5.36×10−27	2.10×10−26	4.81×10−24	4.98×10−24	5.87×10−15	1.85×10−17	3.35×10−16	2.62×10−17	3.58×10−16	1.00×10−16

综上所述, 本文提出的NCCS算法的综合重构性能要比实验的对照算法好, 所以利用NCCS算法能更好地恢复稀疏信号向量.

6. 结论

本文运用MM优化、外点罚函数法、共轭梯度法等方法, 借鉴文献[23]提出的SCSA算法思想, 提出了一种新的稀疏信号重构算法, 称为NCCS算法. 该算法利用逼近性能更优的非凸复合指数函数实现对 $L_{0}$ 范数的逼近. 仿真实验表明, 本文所提出的NCCS算法不仅有效可行, 而且在无噪声情况下, 对比SL0、IRLS、BP、SCSA这4种算法, NCCS算法在稀疏信号恢复实验中表现出更优越的重构性能.

关键字：优秀论文

上一篇：基于改进动态系统稳定估计器的机器人技能学习方法
下一篇：三相结构离子交换膜的构筑及应用研究

栏目分类

热门排行

推荐信息

期刊知识