论文解读者：陈康明，赵田田，李朋

编者按：​

对于大多数一阶算法，我们会在收敛性分析时假设函数是凸的，且梯度满足全局 Lipschitz 条件。而本文中，对于某一类特殊函数。我们不仅不要求函数是凸的，也不再要求梯度满足全局 Lipschitz 条件。

考虑复合优化问题
( P ) min ⁡ { Ψ ( x ) = f ( x ) + g ( x ) : x ∈ C ˉ } , \begin{equation}\nonumber (\mathcal{P})\quad \min \{\Psi(x)=f(x)+g(x): x\in\bar{C}\}, \end{equation} (P)min{Ψ(x)=f(x)+g(x):x∈Cˉ},​
其中$\bar{C}$是$C$的闭包，$C$是 ${\mathbb{R}}^d$的非空开子集。对于大多数一阶算法，我们会在收敛性分析时假设$f$和$g$都是凸函数，且$g$的梯度满足全局 Lipschitz 条件。而本文中，我们不仅不要求函数$f$和$g$是凸函数，也不再要求$g$的梯度的满足全局 Lipschitz 条件，而是使用适应函数g几何形状的凸性条件代替。我们重点研究了一种基于 Bregman 距离而非欧式距离的近端梯度法，该方法涵盖了标准的近端梯度法，并且在一定的合理假设下，证明了该方法全局收敛到临界点。为了展示我们的成果的潜力，我们考虑了一类具有稀疏性约束的二次逆问题，这类问题在许多基础应用中经常出现。并且应用我们的方法推导出了该类问题的新的收敛方案，这是这类重要问题的第一个全局收敛的算法。

第一部分：预备知识​

1.1 Bregman 距离

首先我们给出 kernel generating distance 的定义：

定义1.1 (kernel generating distance). 让$C$是${\mathbb{R}}^d$的凸的非空开集，如果函数$h:{\mathbb{R}}^d\to (-\infty ,+\infty ]$满足下面的条件，那么它被称为 kernel generating distance :
(i) $h$是适当的，下半连续的凸函数， 并且$\mathrm{dom}h\subset \bar{C}$,$\mathrm{dom}\partial h=C$。
(ii) 在$\mathrm{dom}h\equiv C$上，$h$是$C^1$的。

我们用$\mathcal{G}(C)$表示这类 kernel generating distance。
给定$h\in \mathcal{G}(C)$，我们可以通过以下方式定义一个近似度量$D_h:\mathrm{dom}h\times \mathrm{int}\mathrm{dom}h\to {\mathbb{R}}_{+}$:

$$D_h(x,y):=h(x)-[h(y)+⟨\nabla h(y),x-y⟩]$$

这个近似度量$D_h$就被称为 Bregman 距离，它衡量了$x$和$y$的接近程度。

由于梯度不等式，对于所有的$x\in \mathrm{dom}h,y\in \mathrm{int}\mathrm{dom}h$，$h$是凸的当且仅当$D_h(x,y)\ge 0$。并且如果$h$是严格凸的，当且仅当$x=y$时，等号成立。值得注意的是，一般情况下$D_h$ 不是对称的，除非$h=|\cdot {|}^2$，这样得到的就是经典欧式距离的平方。

另外，当$h$不是凸函数时，$D_h$的结构形式也是有用的。它衡量了在给定点$x\in \mathrm{dom}h$处$h$的值与其在$y\in \mathrm{int}\mathrm{dom}h$附近的线性近似之间的差异或者说误差。在这种情况下，前面提到的$D_h(x,y)\ge 0$和$D_h(x,y)=0$当且仅当$x=y$都不再成立。然而，$D_h$仍然具有两个简单但显著的性质，这些性质可以从基本的代数运算中得出：

三点恒等式：对于任意$$y,z\in \mathrm{int}\mathrm{dom}$$和$$x\in \mathrm{dom}h$$，我们有$$D_h(x,z)-D_h(x,y)-D_h(y,z)=⟨\nabla h(y)-\nabla h(z),x-y⟩$$

