有了前面标准Lagrangian乘子法与对偶上升法和增广Lagrangian法的基础,理解ADMM就容易了很多。本文主要来自张贤达《矩阵分析与优化(第二版)》4.7.4节。
ADMM算法
ADMM认为,在统计学与机器学习中,经常会遇到大尺度的等式约束优化问题,即
x
∈
R
n
x\in \mathbb{R}^n
x∈Rn的维数
n
n
n很大。如果
x
x
x可以分解为几个子向量,即
x
=
(
x
1
,
⋯
,
x
r
)
x=(x_1,\cdots,x_r)
x=(x1,⋯,xr),其目标函数也可以分解为:
f
(
x
)
=
∑
i
=
1
r
f
i
(
x
)
x
i
∈
R
n
i
,
∑
i
=
1
r
n
i
=
n
f(x)=\sum_{i=1}^r f_i(x) \\ x_i\in \mathbb{R}^{n_i},\sum_{i=1}^r n_i=n
f(x)=i=1∑rfi(x)xi∈Rni,i=1∑rni=n
则大尺度的优化问题可以转化为分布式优化问题。相应的,等式约束矩阵
A
x
=
b
Ax=b
Ax=b也分块为:
A
=
[
A
1
,
⋯
,
A
r
]
,
A
x
=
∑
i
=
1
r
A
i
x
i
=
b
A=[A_1,\cdots,A_r], Ax=\sum_{i=1}^r A_ix_i=b
A=[A1,⋯,Ar],Ax=i=1∑rAixi=b
于是增广Lagrangian目标函数
L
ρ
(
x
,
λ
)
L_\rho(x,\lambda)
Lρ(x,λ)可以写作:
L
ρ
(
x
,
λ
)
=
∑
i
=
1
r
(
f
i
(
x
i
)
+
λ
T
A
i
x
i
)
−
λ
T
b
+
ρ
2
∥
∑
i
=
1
r
(
A
i
x
i
)
−
b
∥
2
2
L_\rho(x,\lambda)=\sum_{i=1}^r(f_i(x_i)+\lambda^TA_ix_i)-\lambda^Tb+\frac{\rho}{2}\|\sum_{i=1}^r(A_ix_i)-b\|_2^2
Lρ(x,λ)=i=1∑r(fi(xi)+λTAixi)−λTb+2ρ∥i=1∑r(Aixi)−b∥22
所取的罚函数与增广Lagrangian乘子法中的仍相同。再采用对偶上升法,即可得到能进行并行运算的分散算法:
x
i
k
+
1
=
arg min
x
i
∈
R
n
i
L
i
(
x
i
,
λ
k
)
,
i
=
1
,
⋯
,
r
λ
k
+
1
=
λ
k
+
ρ
k
(
∑
i
=
1
r
A
i
x
i
k
+
1
−
b
)
x_i^{k+1}=\argmin_{x_i\in \mathbb{R}^{n_i}} L_i(x_i,\lambda_k),i=1,\cdots,r \\ \lambda_{k+1}=\lambda_k+\rho_k (\sum_{i=1}^r A_i x_i^{k+1}-b)
xik+1=xi∈RniargminLi(xi,λk),i=1,⋯,rλk+1=λk+ρk(i=1∑rAixik+1−b)
这里的
x
i
x_i
xi是可以独立更新的。由于
x
i
x_i
xi以一种交替的或序贯的方式进行更新,所以称为“交替方向”乘子法(ADMM算法)。
举个 r = 2 r=2 r=2的例子
r
=
2
r=2
r=2,则目标函数为:
min
f
(
x
)
+
g
(
z
)
s
.
t
.
A
x
+
B
z
=
c
\min f(x)+g(z)\\ s.t.\ Ax+Bz=c
minf(x)+g(z)s.t. Ax+Bz=c
上式中,
x
∈
R
n
,
z
∈
R
m
,
A
∈
R
p
×
n
,
B
∈
R
p
×
m
,
c
∈
R
p
x\in \mathbb{R}^n,z\in \mathbb{R}^m,A\in \mathbb{R}^{p\times n},B\in \mathbb{R}^{p\times m},c\in \mathbb{R}^{p}
x∈Rn,z∈Rm,A∈Rp×n,B∈Rp×m,c∈Rp。则增广Lagrangian目标函数为:
L
ρ
(
x
,
z
,
λ
)
=
f
(
x
)
+
g
(
z
)
+
λ
T
(
A
x
+
B
z
−
c
)
+
ρ
2
∥
A
x
+
B
z
−
c
∥
2
2
(1)
L_\rho(x,z,\lambda)=f(x)+g(z)+\lambda^T(Ax+Bz-c)+\frac{\rho}{2}\|Ax+Bz-c\|_2^2 \\ \tag{1}
Lρ(x,z,λ)=f(x)+g(z)+λT(Ax+Bz−c)+2ρ∥Ax+Bz−c∥22(1)
上式的交替方向乘子法的更新公式为:
x
k
+
1
=
arg min
x
∈
R
n
L
ρ
(
x
,
z
k
,
λ
k
)
z
k
+
1
=
arg min
z
∈
R
m
L
ρ
(
x
k
+
1
,
z
,
λ
k
)
λ
k
+
1
=
λ
k
+
ρ
k
(
A
x
k
+
1
+
B
z
k
+
1
−
c
)
x_{k+1}=\argmin_{x\in \mathbb{R}^n} L_\rho(x,z_k,\lambda_k)\\ z_{k+1}=\argmin_{z\in \mathbb{R}^m} L_\rho(x_{k+1},z,\lambda_k)\\ \lambda_{k+1}=\lambda_k+\rho_k(Ax_{k+1}+Bz_{k+1}-c)
xk+1=x∈RnargminLρ(x,zk,λk)zk+1=z∈RmargminLρ(xk+1,z,λk)λk+1=λk+ρk(Axk+1+Bzk+1−c)
误差分析与停止条件
公式
(
1
)
(1)
(1)的最优化条件分为原始可行性:
A
x
+
B
z
−
c
=
0
Ax+Bz-c=0
Ax+Bz−c=0
和对偶可行性:
0
∈
∂
f
(
x
)
+
A
T
λ
+
ρ
A
T
(
A
x
+
B
z
−
c
)
=
∂
f
(
x
)
+
A
T
λ
0
∈
∂
f
(
x
)
+
B
T
λ
+
ρ
B
T
(
A
x
+
B
z
−
c
)
=
∂
g
(
z
)
+
B
T
λ
(2)
0\in \partial f(x)+A^T\lambda+\rho A^T(Ax+Bz-c)=\partial f(x)+A^T\lambda \\ 0\in \partial f(x)+B^T\lambda+\rho B^T(Ax+Bz-c)=\partial g(z)+B^T\lambda \\ \tag{2}
0∈∂f(x)+ATλ+ρAT(Ax+Bz−c)=∂f(x)+ATλ0∈∂f(x)+BTλ+ρBT(Ax+Bz−c)=∂g(z)+BTλ(2)
根据我的企业级理解,这原始和对偶可行性关系分别是等式约束成立和偏导为0,是从KKT条件来的,都是必要条件。不过书里没有明确指出是或者不是。
关于推导,这里用的是
0
∈
0\in
0∈而不是
0
=
0=
0=,这是什么企业级逻辑我没弄懂,不过我觉得不影响理解,意思差不多。书上这里求导有问题,疑似纰漏,我改成了公式
(
1
)
(1)
(1)的正确的求导结果,这样也和后文更对的上。要记得,公式
(
1
)
(1)
(1)中向量
A
x
+
B
z
−
c
Ax+Bz-c
Ax+Bz−c二范数平方的一阶导等于向量的二倍,再乘一个系数
A
T
A^T
AT,就可以得到这个结果。求导法则可以参考我的这篇总结。再加上
A
x
+
B
z
−
c
=
0
Ax+Bz-c=0
Ax+Bz−c=0的约束,就能推导下来了。
在迭代的过程中,原始可行性不可能完全满足,设其误差为:
r
k
=
A
x
k
+
B
z
k
−
c
r_k=Ax_k+Bz_k-c
rk=Axk+Bzk−c
称为第
k
k
k次迭代的原始残差(向量),这样Lagrangian乘子向量的更新可以用这个残差重写为:
λ
k
+
1
=
λ
k
+
ρ
k
r
k
+
1
\lambda_{k+1}=\lambda_k+\rho_k r_{k+1}
λk+1=λk+ρkrk+1
同样,对偶可行性也不会完全满足:
0
∈
∂
f
(
x
k
+
1
)
+
A
T
λ
k
+
ρ
A
T
(
A
x
k
+
1
+
B
z
k
−
c
)
=
∂
f
(
x
k
+
1
)
+
A
T
[
λ
k
+
ρ
(
A
x
k
+
1
+
B
z
k
+
1
−
c
)
+
ρ
B
(
z
k
−
z
k
+
1
)
]
=
∂
f
(
x
k
+
1
)
+
A
T
[
λ
k
+
ρ
r
k
+
1
+
ρ
B
(
z
k
−
z
k
+
1
)
]
=
∂
f
(
x
k
+
1
)
+
A
T
λ
k
+
1
+
ρ
A
T
B
(
z
k
−
z
k
+
1
)
0\in \partial f(x_{k+1})+A^T\lambda_k+\rho A^T(Ax_{k+1}+Bz_k-c) \\ =\partial f(x_{k+1})+A^T[\lambda_k+\rho (Ax_{k+1}+Bz_{k+1}-c)+\rho B(z_k-z_{k+1})] \\ =\partial f(x_{k+1})+A^T[\lambda_k+\rho r_{k+1}+\rho B(z_k-z_{k+1})]\\ =\partial f(x_{k+1}) +A^T\lambda_{k+1}+\rho A^TB(z_k-z_{k+1})
0∈∂f(xk+1)+ATλk+ρAT(Axk+1+Bzk−c)=∂f(xk+1)+AT[λk+ρ(Axk+1+Bzk+1−c)+ρB(zk−zk+1)]=∂f(xk+1)+AT[λk+ρrk+1+ρB(zk−zk+1)]=∂f(xk+1)+ATλk+1+ρATB(zk−zk+1)
注意,由于书上公式
(
2
)
(2)
(2)求导是错的,所以这一步更别扭,怎么看都不对劲,这里我也改成了我认为的正确的推导形式。
对照公式
(
2
)
(2)
(2)中的第一个式子可知对偶残差为:
s
k
+
1
=
ρ
A
T
B
(
z
k
−
z
k
+
1
)
s_{k+1}=\rho A^TB(z_k-z_{k+1})
sk+1=ρATB(zk−zk+1)
交替方向乘子法的停止条件就是两个残差都小于阈值:
∥
r
k
+
1
∥
2
≤
ε
p
r
i
a
n
d
∥
s
k
+
1
∥
≤
ε
d
u
a
l
\|r_{k+1}\|_2\le \varepsilon_{pri} \ and \ \|s_{k+1}\|\le \varepsilon_{dual}
∥rk+1∥2≤εpri and ∥sk+1∥≤εdual
缩放形式的ADMM
令
v
=
(
1
/
ρ
)
λ
v=(1/\rho)\lambda
v=(1/ρ)λ为经过
1
/
ρ
1/\rho
1/ρ缩放的Lagrangian乘子向量,则更新公式变为:
x
k
+
1
=
arg min
x
∈
R
n
L
ρ
(
x
,
z
k
,
v
k
)
z
k
+
1
=
arg min
z
∈
R
m
L
ρ
(
x
k
+
1
,
z
,
v
k
)
v
k
+
1
=
v
k
+
A
x
k
+
1
+
B
z
k
+
1
−
c
=
v
k
+
r
k
+
1
x_{k+1}=\argmin_{x\in \mathbb{R}^n} L_\rho(x,z_k,v_k)\\ z_{k+1}=\argmin_{z\in \mathbb{R}^m} L_\rho(x_{k+1},z,v_k)\\ v_{k+1}=v_k+Ax_{k+1}+Bz_{k+1}-c=v_k+r_{k+1}
xk+1=x∈RnargminLρ(x,zk,vk)zk+1=z∈RmargminLρ(xk+1,z,vk)vk+1=vk+Axk+1+Bzk+1−c=vk+rk+1
其第
k
k
k次迭代的残差
r
k
r_k
rk为:
r
k
=
A
x
k
+
B
z
k
−
c
=
v
0
+
∑
i
=
1
k
r
i
r_k=Ax_k+Bz_k-c=v_0+\sum_{i=1}^kr_i
rk=Axk+Bzk−c=v0+i=1∑kri
即,第
k
k
k次迭代的缩放对偶向量是所有
k
k
k次迭代的原始残差之和。这种方法称为有缩放的交替方向乘子法。
最后,个人认为,不论是对偶上升法,增广Lagrangian乘子法,还是ADMM算法,核心思想都相似,而且具体使用时都要与其他最优化方法结合,因为 arg min L ( x , z , λ ) \argmin L(x,z,\lambda) argminL(x,z,λ)的求解是还需要别的方法的,停止条件需要根据使用环境具体再去确定。