凸差算法(DC Algorithm)
一、基本概念
1.1凸差问题
为了解决非凸规划问题，1975年引入，一般形式如下：
i
n
f
{
f
(
x
)
=
g
(
x
)
−
h
(
x
)
,
x
∈
R
n
}
inf\{f(x)=g(x)-h(x),x\in \mathbb R^n\}
inf{f(x)=g(x)−h(x),x∈Rn}
其对偶形式为：
i
n
f
{
f
∗
(
y
)
=
h
∗
(
y
)
−
g
∗
(
y
)
:
y
∈
R
n
}
inf\{f^*(y)=h^*(y)-g^*(y):y \in \mathbb R^n\}
inf{f∗(y)=h∗(y)−g∗(y):y∈Rn}
1.2共轭函数(conjugate function)
这个概念和复数的共轭没有一点关系，事实上从牛津词典上conjugate一词只是表示(动词)的变化形式，数学维基百科里也有四种共轭函数的形式，第四种才是凸分析中的共轭函数。函数
f
:
X
→
R
f:X\rightarrow R
f:X→R的共轭函数为
g
∗
(
y
)
=
s
u
p
x
∈
X
{
<
x
,
y
>
−
f
(
x
)
}
g^*(y)=\mathop{sup}\limits_{x\in X}\{<x,y>-f(x)\}
g∗(y)=x∈Xsup{<x,y>−f(x)}
共轭函数的自变量变成了
y
y
y，
x
x
x变成了它的一个参数。它有一个很好的性质：一定是凸函数，甚至是多面凸函数（下面将讲到什么是多面凸函数）。因为它是一系列关于
y
y
y的仿射函数(线性函数)中取的最大值，所以斜率肯定是分段增大的，可以从下面这张图很清晰地看出这一点。另一方面也可以直接对
y
y
y求导。
这个共轭函数通过对 x x x求导并令其等于0可以求出确切解，想进一步了解的可以阅读参考文献[1] 中的3.3节，或者它的译著参考文献[3] 。
1.3次微分(subdifferntial)和次梯度(subgradient)
函数
φ
\varphi
φ在点
x
0
x_0
x0处的次梯度定义如下：
∂
φ
(
x
0
)
:
=
{
y
0
∈
R
n
:
φ
(
x
)
≥
φ
(
x
0
)
+
<
x
−
x
0
,
y
0
>
,
∀
x
∈
R
n
}
\partial \varphi(x_0):=\{y_0 \in \mathbb{R}^n:\varphi(x)\geq \varphi(x_0)+<x-x_0,y_0>,\forall x \in \mathbb{R}^n\}
∂φ(x0):={y0∈Rn:φ(x)≥φ(x0)+<x−x0,y0>,∀x∈Rn}
一般来说，次梯度是闭凸集，是导数概念在凸函数上的扩展。次微分中的一个元素
y
0
∈
∂
φ
(
x
0
)
y_0 \in \partial \varphi(x_0)
y0∈∂φ(x0)称作函数
φ
\varphi
φ在点
x
0
x_0
x0处的次梯度。当且仅当函数
φ
\varphi
φ在
x
0
x_0
x0处可微，那么闭集
∂
φ
(
x
0
)
\partial \varphi(x_0)
∂φ(x0)变为导数
{
∇
φ
(
x
0
)
}
\{\nabla \varphi(x_0)\}
{∇φ(x0)}。
比如，绝对值函数 φ ( x ) = ∣ x ∣ \varphi(x)=|x| φ(x)=∣x∣： x < 0 x<0 x<0时， ∂ φ ( x ) = − 1 \partial \varphi(x)=-1 ∂φ(x)=−1; x > 0 x>0 x>0时， ∂ φ ( x ) = 1 \partial \varphi(x)=1 ∂φ(x)=1；当 x = 0 x=0 x=0时， ∂ φ = [ − 1 , 1 ] \partial \varphi=[-1,1] ∂φ=[−1,1]这个闭区间。
1.4 ρ − c o n v e x \rho -convex ρ−convex函数
如果定义在凸集
C
C
C上的函数
φ
\varphi
φ对于任意
x
,
y
∈
C
,
λ
∈
[
0
,
1
]
x,y \in C,\lambda\in[0,1]
x,y∈C,λ∈[0,1]，对于某个
ρ
≥
0
\rho \geq0
ρ≥0满足以下不等式
φ
(
λ
x
+
(
1
−
λ
)
y
)
≤
λ
φ
(
x
)
+
(
1
−
λ
)
φ
(
y
)
−
ρ
2
λ
(
1
−
λ
)
∣
∣
x
−
y
∣
∣
2
\varphi(\lambda x+(1-\lambda)y)\leq\lambda \varphi(x)+(1-\lambda)\varphi(y)-\frac{\rho}{2}\lambda(1-\lambda)||x-y||^2
φ(λx+(1−λ)y)≤λφ(x)+(1−λ)φ(y)−2ρλ(1−λ)∣∣x−y∣∣2
那么称
φ
\varphi
φ为凸集
C
C
C上的
ρ
−
c
o
n
v
e
x
\rho -convex
ρ−convex函数。
1.5局部最优条件
原DC规划的充要局部最优条件为：（对偶DC规划问题类似，就不提了）
∂
h
(
x
∗
)
⋂
∂
g
(
x
∗
)
≠
∅
\partial h(x^*) \bigcap \partial g(x^*) \neq \empty
∂h(x∗)⋂∂g(x∗)=∅
这个点就称作临界点，或者说上式是原规划问题的一般KKT条件。并且
∅
≠
∂
h
(
x
∗
)
⊂
∂
g
(
x
∗
)
\empty \neq \partial h(x^*)\subset\partial g(x^*)
∅=∂h(x∗)⊂∂g(x∗)
在很多重要的DC规划问题中也是局部最优的充分条件。
1.6多面凸函数(polyhedral function)
如果一个函数是有限个仿射函数集合中的最大值组成的，那么这个函数一定是多面凸的。也就是说，共轭函数就是多面凸函数。
多面凸差规划是指凸差问题中的函数 g g g或 h h h至少有一个是多面凸的凸差问题。多面凸规划问题在非凸优化和全局优化上发挥了核心作用，是凸差规划和凸差算法的基础，从理论和算法的角度看在局部最优条件和DCA算法的收敛限度上都有非常让人感兴趣的特性，有如下两点：
- 在 h h h是多面凸函数的多面凸问题中，如果 h h h在临界点 x ∗ x^* x∗处可导，那么 x ∗ x^* x∗实际上就是原DC问题的局部最小值。实际上多面凸函数一般处处可导，所以找到了临界点一般就等价于找到了原DC问题的局部极小值。
- 函数 f f f在临界点 x ∗ x^* x∗处是局部凸的，所以就可以用解决凸优化问题的方法来解决原DC问题。这一点的原因是因为 f = g − h f=g-h f=g−h在 h h h可导的地方都是局部凸的。
二、凸差算法(DC Algorithm)
2.1迭代格式
基于凸差规划中的局部最优条件和对偶性，凸差算法(DCA)由分别构建原问题和对偶问题的尝试解的两个序列
x
k
{x_k}
xk和
y
k
y_k
yk组成，这两个序列满足
g
(
x
k
)
−
h
(
x
k
)
g(x^k)-h(x^k)
g(xk)−h(xk)和
h
∗
(
y
k
)
−
g
∗
(
y
k
)
h^*(y^k)-g^*(y^k)
h∗(yk)−g∗(yk)是递减的，并且都收敛至局部最优条件，以及
x
∗
∈
∂
g
∗
(
y
∗
)
,
y
∗
∈
∂
h
(
x
∗
)
x^* \in \partial g^*(y^*),y^* \in \partial h(x^*)
x∗∈∂g∗(y∗),y∗∈∂h(x∗)
这两个序列是分别以
x
k
+
1
x^{k+1}
xk+1、
y
k
+
1
y^{k+1}
yk+1是凸规化问题
P
k
P_k
Pk和
D
k
+
1
D_{k+1}
Dk+1的解的形式计算的，初始点
x
0
∈
d
o
m
∂
h
x^0 \in dom \partial h
x0∈dom∂h，
y
0
∈
∂
h
(
x
0
)
y^0 \in \partial h(x_0)
y0∈∂h(x0)。
(
P
k
)
i
n
f
{
g
(
x
)
−
[
h
(
x
k
)
+
<
x
−
x
k
,
y
k
>
]
:
x
∈
R
n
}
(P_k) \quad inf\{g(x)-[h(x^k)+<x-x^k,y^k>]:x \in \mathbb{R}^n\}
(Pk)inf{g(x)−[h(xk)+<x−xk,yk>]:x∈Rn}
( D k + 1 ) i n f { h ∗ ( y ) − [ g ∗ ( y k ) + < y − y k , x k + 1 > ] : y ∈ R n } (D_{k+1})\quad inf\{h^*(y)-[g^*(y^k)+<y-y^k,x^{k+1}>]:y \in \mathbb{R}^n\} (Dk+1)inf{h∗(y)−[g∗(yk)+<y−yk,xk+1>]:y∈Rn}
DCA算法有一个非常简单的理解：在第k次迭代中，我们把初始凸差问题中的第二个部分
h
h
h替换成它的仿射较小值(affine minorization)
h
(
x
k
)
+
<
x
−
x
k
,
y
k
>
h(x^k)+<x-x^k,y^k>
h(xk)+<x−xk,yk>。这个仿射较小值是通过函数h在
x
k
x^k
xk点处的次梯度
y
k
y^k
yk来定义的，目的是为了生成第
k
k
k个原始凸规化
(
P
k
)
(P_k)
(Pk)，这个规划问题的解恰恰就是
∂
g
∗
(
y
∗
)
\partial g^*(y^*)
∂g∗(y∗)。对偶地，
(
P
k
)
(P_k)
(Pk)的一个解
x
k
+
1
x^{k+1}
xk+1进而就被用来生成凸规化问题
(
D
k
+
1
)
(D_{k+1})
(Dk+1)——通过把对偶规划问题
(
D
d
c
)
(D_{dc})
(Ddc)中的第二项
g
∗
g^*
g∗换成其仿射最小值
(
g
∗
)
(
k
)
(
y
)
:
=
g
∗
(
y
k
)
+
<
y
−
y
k
,
x
k
+
1
>
(g^*)^{(k)}(y):=g^*(y^k)+<y-y^k,x^{k+1}>
(g∗)(k)(y):=g∗(yk)+<y−yk,xk+1>。这个对偶规划问题的解集也恰恰就是
h
h
h函数在点
x
k
+
1
x^{k+1}
xk+1处的微分
∂
h
(
x
k
+
1
)
\partial h(x^{k+1})
∂h(xk+1)。这个过程一直重复直至收敛。DCA借助函数
h
h
h和
g
∗
g*
g∗的次梯度执行了一次双重线性化，从而生成下一个迭代凸规化问题：
y
k
∈
∂
h
(
x
k
)
;
x
k
+
1
∈
∂
g
∗
(
y
k
)
,
∀
k
≥
0.
s
t
a
r
t
i
n
g
f
r
o
m
x
0
∈
d
o
m
∂
h
y^k \in \partial h(x^k);\quad x^{k+1} \in \partial g^*(y^k),\forall k\geq0.\\ starting \quad from\quad x^0\in dom \partial h
yk∈∂h(xk);xk+1∈∂g∗(yk),∀k≥0.startingfromx0∈dom∂h
这就是标准的DCA算法迭代格式，中心思想就是把第二项通过线性化替换成局部极小值，然后根据多面凸规化的特性，原问题（对偶问题）就变成了标准的凸优化问题。
2.2收敛特性
DCA是一个不带线性搜索的下降算法，但是具有全局收敛性，享有下面的特点：
(C和D是两个分别包含序列 { x k } \{x^k\} {xk} { y k } \{y^k\} {yk}在 { R n } \{\mathbb R^n\} {Rn}的凸集)
-
序列 { g ( x k ) − h ( x k ) } \{g(x^k)-h(x^k)\} {g(xk)−h(xk)}和 { h ∗ ( y k ) − g ∗ ( y k ) } \{h^*(y^k)-g^*(y^k)\} {h∗(yk)−g∗(yk)}是递减的并且
-
g ( x k + 1 ) − h ( x k + 1 ) = g ( x k ) − h ( x k ) g(x^{k+1})-h(x^{k+1})=g(x^k)-h(x^k) g(xk+1)−h(xk+1)=g(xk)−h(xk)当且仅当 y k ∈ ∂ g ( x k ) ⋂ ∂ h ( x k ) y^k \in \partial g(x^k)\bigcap \partial h(x^k) yk∈∂g(xk)⋂∂h(xk), y k ∈ ∂ g ( x k + 1 ) ⋂ ∂ h ( x k + 1 ) y^k \in \partial g(x^{k+1})\bigcap \partial h(x^{k+1}) yk∈∂g(xk+1)⋂∂h(xk+1)，并且 [ ρ ( g , C ) + ρ ( h , C ) ] ∣ ∣ x k + 1 − x k ∣ ∣ = 0 [\rho(g,C)+\rho(h,C)]||x^{k+1}-x^k||=0 [ρ(g,C)+ρ(h,C)]∣∣xk+1−xk∣∣=0.而且如果 g g g和 h h h在 C C C上是严格凸的，那么 x k = x k + 1 x^k=x^{k+1} xk=xk+1.
此时DCA算法就在第 k k k次迭代终止了(DCA算法的有限收敛)
-
对偶地有， h ∗ ( y k + 1 ) − g ∗ ( y k + 1 ) = h ∗ ( y k ) − g ∗ ( y k ) h^*(y^{k+1})-g^*(y^{k+1})=h^*(y^k)-g^*(y^k) h∗(yk+1)−g∗(yk+1)=h∗(yk)−g∗(yk)当且仅当 x k + 1 ∈ ∂ g ∗ ( x k ) ⋂ ∂ h ∗ ( y k ) x^{k+1} \in \partial g^*(x^k)\bigcap \partial h^*(y^k) xk+1∈∂g∗(xk)⋂∂h∗(yk), x k + 1 ∈ ∂ g ∗ ( y k + 1 ) ⋂ ∂ h ∗ ( y k + 1 ) x^{k+1} \in \partial g^*(y^{k+1})\bigcap \partial h^*(y^{k+1}) xk+1∈∂g∗(yk+1)⋂∂h∗(yk+1)，并且 [ ρ ( g ∗ , D ) + ρ ( h ∗ , D ) ] ∣ ∣ y k + 1 − y k ∣ ∣ = 0 [\rho(g^*,D)+\rho(h^*,D)]||y^{k+1}-y^k||=0 [ρ(g∗,D)+ρ(h∗,D)]∣∣yk+1−yk∣∣=0.而且如果 g ∗ g^* g∗和 h ∗ h^* h∗在 D D D上是严格凸的，那么 y k + 1 = y k y^{k+1}=y^k yk+1=yk。
此时DCA算法就在第 k k k次迭代终止了(DCA算法的有限收敛)
-
-
如果 ρ ( g , C ) + ρ ( h , C ) > 0 \rho(g,C)+\rho(h,C)>0 ρ(g,C)+ρ(h,C)>0(对应地， ρ ( g ∗ , D ) + ρ ( h ∗ , D ) > 0 \rho(g^*,D)+\rho(h^*,D)>0 ρ(g∗,D)+ρ(h∗,D)>0)),那么序列 ∣ ∣ x k + 1 − x k ∣ ∣ 2 {||x^{k+1}-x^k||^2} ∣∣xk+1−xk∣∣2（对应地, ∣ ∣ y k + 1 − y k ∣ ∣ ||y^{k+1}-y^k|| ∣∣yk+1−yk∣∣)收敛。
-
如果原DC问题的最优值是有限个的并且无穷序列 x k x^k xk和 y k y^k yk是有界的，那么每个无穷序列 x k x^k xk和 y k y^k yk的极限点 x ~ \tilde x x~(对应地, y ~ \tilde y y~)都是 g − h g-h g−h的临界点(对应地, h ∗ − g ∗ h^*-g^* h∗−g∗)。
-
对于一般地凸差规划，DCA线性收敛。
-
对于多面凸差问题，DCA有限收敛。
2.3相关问题
凸差分解，DCA中凸规化的求解，DCA初始点的选取，DCA多启动，凸差问题的最近分解技术等都是与凸差规划相关的问题，想进一步研究的可以去阅读参考文献[2].
2.4算法流程
算法：DCA
输入:
T:最大迭代次数
ϵ
\epsilon
ϵ:收敛误差
算法流程：
1.令t=0,初始化
x
t
x_t
xt
2.while (t<T) do
3. 计算
y
t
∈
∂
h
(
x
t
)
y_t \in \partial h(x_t)
yt∈∂h(xt);
4. 求解凸优化问题
a
r
g
m
i
n
{
g
(
x
)
−
[
h
(
x
t
)
+
<
x
−
x
t
,
y
t
>
]
}
argmin \{ g(x)-[h(x_t)+<x-x_t,y_t>]\}
argmin{g(x)−[h(xt)+<x−xt,yt>]}
得到
x
t
+
1
x_{t+1}
xt+1
5. if
∣
x
t
+
1
−
x
t
∣
≤
ϵ
|x_{t+1}-x_t|\leq\epsilon
∣xt+1−xt∣≤ϵ then
6. 算法收敛，跳出循环
7. end if
8. 令t=t+1
9.end while
10.return
x
t
x_t
xt
有一个小细节值得一提的是，步骤4中的凸优化问题中的一项 h ( x t ) h(x_t) h(xt)一项是一个常数，实际算的时候可以省去。
三、小案例
考虑以下非凸问题：
i
n
f
{
f
(
x
)
=
x
4
−
3
x
2
−
x
,
x
∈
R
}
inf\{f(x)=x^4-3x^2-x,x\in \mathbb R\}
inf{f(x)=x4−3x2−x,x∈R}
将其DC分解为两个凸函数：
f
(
x
)
=
g
(
x
)
−
h
(
x
)
w
h
e
r
e
g
(
x
)
=
x
4
,
h
(
x
)
=
3
x
2
+
x
f(x)=g(x)-h(x)\\ where\quad g(x)=x^4,h(x)=3x^2+x
f(x)=g(x)−h(x)whereg(x)=x4,h(x)=3x2+x
x0=input('请输入初始点位置\n');
ah0=6*x0+1; %h的导数*
Iter=10;it=1;
x=zeros(Iter,1);
y=zeros(Iter,1);
x(1)=x0;y(1)=ah0;
while it<=Iter
%求取g的导数
[x(it+1),~]=fminunc(@(a)a^4-(a-x(it))*y(it),0);*
y(it+1)=6*x(it+1)+1;
it=it+1;
end
disp([x,y])
当选取初始点 x = 0 x=0 x=0时，序列 x k x_k xk经过9次迭代收敛至1.3008：x=[ 0 0.6300 1.0612 1.2258 1.2783 1.2941 1.2989 1.3003 1.3007 1.3008 1.3008]。在 x ∗ = 1.3008 x^*=1.3008 x∗=1.3008这个临界点有 ∂ h ( x ∗ ) = 6 x ∗ + 1 = ∂ g ( x ∗ ) = 4 x ∗ 3 = 8.8048 \partial h(x*)=6x^*+1=\partial g(x^*)=4x^{*3}=8.8048 ∂h(x∗)=6x∗+1=∂g(x∗)=4x∗3=8.8048,也就是满足了凸差规划问题的局部最优条件。
参考文献
[1]Convex Analysis[2004],Stephen Boyd,et al
[2]Recent Advances in DC Programming and DCA [2014],Tao Pham Dinh
[3]凸优化[2013]，王书宁等译
[4]走进中神通Fenchel
3677
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