线性可加性：对于任意$\alpha ,\beta \in \mathbb{R}$，以及任意函数$h_1$和$h_2$，我们有$$D_{\alpha h_1+\beta h_2}(x,y)=\alpha D_{h_1}(x,y)+\beta D_{h_2}(x,y)$$
对于所有$x,y\in \mathrm{dom}{h}_1\cap \mathrm{dom}{h}_2$，使得$h_1$和$h_2$在$y$处可导。

1.2 L-smooth adaptable 条件我们想要选择合适的函数$h\in \mathcal{G}(C)$，并用对应的 Bregman 函数$D_h$来代替近似点梯度法中的欧氏距离平方项。注意，本文所考虑的函数$f$和$g$未必是凸函数。其中$g$满足假设：

$g:{\mathbb{R}}^d\to (-\infty ,+\infty ]$是适当的下半连续函数，定义域满足$\text{dom}h\subset\text{dom}g $,且$g$在$C$上连续可微。

基于上述$g$有关假设, 我们可以给出 L-smooth adaptable 的定义如下：

定义1.2 函数对$(g,h)$在$C$上满足 L-smooth adaptable 条件，当且仅当存在$L>0$使得$Lh+g$和$Lh-g$在$C$上都是凸函数。

结合1.1节中 Bregman 函数的定义，容易得到它的一个等价定义：

定义1.2’ 函数对$(g,h)$在$C$上满足 L-smooth adaptable 条件，当且仅当存在$L>0$使得KaTeX parse error: {equation} can be used only in display mode.

上述定义可看作是 L-smooth 条件的推广。如果取$C={\mathbb{R}}^d$, $h=\frac{1}{2}\Vert \cdot {\Vert }^2$, 则对应的不等式可写为
$$\begin{array}{c}\left|D_g(x,y)\right|=|g(x)-g(y)-⟨\nabla g(y),x-y⟩|\le \frac{L}{2}\Vert x-y{\Vert }^2,\forall x,y\in {\mathbb{R}}^d,\end{array}$$

相当于$g$满足 L-smooth条件。

另外，第二节的证明只需要$Lh-g$是凸函数这个条件。我们把它记作L-smad 条件。

第二部分：BPG 算法

2.1 BPG 算法

根据第一节的分析，我们可以作出以下初步假设：

假设2.1 (1)$h\in \mathcal{G}(C)$, 且$\overline{C}=\overline{\text{dom}h}$;

(2) $f:{\mathbb{R}}^d\to (-\infty ,+\infty ]$是适当的下半连续函数，定义域满足$\text{dom}f\cap C\ne \varnothing$;

(3)$g:{\mathbb{R}}^d\to (-\infty ,+\infty ]$是适当的下半连续函数，定义域满足$\text{dom}h\subset \text{dom}g$,
且$g$在$C$上连续可微;

(4) $(h,g)$满足 L-smad 条件；

(5) v ( P ) = inf ⁡ { Ψ ( x ) : x ∈ C ‾ } > − ∞ v(\mathcal{P})=\inf\{\Psi(x):x\in\overline{C}\}>-\infty v(P)=inf{Ψ(x):x∈C}>−∞.

基于以上假设，我们可以利用函数$h$，构造求解问题$\mathcal{P}$的 BPG 算法如下：

不妨记 $T_{\lambda }(x):=\underset{u\in {\mathbb{R}}^d}{\mathrm{argmin}}\left\{f(u)+⟨\nabla g(x),u-x⟩+\frac{1}{\lambda }{D}_h(u,x)\right\}$

为了保证算法中的 (3.4) 式能够顺利求解，我们需要添加如下假设：

假设2.2 对任意的$\lambda >0$，都有$\underset{\Vert u\Vert \to \infty }{\lim }\frac{h(u)+\lambda f(u)}{\Vert u\Vert }=+\infty .$

假设2.3 对任意的$x\in C$，都有 $T_{\lambda }(x)\subset C$.

这两条假设都是易于实现的${}^{[1]}$. 可以证明，在假设2.1—2.3之下，对任意的$x\in \text{intdom}h$和$x\in \text{intdom}h$,$T_{\lambda }(x)$是$C$的非空紧子集。此时，我们认为求解 (3.4) 这一步确实是可行的。

2.2 充分下降性质

在假设2.1—2.3之下，可证明算法具有充分下降性质：

引理2.1 对于任意$x\in \text{intdom}h$，$\lambda >0$以及$x^{+}\in T_{\lambda }(x)$, 都有不等式$$\begin{array}{c}\lambda \Psi (x^{+})\le \lambda \Psi (x)-(1-\lambda L)D_h(x^{+},x).\end{array}$$

由$h$的凸性可知$D_h$是非负函数。结合引理2.1，可得如下定理：

定理2.1 如果假设2.1—2.3成立，$0<\lambda L<1$, { x k } \{x^k\} {xk}是 BPG 算法生成的序列，则有以下结论：

(1) 序列 { Ψ ( x k ) } \{\Psi(x^k)\} {Ψ(xk)}单调不增；

(2) ${\sum }_{k=0}^{+\infty }{D}_h(x^k,x^{k-1})<\infty$, 因此有$D_h(x^k,x^{k-1})\to 0(k\to \infty )$.

(3) ${\min }_{1\le k\le n}{D}_h(x^k,x^{k-1})\le \frac{\lambda }{n}(\frac{\Psi (x^0)-{\Psi }_{\ast }}{1-\lambda L})$，其中${\Psi }_{\ast }=v(\mathcal{P})>-\infty$.

实际上我们不难看出，如果函数$h$满足假设2.1—2.3，那么$h+\frac{\sigma }{2}\Vert \cdot {\Vert }^2$一定也满足假设，其中$\sigma >0$. 因此不妨设$h$是强凸函数，对应的强凸系数为$\sigma$. 此时定理2.1中的 (3) 可推出${\min }_{1\le k\le n}\Vert x^k-x^{k-1}{\Vert }^2\le \frac{\lambda }{n}\frac{\Psi (x^0)-{\Psi }_{\ast }}{\sigma (1-\lambda L)}$.

2.3 收敛速度

为了证明算法的全局收敛性，本节我们设$C={\mathbb{R}}^d$, 并添加了如下假设：

假设2.4 (1) $\text{dom}h={\mathbb{R}}^d$, 且$h$在 ${\mathbb{R}}^d$上是$\sigma -$强凸的；

(2) $\nabla h$和$\nabla g$在${\mathbb{R}}^d$上都是局部 Lipschitz 连续的。

在假设2.1—2.4之下，可证明算法生成的序列 { x k } \{x^k\} {xk}是极小化$\Psi$的一个类梯度下降序列。其定义如下：

定义1.3 记$F:{\mathbb{R}}^d\to (-\infty ,+\infty ]$是适当的下半连续函数。我们称 { x k } \{x^k\} {xk}是极小化$F$的一个类梯度下降序列，当且仅当以下三个条件成立：

(1) 存在${\rho }_1>0$, 使得${\rho }_1\Vert x^k-x^{k-1}{\Vert }^2\le F(x^k)-F(x^{k-1})$对所有$k$成立；

(2) 存在${\rho }_2>0$，使得对任意的$k$都存在${\omega }^{k+1}\in \partial F(x^{k+1})$ ,
满足$\Vert {\rho }_{k+1}\Vert \le {\rho }_2\Vert x^{k+1}-x^k\Vert$；

(3) 对于 { x k } \{x^k\} {xk}的聚点$\bar{x}$,
不妨设$\underset{k\to \infty ,k\in \mathcal{K}}{\lim }{x}^k=\bar{x}$.
此时有${\mathrm{limsup}}_{k\to \infty ,k\in \mathcal{K}}F(x^k)\le F(\bar{x})$.

利用类梯度下降序列的性质，我们可以证明算法的全局收敛性。记$\Psi$的稳定点集合为 crit Ψ = { x ∈ R d : 0 ∈ ∂ Ψ ( x ) = ∂ f ( x ) + ∇ g ( x ) } , \begin{equation}\nonumber \text{crit}\Psi=\{x\in\mathbb{R}^d:0\in\partial\Psi(x)=\partial{f}(x)+\nabla{g}(x)\}, \end{equation} critΨ={x∈Rd:0∈∂Ψ(x)=∂f(x)+∇g(x)},​

序列 { x k } \{x^k\} {xk}所有聚点构成的集合为$\omega (x^0)$. 对于满足定义1.3的序列 { x k } \{x^k\} {xk}和对应的函数$F$, 可证明$\omega (x^0)$是$\text{crit}F$的非空紧子集，且$F$在$\omega (x^0)$中每点的取值是相同的。进一步，我们可得到如下结论：

定理2.2 如果假设2.1—2.4成立，且$0<\lambda L<1$, 则有：

(1) { x k } \{x^k\} {xk}任意聚点都是$\Psi$的稳定点；

(2) 如果$\Psi$满足 KL 性质，那么$\sum \Vert x^{k+1}-x^k\Vert <\infty$且 { x k } \{x^k\} {xk}收敛到某一个稳定点。

第三部分：数值算例

3.1 问题模型 (SQIP)

为证明算法的有效性，作者用提出的算法近似求解一个二次方程问题，问题的目标是近似寻找一个 $x\in {\mathbb{R}}^d$满足下面的一系列方程
$$\begin{array}{c}{x}^T{A}_{i}x\approx b_i,i=1,2,\dots ,m\end{array}$$

其中 $$\begin{gathered} A_i\in \\ {R}^d \end{gathered}$$是对称矩阵， $$\begin{gathered} b_i\in \\ R \end{gathered}$$是包含噪声的测量值。

通常，研究的系统是欠定的，因此一般利用正则项把原始信号的一些先验信息包含进模型。正则项通常用一个函数$f$表示，这个函数可能是非凸、非光滑、扩展值函数 (为包含约束)。当用最小平方模型来描述测量误差，那么问题能够重新描述为
$$\begin{gathered} \begin{array}{c}\text{(QIP)}\min \{\Psi (x):=\frac{1}{4}\sum _{i=1}^m(x^T{A}_{i}x-b_i{)}^2+\theta f(x):x\in \\ {R}^d\}\end{array} \end{gathered}$$
其中 $\theta >0$是一个惩罚参数，主要对数据的真实性和正则项 $f$之间进行平衡。
定义非凸函数 $$\begin{gathered} g: \\ {R}^d\to \\ R \end{gathered}$$
$g(x)=\frac{1}{4}{\sum }_{i=1}^m(x^T{A}_{i}x-b_i{)}^2.$
函数 $g$$在 ${\mathbb{R}}^d$是连续可微的，但是它的梯度不是全局利普希茨连续的，因此不能够采用经典的近端梯度法求解问题(QIP)。

3.2 算法求解

在这一部分，基本空间是 $C\equiv {\mathbb{R}}^d$，非凸函数 $g:{\mathbb{R}}^d\to \mathbb{R}$被定义为
$$g(x)=\frac{1}{4}\sum _{i=1}^m(x^T{A}_{i}x-b_i{)}^2.$$
对于非凸模型，我们考虑下面两种情况：

(a) 凸 $l_1$范数正则项，即 $f:{\mathbb{R}}^d\to \mathbb{R}$，其中$f(x)=\Vert x{\Vert }_1$
(b)非凸 $l_0$球约束。 $f:{\mathbb{R}}^d\to \mathbb{R}$，其中$f(x)={\delta }_{{\mathbb{B}}_0^s}(x)$，$l_0$球上的指示函数，正整数$s<d$，
B 0 s = { x : ∥ x ∥ 0 ≤ s } , \mathbb{B}_{0}^{s}=\{x: \|x\|_{0}\leq s\}, B0s​={x:∥x∥0​≤s},
$\Vert x{\Vert }_0$ 是$l_0$范数，表示向量$x$的非零元素个数。

为了把我们的方法应用到问题(a)和(b)中，我们首先需要选择一个合适的函数 $$\begin{gathered} h\in \mathcal{G}( \\ {R}^d) \end{gathered}$$使得对于$(g,h)$，$\text{L-smad}$成立。这里，我们采用的 $h:{\mathbb{R}}^d\to \mathbb{R}$为
$$h(x)=\frac{1}{4}\Vert x{\Vert }_2^4+\frac{1}{2}\Vert x{\Vert }_2^2$$

现在，我们证明 $\text{L-smad}$成立，即存在 $L>0$使得 $Lh-g$在 ${\mathbb{R}}^d$上为凸。

引理3.1 假设 $g$和 $h$是上面定义的函数，那么对任意 $L$满足 $$L\ge \sum _{i=1}^{m}3\Vert A_i{\Vert }^2+\Vert A_i\Vert |b_i|,$$
函数 $Lh-g$在 ${\mathbb{R}}^d$上为凸函数。

为了把 $2.2$节的结果应用到问题(a)和(b)中，我们观察到上面的函数 $h$在 ${\mathbb{R}}^d$上是 $1-$强凸，很容易看出假设 $2.1-2.4$是成立的。另外，$g$是实多项式函数，因此是半代数函数。函数 $\Vert x{\Vert }_0$和 $\Vert x{\Vert }_1$也是半代数函数([4] 附录5)。因此，由于半代数函数的和是半代数函数，可得模型(a)和(b)的目标函数 $\Psi$是半代数函数，因此提出的BPG算法能够应用到模型(QIP) (a)和(b)，且能够产生一个全局收敛序列收敛到 $\Psi$的临界点。另外，对于模型(a)和(b)，全局收敛策略具有一个简明的显式迭代步，接下来会详细进行介绍。

在BPG算法中，我们需要计算Bregman近似梯度映射：
$$T_{\lambda }(x)=\arg \min \{f(u)+⟨\nabla g(x),u-x⟩+\frac{1}{\lambda }{D}_h(u,x):u\in {\mathbb{R}}^d\}(\lambda >0).$$

对于模型 (a)和(b)，我们将给出这一迭代步能够产生一个显式的解析解。

在描述之前，我们首先介绍一些简便的符号和一些余下章节将用到的著名算子。令 $\lambda >0$并固定任意 $x\in {\mathbb{R}}^d$ 。定义
$$\begin{array}{c}p\equiv p_{\lambda }(x)=\lambda \nabla g(x)-\nabla h(x)\text{(为了简便，通常省略}\lambda \text{和}x)\end{array}\text{\hspace{1em}(3.1)}$$

对于 $(g,h)$，它们梯度的直接计算结果是 $p_{\lambda }(x)$。现在，忽略掉表达式 $T_{\lambda }$中的常数项，可得
$$\begin{array}{c}{T}_{\lambda }(x)=\arg \min \{\lambda f(u)+⟨p_{\lambda }(x),u⟩+h(u):u\in {\mathbb{R}}^d\}.\end{array}\text{\hspace{1em}(3.2)}$$

接下来，我们介绍两个非常著名的算子，它们会用于计算$T_{\lambda }$。
具有参数$\tau$的软阈值算子。对任意$y\in {\mathbb{R}}^d$，

S τ ( y ) = arg ⁡ min ⁡ x ∈ R d { τ ∥ x ∥ 1 + 1 2 ∥ x − y ∥ 2 } = max ⁡ { ∣ y ∣ − τ ， 0 } sgn ( y ) , \begin{equation}\tag{3.3} S_{\tau}(y)=\arg\min_{x\in\R^{d}}\Big\{\tau\|x\|_{1}+\frac{1}{2}\|x-y\|^{2}\Big\}=\max\{|y|-\tau，0\}\text{sgn}(y), \end{equation} Sτ​(y)=argx∈Rdmin​{τ∥x∥1​+21​∥x−y∥2}=max{∣y∣−τ，0}sgn(y),​(3.3)​

其中绝对值按照分量进行计算。具有参数 $\tau$的硬阈值算子。对任意 $y\in {\mathbb{R}}^d$，
$$\begin{array}{c}{H}_{\tau }(y)=\arg \underset{x\in {\mathbb{R}}^d}{\min }\{\Vert x-y{\Vert }^2:x\in {\mathbb{B}}_0^{\tau }\}=\{\begin{array}{l}{y}_i,i\le \tau ,\\ 0,\text{否则,}\end{array}\end{array}\text{\hspace{1em}(3.4)}$$

对于问题(a)和(b)，我们分别建立 $T_{\lambda }$的显式表达式。

命题3.1 ( $l_1$范数正则化的Bregman近似公式) 令 $f=\Vert \cdot {\Vert }_1$ 且对
$x\in {\mathbb{R}}^d$，令
$v(x):=S_{\lambda \theta }(p_{\lambda }(x))$。那么$x^{+}=T_{\lambda }(x)$为
$$x^{+}=-t^{\ast }v(x)=t^{\ast }{S}_{\lambda \theta }(\nabla h(x)-\lambda \nabla g(x)),$$ 是显示公式，其中 $t^{\ast }$是下面方程的唯一正实根，且具有显式公式形式。
$t^3\Vert v(x){\Vert }_2^2+t-1=0$

接下来，我们考虑 $l_0$范数约束的稀疏二次逆问题。首先，我们回顾下下面的结果[5，命题4.3，79页]。

引理3.2 如果 $0\ne a\in {\mathbb{R}}^d$和正整数 $s<d$，可得$$\max \{⟨a,z⟩:\Vert z{\Vert }_2=1,\Vert z{\Vert }_0\le s\}=\Vert {\mathcal{H}}_s(a){\Vert }_2,$$
其中最优解为 $z^{\ast }={\mathcal{H}}_s(a)/\Vert {\mathcal{H}}_s(a){\Vert }_2$。

命题3.2 ( $l_0$范数正则化的Bregman近似公式) 令 $f={\delta }_{{\mathbb{B}}_0^s}$， $x\in {\mathbb{R}}^d$。那么
$x^{+}=T_{\lambda }(x)$为
$$x^{+}=-\sqrt{t^{\ast }}\Vert {\mathcal{H}}_s(p_{\lambda }(x)){\Vert }_2^{-1}{\mathcal{H}}_s(p_{\lambda }(x))$$
其中 $\sqrt{t^{\ast }}\equiv {\eta }^{\ast }$是下面立方方程的唯一非负实根
$$\begin{array}{c}{\eta }^3+\eta -\Vert {\mathcal{H}}_s(p_{\lambda }(x)){\Vert }_2=0.\end{array}\text{\hspace{1em}(3.5)}$$

参考文献：
[1] Bolte, J., Sabach, S., Teboulle, M., & Vaisbourd, Y. (2018). First order methods beyond convexity and Lipschitz gradient continuity with applications to quadratic inverse problems. SIAM Journal on Optimization, 28(3), 2131-2151.

[2] Bauschke, H. H., Bolte, J., & Teboulle, M. (2017). A descent lemma beyond Lipschitz gradient continuity: first-order methods revisited and applications. Mathematics of Operations Research, 42(2), 330-348.

[3] Geiping, J., & Moeller, M. (2018). Composite optimization by nonconvex majorization-minimization. SIAM Journal on Imaging Sciences, 11(4), 2494-2528.

[4] Bolte, Jérôme, Sabach, S. , & Teboulle, M. (2014). Proximal alternating linearized minimization for nonconvex and nonsmooth problems. Mathematical Programming, 146(1-2), 459-494.

[5] Luss, R. , & Teboulle, M. . (2012). Conditional gradient algorithms for rank-one matrix approximations with a sparsity constraint. Siam Review, 55(1), 65-98.