随机优化是机器学习中至关重要的组成部分,其核心是随机梯度下降算法(SGD),这种方法自60多年前首次提出以来一直被广泛使用。近八年来,我们见证了一个激动人心的新进展:随机优化方法的方差降低技术。这些方差降低的方法(VR方法)在允许多次迭代训练数据的场景下表现出色,无论是在理论上还是在实践中,它们都显示出比SGD更快的收敛速度。这种速度的提高凸显了VR方法的日益增长的兴趣和这一领域迅速积累的研究成果。本文综述了VR方法在有限数据集优化中的关键原则和主要进展,旨在为非专家读者提供信息。我们主要集中讨论了凸优化环境,并为那些对非凸函数最小化扩展感兴趣的读者提供了参考。
关键词 | 机器学习;优化;方差降低
1.引言
在机器学习的研究领域中,一个基础而重要的问题是如何将模型适配到庞大的数据集上。例如,我们可以考虑线性最小二乘模型的典型案例:
x ∗ ∈ arg min x ∈ R d 1 n ∑ i = 1 n ( a i T x − b i ) 2 x^* \in \arg\min_{x \in \mathbb{R}^d} \frac{1}{n} \sum_{i=1}^{n} (a_i^T x - b_i)^2 x∗∈argx∈Rdminn1i=1∑n(aiTx−bi)2
在这个模型中,我们有 d d d 个参数,它们由向量 x ∈ R d x \in \mathbb{R}^d x∈Rd 给出。同时,我们手头有 n n n 个数据点,包括特征向量 a i ∈ R d a_i \in \mathbb{R}^d ai∈Rd 和目标值 b i ∈ R b_i \in \mathbb{R} bi∈R。模型的适配过程就是调整这些参数,以使得模型的预测输出 a i T x a_i^T x aiTx 平均上尽可能接近目标值 b i b_i bi。
更广泛地说,我们可能会使用损失函数 f i ( x ) f_i(x) fi(x) 来衡量模型预测与第 i i i 个数据点的接近程度:
x ∗ ∈ arg min x ∈ R d f ( x ) : = 1 n ∑ i = 1 n f i ( x ) x^* \in \arg\min_{x \in \mathbb{R}^d} f(x) := \frac{1}{n} \sum_{i=1}^{n} f_i(x) x∗∈argx∈Rdminf(x):=n1i=1∑nfi(x)
损失函数 f i ( x ) f_i(x) fi(x) 如果较大,表明模型的预测与数据有较大偏差;如果 f i ( x ) f_i(x) fi(x) 等于零,则表示模型完美地拟合了数据点。函数 f ( x ) f(x) f(x) 反映了模型在整个数据集上的平均损失。
类似上述形式 (2) 的问题不仅适用于线性最小二乘问题,也适用于机器学习中研究的许多其他模型。例如,在逻辑回归模型中,我们解决的是:
x ∗ ∈ arg min x ∈ R d 1 n ∑ i = 1 n log ( 1 + e − b i a i T x ) + λ 2 ∥ x ∥ 2 2 x^* \in \arg\min_{x \in \mathbb{R}^d} \frac{1}{n} \sum_{i=1}^{n} \log(1 + e^{-b_i a_i^T x}) + \frac{\lambda}{2} \|x\|_2^2 x∗∈argx∈Rdminn1i=1∑nlog(1+e−biaiTx)+2λ∥x∥22
这里,我们处理的是 b i ∈ { − 1 , + 1 } b_i \in \{-1, +1\} bi∈{−1,+1} 的二元分类问题,预测是基于 a i T x a_i^T x aiTx 的符号来进行的。公式中还引入了正则化项 λ 2 ∥ x ∥ 2 2 \frac{\lambda}{2} \|x\|_2^2 2λ∥x∥22 来避免对数据的过拟合,其中 ∥ x ∥ 2 2 \|x\|_2^2 ∥x∥22 表示 x x x 的欧几里得范数的平方。
在大多数监督学习模型中,训练过程可以表示为形式 (2),包括 L1 正则化最小二乘、支持向量机 (SVM)、主成分分析、条件随机场和深度神经网络等。
现代问题实例中的一个关键挑战是数据点的数量
n
n
n 可能极其庞大。我们经常处理的数据集大小远远超出了太字节的范围,这些数据可能来自互联网、卫星、远程传感器、金融市场和科学实验等多种来源。为了应对如此庞大的数据集,一种常见的方法是使用随机梯度下降(SGD)算法,该算法在每次迭代中仅使用少量随机选取的数据点。此外,最近对方差减少(VR)的随机梯度方法的兴趣急剧上升,这些方法比传统的随机梯度方法具有更快的收敛速度。
图1. 在基于蘑菇数据集[7]的逻辑回归问题上,将梯度下降(GD)、加速梯度下降(AGD,即[50]中的加速GD)、随机梯度下降(SGD)和ADAM[30]方法与方差减少(VR)方法SAG和SVRG进行了比较,其中n=8124,d=112。
1.1. 梯度与随机梯度下降方法
梯度下降(GD)是一种经典算法,用于解决上述问题 (2),其迭代更新公式如下所示:
x
k
+
1
=
x
k
−
γ
1
n
∑
i
=
1
n
∇
f
i
(
x
k
)
x_{k+1} = x_k - \gamma \frac{1}{n} \sum_{i=1}^{n} \nabla f_i(x_k)
xk+1=xk−γn1i=1∑n∇fi(xk)
这里, γ \gamma γ 是一个大于零的固定步长值。在GD算法的每次迭代过程中,必须对每一个数据点 i i i 计算梯度 ∇ f i ( x k ) \nabla f_i(x_k) ∇fi(xk),这意味着GD需要对所有 n n n 个数据点进行完整的遍历。当数据集的大小 n n n 变得非常大时,GD算法的每次迭代成本会变得非常高,从而限制了其应用。
作为替代,我们可以考虑随机梯度下降(SGD)方法,这是由 Robbins 和 Monro 首次提出的,其迭代更新公式如下:
x
k
+
1
=
x
k
−
γ
∇
f
i
k
(
x
k
)
x_{k+1} = x_k - \gamma \nabla f_{i_k}(x_k)
xk+1=xk−γ∇fik(xk)
SGD算法通过在每次迭代中仅使用一个随机选取的数据点的梯度 ∇ f i k ( x k ) \nabla f_{i_k}(x_k) ∇fik(xk) 来降低每次迭代的成本。在图 1 中,我们可以看到SGD在优化过程的初期阶段比GD(包括加速的GD方法)取得了更显著的进步。该图根据轮次(epoch)来展示优化的进展,轮次定义为计算所有 n n n 个训练样本的梯度的次数。GD算法在每个轮次进行一次迭代,而SGD算法则在每个轮次进行 n n n 次迭代。我们以轮次作为比较SGD和GD的依据,因为在假设 n n n 非常大的情况下,两种方法的主要成本都集中在梯度 ∇ f i ( x k ) \nabla f_i(x_k) ∇fi(xk) 的计算上。
1.2.方差问题
让我们考虑随机索引
i
k
i_k
ik 从集合
{
1
,
…
,
n
}
\{1, \ldots, n\}
{1,…,n} 中均匀随机选择的情况,这意味着对于所有
i
i
i,选择
i
k
=
i
i_k = i
ik=i 的概率
P
[
i
k
=
i
]
P[i_k = i]
P[ik=i] 等于
1
n
\frac{1}{n}
n1。在这种情况下,
∇
f
i
k
(
x
k
)
\nabla f_{i_k}(x_k)
∇fik(xk) 作为
∇
f
(
x
k
)
\nabla f(x_k)
∇f(xk) 的估计量是无偏的,因为根据期望的定义,我们有:
E
[
∇
f
i
k
(
x
k
)
∣
x
k
]
=
1
n
∑
i
=
1
n
∇
f
i
(
x
k
)
=
∇
f
(
x
k
)
(
6
)
E[\nabla f_{i_k}(x_k) | x_k] = \frac{1}{n} \sum_{i=1}^{n} \nabla f_i(x_k) = \nabla f(x_k) \quad (6)
E[∇fik(xk)∣xk]=n1i=1∑n∇fi(xk)=∇f(xk)(6)
尽管SGD(随机梯度下降)方法在每次迭代中不保证函数 f f f 的值会减少,但平均而言,它朝着负的完整梯度方向移动,这代表了下降方向。
然而,拥有一个无偏梯度估计量并不足以确保SGD迭代的收敛性。为了说明这一点,图 2(左侧)展示了使用常数步长对LIBSVM [7] 提供的四类数据集应用逻辑回归函数时SGD的迭代轨迹。图中的同心椭圆代表了函数的等高线,即函数值 f ( x ) = c f(x) = c f(x)=c 对应的点 x x x 集合, c c c 是实数集中的特定常数。不同的常数值 c c c 对应不同的椭圆。
SGD的迭代轨迹并没有收敛到最优解(图中以绿色星号表示),而是在最优解周围形成了一个点云。与此相反,我们在图 2中使用相同的常数步长展示了一种方差减少(VR)方法——随机平均梯度(SAG)的迭代轨迹,我们将在后文中介绍这种方法。SGD在这个例子中未能收敛的原因是随机梯度本身没有收敛到零,因此,常数步长的SGD方法(5)永远不会停止。这与梯度下降(GD)方法形成鲜明对比,GD方法会自然停止,因为随着
x
k
x_k
xk 趋近于
x
∗
x^*
x∗,梯度
∇
f
(
x
k
)
\nabla f(x_k)
∇f(xk) 会趋向于零。
图2. 使用固定步长的SGD(左)和SAG(右)迭代方法的二维逻辑回归的水平集图。绿色星号表示x解。
1.3.经典方差减少方法
处理由于 ∇ f i ( x k ) \nabla f_i(x_k) ∇fi(xk) 值的方差导致的非收敛性问题,有几种经典技术。例如,Robbins 和 Monro [64] 通过使用一系列递减的步长 γ k \gamma_k γk 来解决方差问题,确保乘积 γ k ∇ f i k ( x k ) \gamma_k \nabla f_{i_k}(x_k) γk∇fik(xk) 能够收敛到零。然而,调整这个递减步长序列以避免过早或过晚停止算法是一个难题。
另一种减少方差的经典技术是每次迭代中使用多个
∇
f
i
(
x
k
)
\nabla f_i(x_k)
∇fi(xk) 的平均值,以获得对完整梯度
∇
f
(
x
)
\nabla f(x)
∇f(x) 的更准确估计。这种方法称为小批量处理(minibatch),尤其适用于可以并行评估多个梯度的情况。这导致迭代形式如下:
x
k
+
1
=
x
k
−
γ
1
∣
B
k
∣
∑
i
∈
B
k
∇
f
i
(
x
k
)
(
7
)
x_{k+1} = x_k - \gamma \frac{1}{|B_k|} \sum_{i \in B_k} \nabla f_i(x_k) \quad (7)
xk+1=xk−γ∣Bk∣1i∈Bk∑∇fi(xk)(7)
其中
B
k
B_k
Bk 是一个随机索引集合,
∣
B
k
∣
|B_k|
∣Bk∣ 表示
B
k
B_k
Bk 的大小。如果
B
k
B_k
Bk 以有放回的方式均匀采样,那么这个梯度估计的方差与“批量大小”
∣
B
k
∣
|B_k|
∣Bk∣ 成反比,因此可以通过增加批量大小来降低方差。
但是,这种迭代的成本与批量大小成正比,因此这种方差减少形式是以增加计算成本为代价的。
另一个常见的减少方差并提高SGD经验性能的策略是添加“动量”,这是一个基于过去步骤中使用的方向的额外项。特别是,带动量的SGD形式如下:
x
k
+
1
=
x
k
−
γ
m
k
(
9
)
x_{k+1} = x_k - \gamma m_k \quad (9)
xk+1=xk−γmk(9)
其中动量参数
β
\beta
β 位于 (0, 1) 范围内。如果初始动量
m
0
=
0
m_0 = 0
m0=0,并在 (8) 中展开
m
k
m_k
mk 的更新,我们得到
m
k
m_k
mk 是之前梯度的加权平均:
m
k
=
∑
t
=
0
k
β
k
−
t
∇
f
i
t
(
x
t
)
(
10
)
m_k = \sum_{t=0}^{k} \beta^{k-t} \nabla f_{i_t}(x_t) \quad (10)
mk=t=0∑kβk−t∇fit(xt)(10)
因此,
m
k
m_k
mk 是随机梯度的加权和。由于
∑
t
=
0
k
β
k
−
t
=
1
−
β
k
+
1
1
−
β
\sum_{t=0}^{k} \beta^{k-t} = \frac{1 - \beta^{k+1}}{1 - \beta}
∑t=0kβk−t=1−β1−βk+1,我们可以将
1
−
β
1
−
β
k
m
k
\frac{1 - \beta}{1 - \beta^k} m_k
1−βk1−βmk 视为随机梯度的加权平均。如果我们将其与完整梯度的表达式
∇
f
(
x
k
)
=
1
n
∑
i
=
1
n
∇
f
i
(
x
k
)
\nabla f(x_k) = \frac{1}{n} \sum_{i=1}^{n} \nabla f_i(x_k)
∇f(xk)=n1∑i=1n∇fi(xk) 进行比较,我们可以将
1
−
β
1
−
β
k
m
k
\frac{1 - \beta}{1 - \beta^k} m_k
1−βk1−βmk(以及
m
k
m_k
mk)解释为对完整梯度的估计。这种加权和虽然减少了方差,但也带来了关键问题。由于加权和(10)对最近采样的梯度赋予了更多权重,它不会收敛到完整梯度
∇
f
(
x
k
)
\nabla f(x_k)
∇f(xk),后者是一个简单平均。我们将在第二节A中看到的第一种方差减少方法通过使用简单平均而不是任何加权平均来解决这个问题。
1.4.现代方差减少方法
与经典方法不同,它们直接使用一个或多个
∇
f
i
(
x
k
)
\nabla f_i(x_k)
∇fi(xk) 作为
∇
f
(
x
k
)
\nabla f(x_k)
∇f(xk) 的近似值,现代方差减少(VR)方法采用了一种不同的策略。这些方法利用
∇
f
i
(
x
k
)
\nabla f_i(x_k)
∇fi(xk) 来更新梯度的估计值
g
k
g_k
gk,其目标是让
g
k
g_k
gk 逼近
∇
f
(
x
k
)
\nabla f(x_k)
∇f(xk)。具体来说,我们希望
g
k
g_k
gk 能够满足
g
k
≈
∇
f
(
x
k
)
g_k \approx \nabla f(x_k)
gk≈∇f(xk)。基于这样的梯度估计,我们接着执行形式如下的近似梯度步骤:
x
k
+
1
=
x
k
−
γ
g
k
(
11
)
x_{k+1} = x_k - \gamma g_k \quad (11)
xk+1=xk−γgk(11)
这里的
γ
>
0
\gamma > 0
γ>0 是步长参数。
为了确保使用常数步长
γ
\gamma
γ 时迭代 (11) 能够收敛,我们需要保证梯度估计
g
k
g_k
gk 的方差趋向于零。数学上,这可以表达为:
E
[
∥
g
k
−
∇
f
(
x
k
)
∥
2
]
→
0
as
k
→
∞
(
12
)
E\left[ \| g_k - \nabla f(x_k) \|^2 \right] \rightarrow 0 \quad \text{as } k \rightarrow \infty \quad (12)
E[∥gk−∇f(xk)∥2]→0as k→∞(12)
这里的期望
E
E
E 是基于算法中直到第
k
k
k 次迭代的所有随机变量计算的。属性 (12) 确保了 VR 方法在达到最优解时能够停止。我们将此属性视为 VR 方法的一个标志性特征,因此称之为 VR 属性。值得注意的是,“减少的”方差这个表述可能会引起误解,因为实际上方差是趋向于零的。属性 (12) 是 VR 方法在理论上(在适当的假设条件下)和实践中(如图 1 所展示的)能够实现更快收敛的关键因素。
1.5.第一个方差减少方法的例子:SGD²
一种简单的改进方法可以使SGD递归式(5)在不减小步长的情况下实现收敛,那就是对每个梯度进行平移,具体做法是减去
∇
f
i
(
x
∗
)
\nabla f_i(x^*)
∇fi(x∗),这种方法定义如下:
x
k
+
1
=
x
k
−
γ
(
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
∗
)
)
(
13
)
x_{k+1} = x_k - \gamma (\nabla f_{i_k}(x_k) - \nabla f_{i_k}(x^*)) \quad (13)
xk+1=xk−γ(∇fik(xk)−∇fik(x∗))(13)
这种方法被称为SGD² [22]。虽然我们通常无法确切知道每个
∇
f
i
(
x
∗
)
\nabla f_i(x^*)
∇fi(x∗),但SGD²作为一个例子,能够很好地阐释方差减少方法的基本特性。此外,许多方差减少方法都可以看作是SGD²方法的一种近似形式;这些方法不是依赖于已知的每个
∇
f
i
(
x
∗
)
\nabla f_i(x^*)
∇fi(x∗),而是使用能够逼近
∇
f
i
(
x
∗
)
\nabla f_i(x^*)
∇fi(x∗) 的估计值。
值得注意的是,SGD²使用的是对完整梯度的无偏估计。因为
∇
f
(
x
∗
)
=
0
\nabla f(x^*) = 0
∇f(x∗)=0,所以有:
E
[
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
∗
)
]
=
∇
f
(
x
k
)
−
∇
f
(
x
∗
)
=
∇
f
(
x
k
)
E[\nabla f_{i_k}(x_k) - \nabla f_{i_k}(x^*)] = \nabla f(x_k) - \nabla f(x^*) = \nabla f(x_k)
E[∇fik(xk)−∇fik(x∗)]=∇f(xk)−∇f(x∗)=∇f(xk)
另外,当SGD²达到最优解时,它自然会停止,因为对于任意
i
i
i,有:
(
∇
f
i
(
x
)
−
∇
f
i
(
x
∗
)
)
∣
x
=
x
∗
=
0
(\nabla f_i(x) - \nabla f_i(x^*)) \bigg|_{x=x^*} = 0
(∇fi(x)−∇fi(x∗))
x=x∗=0
进一步观察,随着
x
k
x_k
xk 接近
x
∗
x^*
x∗(对于连续的
∇
f
i
\nabla f_i
∇fi),SGD²满足方差减少属性(12),因为:
E
[
∥
g
k
−
∇
f
(
x
k
)
∥
2
]
=
E
[
∥
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
∗
)
−
∇
f
(
x
k
)
∥
2
]
≤
E
[
∥
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
∗
)
∥
2
]
E\left[ \| g_k - \nabla f(x_k) \|^2 \right] = \\E\left[ \| \nabla f_{i_k}(x_k) - \nabla f_{i_k}(x^*) - \nabla f(x_k) \|^2 \right] \leq E\left[ \| \nabla f_{i_k}(x_k) - \nabla f_{i_k}(x^*) \|^2 \right]
E[∥gk−∇f(xk)∥2]=E[∥∇fik(xk)−∇fik(x∗)−∇f(xk)∥2]≤E[∥∇fik(xk)−∇fik(x∗)∥2]
这里我们使用了引理2,令
X
=
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
∗
)
X = \nabla f_{i_k}(x_k) - \nabla f_{i_k}(x^*)
X=∇fik(xk)−∇fik(x∗),并利用了
E
[
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
∗
)
]
=
∇
f
(
x
k
)
E[\nabla f_{i_k}(x_k) - \nabla f_{i_k}(x^*)] = \nabla f(x_k)
E[∇fik(xk)−∇fik(x∗)]=∇f(xk) 的性质。这个属性表明SGD²具有比传统SGD方法更快的收敛速度,我们已在附录B中对此进行了详细说明。
1.6.方差减少方法的快速收敛性
本节我们将介绍两个标准假设,这些假设用于分析方差减少(VR)方法,并讨论在这些假设下相比于传统SGD方法所能够实现的加速效果。首先,我们假设梯度具有Lipschitz连续性,这表示梯度的变化速度是有限的。
假设1(Lipschitz连续性)
我们假设函数
f
f
f是可微的并且是
L
L
L-平滑的,对于所有
x
x
x 和
y
y
y 以及某个
0
<
L
<
∞
0 < L < \infty
0<L<∞,满足以下条件:
∥
∇
f
(
x
)
−
∇
f
(
y
)
∥
≤
L
∥
x
−
y
∥
(
14
)
\|\nabla f(x) - \nabla f(y)\| \leq L\|x - y\| \quad (14)
∥∇f(x)−∇f(y)∥≤L∥x−y∥(14)
这意味着每个
f
i
:
R
d
→
R
fi: \mathbb{R}^d \rightarrow \mathbb{R}
fi:Rd→R 是可微的,
L
i
L_i
Li-平滑的,我们定义
L
max
L_{\text{max}}
Lmax 为
max
{
L
1
,
.
.
.
,
L
n
}
\max\{L_1, . . . , L_n\}
max{L1,...,Ln}。
虽然这通常被认为是一个较弱的假设,但在后续章节中,我们将讨论适用于非光滑问题的VR方法。对于两次可微的单变量函数, L L L-平滑性可以直观理解为:它等同于假设二阶导数被 L L L 上限,即 ∣ f ′ ′ ( x ) ∣ ≤ L |f''(x)| \leq L ∣f′′(x)∣≤L 对于所有 x ∈ R d x \in \mathbb{R}^d x∈Rd。对于多变量的两次可微函数,它等同于假设Hessian矩阵 ∇ 2 f ( x ) \nabla^2 f(x) ∇2f(x) 的奇异值被 L L L 上限。
假设2(强凸性)
我们考虑的第二个假设是函数 (f) 是 μ \mu μ-强凸的,这意味着对于某个 μ > 0 \mu > 0 μ>0,函数 x ↦ f ( x ) − μ 2 ∥ x ∥ 2 x \mapsto f(x) - \frac{\mu}{2}\|x\|^2 x↦f(x)−2μ∥x∥2 是凸的。此外,对于每个 i = 1 , . . . , n i = 1, . . . , n i=1,...,n, f i : R d → R fi: \mathbb{R}^d \rightarrow \mathbb{R} fi:Rd→R 是凸的。
这是一个较强的假设。在最小二乘问题中,每个 (fi$ 是凸的,但总体函数 (f) 只有在设计矩阵 A : = [ a 1 , . . . , a n ] A := [a_1, . . . , a_n] A:=[a1,...,an] 具有完全行秩时才是强凸的。L2正则化的逻辑回归问题由于正则化项的存在,满足这个假设,其中 μ ≥ λ \mu \geq \lambda μ≥λ。
满足这些假设的一个重要问题类别是形式如下的优化问题:
x
∗
∈
arg
min
x
∈
R
d
f
(
x
)
=
1
n
∑
i
=
1
n
ℓ
i
(
a
i
T
x
)
+
λ
2
∥
x
∥
2
(
15
)
x^* \in \arg\min_{x \in \mathbb{R}^d} f(x) = \frac{1}{n} \sum_{i=1}^{n} \ell_i(a_i^Tx) + \frac{\lambda}{2}\|x\|^2 \quad (15)
x∗∈argx∈Rdminf(x)=n1i=1∑nℓi(aiTx)+2λ∥x∥2(15)
其中每个“损失”函数
ℓ
i
:
R
→
R
\ell_i: \mathbb{R} \rightarrow \mathbb{R}
ℓi:R→R 是两次可微的,并且其二阶导数
ℓ
i
′
′
\ell_i''
ℓi′′ 被限制在0和某个上界
M
M
M 之间。这包括了机器学习中带有L2正则化的多种损失函数,例如最小二乘、逻辑回归、probit回归、Huber稳健回归等。在这种情况下,对于所有
i
i
i,我们有
L
i
≤
M
∥
a
i
∥
2
+
λ
L_i \leq M\|a_i\|^2 + \lambda
Li≤M∥ai∥2+λ 并且
μ
≥
λ
\mu \geq \lambda
μ≥λ。
在这些假设下,梯度下降(GD)方法的收敛速率由条件数 κ : = L / μ \kappa := L/\mu κ:=L/μ 决定。条件数总是大于或等于1,当它显著大于1时,函数的等值线变得非常椭圆形,导致GD方法的迭代产生振荡。相反,当 κ \kappa κ 接近1时,GD方法收敛得更快。
在假设1和假设2下,VR方法以线性速率收敛。我们说一个随机方法的函数值 ({f(x_k)}) 以
0
<
ρ
≤
1
0 < \rho \leq 1
0<ρ≤1 的速率线性收敛(在期望下),如果存在一个常数
C
>
0
C > 0
C>0 使得:
E
[
f
(
x
k
)
]
−
f
(
x
∗
)
≤
(
1
−
ρ
)
k
C
=
O
(
exp
(
−
k
ρ
)
)
∀
k
(
16
)
E[f(x_k)] - f(x^*) \leq (1 - \rho)^k C = O(\exp(-k\rho)) \quad \forall k \quad (16)
E[f(xk)]−f(x∗)≤(1−ρ)kC=O(exp(−kρ))∀k(16)
这与每次迭代仅依赖于梯度无偏估计的经典SGD方法形成对比,后者在这些假设下只能获得次线性速率:
E
[
f
(
x
k
)
]
−
f
(
x
∗
)
≤
O
(
1
/
k
)
E[f(x_k)] - f(x^*) \leq O(1/k)
E[f(xk)]−f(x∗)≤O(1/k)
满足这个不等式的最小
k
k
k 称为算法的迭代复杂度。以下是GD、SGD和VR方法的基本变体的迭代复杂度和一次迭代的成本:
算法 | 迭代次数 | 一次迭代的成本 |
---|---|---|
GD | O ( κ log ( 1 / ϵ ) ) O(\kappa \log(1/\epsilon)) O(κlog(1/ϵ)) | O ( n ) O(n) O(n) |
SGD | O ( κ max max ( 1 / ϵ ) ) O(\kappa_{\text{max}} \max(1/\epsilon)) O(κmaxmax(1/ϵ)) | O ( 1 ) O(1) O(1) |
VR | O ( ( κ max + n ) log ( 1 / ϵ ) ) O((\kappa_{\text{max}} + n) \log(1/\epsilon)) O((κmax+n)log(1/ϵ)) | O ( 1 ) O(1) O(1) |
算法的总运行时间由迭代复杂度和迭代运行时间的乘积决定。这里使用了 κ max : = max i L i / μ \kappa_{\text{max}} := \max_i L_i/\mu κmax:=maxiLi/μ。注意 κ max ≥ κ \kappa_{\text{max}} \geq \kappa κmax≥κ;因此,GD的迭代复杂度小于VR方法。
然而,由于GD的每次迭代成本是VR方法的 n n n 倍,VR方法在总运行时间方面更优越。
经典SGD方法的优势在于它们的运行时间和收敛速率不依赖于 n n n,但它对容差 ϵ \epsilon ϵ 的依赖性要差得多,这解释了当容差很小时SGD的性能较差。
在附录B中,我们提供了一个简单的证明,表明SGD²方法具有与VR方法相同的迭代复杂度。
2.基础方差减少方法
方差减少(VR)方法的发展经历了几个阶段,最初的一批方法使得收敛速率得到了显著提升。这一系列方法的开端是SAG算法。随后,随机对偶坐标上升(SDCA)算法、MISO算法、随机方差减少梯度(SVRG/S2GD)算法,以及SAGA(意为“改进的”SAG)算法相继问世。
在本章中,我们将详细介绍这些开创性的VR方法。而在第四章,我们会探讨一些更新的方法,它们在特定的应用场景中相比这些基础方法展现出了更优越的特性。
2.1.随机平均梯度方法(SAG)
我们对第一种方差减少(VR)方法的探索,始于对完整梯度结构的模仿。既然完整梯度 ∇ f ( x ) \nabla f(x) ∇f(x) 是所有 ∇ f i ( x ) \nabla f_i(x) ∇fi(x) 梯度的简单平均,那么我们对完整梯度的估计 g k g_k gk 也应该是这些梯度估计的平均。这种思想催生了我们的第一种VR方法:随机平均梯度(SAG)方法。
SAG方法[37],[65]是早期增量聚合梯度(IAG)方法[4]的随机化版本。SAG的核心思想是对每个数据点
i
i
i 维护一个估计值
v
i
k
≈
∇
f
i
(
x
k
)
v_{ik} \approx \nabla f_i(x_k)
vik≈∇fi(xk)。然后,用这些
v
i
k
v_{ik}
vik 值的平均来作为对完整梯度的估计,即:
g
ˉ
k
=
1
n
∑
j
=
1
n
v
j
k
≈
1
n
∑
j
=
1
n
∇
f
j
(
x
k
)
=
∇
f
(
x
k
)
(
18
)
\bar{g}_k = \frac{1}{n} \sum_{j=1}^{n} v_{jk} \approx \frac{1}{n} \sum_{j=1}^{n} \nabla f_j(x_k) = \nabla f(x_k) \quad (18)
gˉk=n1j=1∑nvjk≈n1j=1∑n∇fj(xk)=∇f(xk)(18)
在SAG的每次迭代中,从集合
{
1
,
…
,
n
}
\{1, \ldots, n\}
{1,…,n} 中抽取一个索引
i
k
i_k
ik,然后根据以下规则更新
v
j
k
v_{jk}
vjk:
v
j
k
k
+
1
=
{
∇
f
i
k
(
x
k
)
,
if
j
=
i
k
v
j
k
k
,
if
j
≠
i
k
(
19
)
v_{jk}^{k+1} = \begin{cases} \nabla f_{i_k}(x_k), & \text{if } j = i_k \\ v_{jk}^k, & \text{if } j \neq i_k \end{cases} \quad (19)
vjkk+1={∇fik(xk),vjkk,if j=ikif j=ik(19)
其中,每个
v
0
i
v_{0i}
v0i 可以初始化为零或
∇
f
i
(
x
0
)
\nabla f_i(x_0)
∇fi(x0) 的近似值。随着解
x
∗
x^*
x∗ 的逼近,每个
v
i
k
v_{ik}
vik 会逐渐收敛到
∇
f
i
(
x
∗
)
\nabla f_i(x^*)
∇fi(x∗),从而满足VR属性(12)。
为了高效实现SAG,我们需要注意在计算
g
ˉ
k
\bar{g}_k
gˉk 时避免每次都从头开始求和
n
n
n 个向量,因为这在
n
n
n 很大时成本很高。幸运的是,由于每次迭代只有一个
v
i
k
v_{ik}
vik 项会改变,我们可以不必每次都重新计算整个和。具体来说,假设在迭代
k
k
k 中抽取了索引
i
k
i_k
ik,则有:
g
ˉ
k
=
1
n
∑
j
=
1
j
≠
i
k
n
v
j
k
+
1
n
v
i
k
k
=
g
ˉ
k
−
1
−
1
n
v
i
k
k
−
1
+
1
n
v
i
k
k
(
20
)
\bar{g}_k = \frac{1}{n} \sum_{\substack{j=1 \\ j \neq i_k}}^{n} v_{jk} + \frac{1}{n} v_{i_k}^k = \bar{g}_{k-1} - \frac{1}{n} v_{i_k}^{k-1} + \frac{1}{n} v_{i_k}^k \quad (20)
gˉk=n1j=1j=ik∑nvjk+n1vikk=gˉk−1−n1vikk−1+n1vikk(20)
由于除了 v i k v_{i_k} vik 之外的所有 v j k v_{jk} vjk 值都保持不变,我们只需存储每个 j j j 对应的一个向量 v j v_j vj。算法1展示了SAG方法的具体实现。
SAG是首个实现线性收敛的随机方法,其迭代复杂度为 O ( ( κ max + n ) log ( 1 / ϵ ) ) O((\kappa_{\text{max}} + n) \log(1/\epsilon)) O((κmax+n)log(1/ϵ)),使用步长 γ = O ( 1 / L max ) \gamma = O(1/L_{\text{max}}) γ=O(1/Lmax)。这种线性收敛性可以在图1中观察到。值得注意的是,由于 L max L_{\text{max}} Lmax-平滑函数对于任何 L ′ ≥ L max L' \geq L_{\text{max}} L′≥Lmax 也是 L ′ L' L′-平滑的,SAG方法对于足够小的步长都能获得线性收敛速率,这与经典SGD方法形成鲜明对比,后者只有在难以在实践中调整的递减步长序列下才能获得次线性速率。
在当时,SAG的线性收敛是一个显著的进展,因为它在每次迭代中只计算了一个随机梯度(处理单个数据点)。然而,Schmidt等人[65]提供的收敛证明非常复杂,并且依赖于计算机验证的步骤。SAG难以分析的一个关键原因是 g k g_k gk 是梯度的一个有偏估计。
接下来,我们将介绍SAGA方法,这是SAG的一个变体,它利用协变量的概念来创建一个无偏的SAG方法变体,该变体具有类似的性能但更易于分析。
算法 1:SAG 方法
- 参数:步长 γ > 0 \gamma > 0 γ>0
- 初始化: x 0 x_0 x0, v i = 0 ∈ R d v_i = 0 \in \mathbb{R}^d vi=0∈Rd 对于 i = 1 , … , n i = 1, \ldots, n i=1,…,n
- 对
k
=
1
,
…
,
T
−
1
k = 1, \ldots, T - 1
k=1,…,T−1 执行:
a. 随机抽取 i k ∈ { 1 , … , n } i_k \in \{1, \ldots, n\} ik∈{1,…,n}
b. 计算 g ˉ k = g ˉ k − 1 − 1 n v i k k − 1 \bar{g}_k = \bar{g}_{k-1} - \frac{1}{n} v_{i_k}^{k-1} gˉk=gˉk−1−n1vikk−1
c. 更新 v i k k = ∇ f i k ( x k ) v_{i_k}^k = \nabla f_{i_k}(x_k) vikk=∇fik(xk)
d. 更新梯度估计 g ˉ k = g ˉ k + 1 n v i k k \bar{g}_k = \bar{g}_k + \frac{1}{n} v_{i_k}^k gˉk=gˉk+n1vikk
e. 更新 x k + 1 = x k − γ g ˉ k x_{k+1} = x_k - \gamma \bar{g}_k xk+1=xk−γgˉk - 输出: x T x_T xT
2.2.SAGA方法
一种减少基本无偏梯度估计
∇
f
i
k
(
x
k
)
\nabla f_{i_k}(x_k)
∇fik(xk) 方差的方法是通过使用所谓的协变量(或称控制变量)。对于
i
=
1
,
…
,
n
i = 1, \ldots, n
i=1,…,n,设
v
i
∈
R
d
v_i \in \mathbb{R}^d
vi∈Rd 是一个向量。利用这些向量,我们可以将完整梯度
∇
f
(
x
)
\nabla f(x)
∇f(x) 重写为:
∇
f
(
x
)
=
1
n
∑
i
=
1
n
(
∇
f
i
(
x
)
−
v
i
+
v
i
)
=
1
n
∑
i
=
1
n
∇
f
i
(
x
)
−
v
i
+
1
n
∑
j
=
1
n
v
j
\nabla f(x) = \frac{1}{n} \sum_{i=1}^{n}(\nabla f_i(x) - v_i + v_i) = \frac{1}{n} \sum_{i=1}^{n} \nabla f_i(x) - v_i + \frac{1}{n} \sum_{j=1}^{n} v_j
∇f(x)=n1i=1∑n(∇fi(x)−vi+vi)=n1i=1∑n∇fi(x)−vi+n1j=1∑nvj
:
=
1
n
∑
i
=
1
n
∇
f
i
(
x
,
v
)
(
21
)
:= \frac{1}{n} \sum_{i=1}^{n} \nabla f_i(x, v) \quad (21)
:=n1i=1∑n∇fi(x,v)(21)
其中定义
∇
f
i
(
x
,
v
)
:
=
∇
f
i
(
x
)
−
v
i
+
1
n
∑
j
=
1
n
v
j
\nabla f_i(x, v) := \nabla f_i(x) - v_i + \frac{1}{n} \sum_{j=1}^{n} v_j
∇fi(x,v):=∇fi(x)−vi+n1∑j=1nvj。现在,我们可以通过随机抽样一个
∇
f
i
(
x
,
v
)
\nabla f_i(x, v)
∇fi(x,v) 来构建完整梯度
∇
f
(
x
)
\nabla f(x)
∇f(x) 的无偏估计,对于
i
∈
{
1
,
…
,
n
}
i \in \{1, \ldots, n\}
i∈{1,…,n},可以应用SGD方法,并使用梯度估计:
g
k
=
∇
f
i
k
(
x
k
,
v
)
=
∇
f
i
k
(
x
k
)
−
v
i
k
+
1
n
∑
j
=
1
n
v
j
(
22
)
g_k = \nabla f_{i_k}(x_k, v) = \nabla f_{i_k}(x_k) - v_{i_k} + \frac{1}{n} \sum_{j=1}^{n} v_j \quad (22)
gk=∇fik(xk,v)=∇fik(xk)−vik+n1j=1∑nvj(22)
为了观察
v
i
v_i
vi 的选择对方差
g
k
g_k
gk 的影响,我们可以将
g
k
=
∇
f
i
k
(
x
k
,
v
)
g_k = \nabla f_{i_k}(x_k, v)
gk=∇fik(xk,v) 代入,并利用
E
i
∼
1
n
[
v
i
]
=
1
n
∑
j
=
1
n
v
j
E_i \sim \frac{1}{n}[v_i] = \frac{1}{n} \sum_{j=1}^{n} v_j
Ei∼n1[vi]=n1∑j=1nvj 来计算期望,得到:
E
[
∥
∇
f
i
(
x
k
)
−
v
i
+
E
i
∼
1
n
[
v
i
−
∇
f
i
(
x
k
)
]
∥
2
]
≤
E
[
∥
∇
f
i
(
x
k
)
−
v
i
∥
2
]
(
23
)
E \left[ \|\nabla f_i(x_k) - v_i + E_i \sim \frac{1}{n}[v_i - \nabla f_i(x_k)]\|^2 \right] \leq E \left[ \|\nabla f_i(x_k) - v_i\|^2 \right] \quad (23)
E[∥∇fi(xk)−vi+Ei∼n1[vi−∇fi(xk)]∥2]≤E[∥∇fi(xk)−vi∥2](23)
这里使用了引理2,其中
X
=
∇
f
i
(
x
k
)
−
v
i
X = \nabla f_i(x_k) - v_i
X=∇fi(xk)−vi。这个界限 (23) 表明,如果
v
i
v_i
vi 随着
k
k
k 的增加接近
∇
f
i
(
x
k
)
\nabla f_i(x_k)
∇fi(xk),我们就能获得VR属性 (12)。这就是为什么我们称
v
i
v_i
vi 为协变量,并且我们可以选择它们来减少方差。
例如,SGD² 方法 (13) 也实现了这种方法,其中
v
i
=
∇
f
i
(
x
∗
)
v_i = \nabla f_i(x^*)
vi=∇fi(x∗)。然而,这在实践中不常用,因为我们通常不知道
∇
f
i
(
x
∗
)
\nabla f_i(x^*)
∇fi(x∗)。一个更实用的选择是
v
i
v_i
vi 作为我们知道的点
x
ˉ
i
∈
R
d
\bar{x}_i \in \mathbb{R}^d
xˉi∈Rd 附近的梯度
∇
f
i
(
x
ˉ
i
)
\nabla f_i(\bar{x}_i)
∇fi(xˉi)。SAGA 对每个函数
f
i
f_i
fi 使用一个参考点
x
ˉ
i
∈
R
d
\bar{x}_i \in \mathbb{R}^d
xˉi∈Rd,并使用协变量
v
i
=
∇
f
i
(
x
ˉ
i
)
v_i = \nabla f_i(\bar{x}_i)
vi=∇fi(xˉi),其中每个
x
ˉ
i
\bar{x}_i
xˉi 将是我们最后一次评估
f
i
f_i
fi 的点。使用这些协变量,我们可以构建梯度估计,按照 (22),给出:
g
k
=
∇
f
i
k
(
x
k
)
−
∇
f
i
k
(
x
ˉ
i
k
)
+
1
n
∑
j
=
1
n
∇
f
j
(
x
ˉ
j
)
(
24
)
g_k = \nabla f_{i_k}(x_k) - \nabla f_{i_k}(\bar{x}_{i_k}) + \frac{1}{n} \sum_{j=1}^{n} \nabla f_j(\bar{x}_j) \quad (24)
gk=∇fik(xk)−∇fik(xˉik)+n1j=1∑n∇fj(xˉj)(24)
为了实现SAGA,我们可以存储梯度 ∇ f i ( x ˉ i ) \nabla f_i(\bar{x}_i) ∇fi(xˉi) 而不是 n n n 个参考点 x ˉ i \bar{x}_i xˉi。也就是说,设 v j = ∇ f j ( x ˉ j ) v_j = \nabla f_j(\bar{x}_j) vj=∇fj(xˉj) 对于 j ∈ { 1 , … , n } j \in \{1, \ldots, n\} j∈{1,…,n},在每次迭代中,我们像SAG一样更新一个随机梯度的 v j v_j vj。
算法 2 SAGA
- 参数:步长 γ > 0 \gamma > 0 γ>0
- 初始化: x 0 x_0 x0, v i = 0 ∈ R d v_i = 0 \in \mathbb{R}^d vi=0∈Rd 对于 i = 1 , … , n i = 1, \ldots, n i=1,…,n
- 进行
k
=
1
,
…
,
T
−
1
k = 1, \ldots, T - 1
k=1,…,T−1 次迭代:
a. 随机抽取 i k ∈ { 1 , … , n } i_k \in \{1, \ldots, n\} ik∈{1,…,n}
b. 保存旧值 v old = v i k v_{\text{old}} = v_{i_k} vold=vik
c. 更新 v i k = ∇ f i k ( x k ) v_{i_k} = \nabla f_{i_k}(x_k) vik=∇fik(xk)
d. 更新 x k + 1 = x k − γ ( v i k − v old + g ˉ k ) x_{k+1} = x_k - \gamma (v_{i_k} - v_{\text{old}} + \bar{g}_k) xk+1=xk−γ(vik−vold+gˉk)
e. 更新梯度估计 g ˉ k = g ˉ k − 1 + 1 n ( v i k − v old ) \bar{g}_k = \bar{g}_{k-1} + \frac{1}{n} (v_{i_k} - v_{\text{old}}) gˉk=gˉk−1+n1(vik−vold) - 输出: x T x_T xT
SAGA方法具有与SAG相同的迭代复杂度 O ( ( κ max + n ) log ( 1 / ϵ ) ) O((\kappa_{\text{max}} + n) \log(1/\epsilon)) O((κmax+n)log(1/ϵ)),使用步长 γ = O ( 1 / L max ) \gamma = O(1/L_{\text{max}}) γ=O(1/Lmax),但证明要简单得多。然而,与SAG一样,SAGA方法需要存储 n n n 个辅助向量 v i ∈ R d v_i \in \mathbb{R}^d vi∈Rd 对于 i = 1 , … , n i = 1, \ldots, n i=1,…,n,这意味着需要 O ( n d ) O(nd) O(nd) 的存储空间。当 d d d 和 n n n 都很大时,这可能是不可行的。在下一节中,我们将详细说明如何为常见模型(如正则化线性模型)减少这种内存需求。
当能够将 n n n 个辅助向量存储在内存中时,SAG和SAGA的表现往往相似。如果这个内存需求太高,我们将在下一节中回顾的SVRG方法是一个不错的选择。SVRG方法实现了相同的收敛速率,并且在实践中通常几乎一样快,但只要求 O ( d ) O(d) O(d) 的内存,对于一般问题。
2.3.SVRG方法
在SAGA方法出现之前,一些早期的工作首次引入了协变量,以解决SAG方法所要求的高内存问题。这些研究构建了基于一个固定参考点 x ˉ ∈ R d \bar{x} \in \mathbb{R}^d xˉ∈Rd 的协变量,我们已经在该点计算了完整的梯度 ∇ f ( x ˉ ) \nabla f(\bar{x}) ∇f(xˉ)。通过存储参考点 x ˉ \bar{x} xˉ 和对应的完整梯度 ∇ f ( x ˉ ) \nabla f(\bar{x}) ∇f(xˉ),我们可以在不存储每个 ∇ f j ( x ˉ ) \nabla f_j(\bar{x}) ∇fj(xˉ) 的情况下,使用 x ˉ j = x ˉ \bar{x}_j = \bar{x} xˉj=xˉ 对所有 j j j 来实现更新 (24)。具体来说,我们不是存储这些向量,而是在每次迭代中利用存储的参考点 x ˉ \bar{x} xˉ 来计算 ∇ f i k ( x ˉ ) \nabla f_{i_k}(\bar{x}) ∇fik(xˉ)。这个方法最初由不同的作者以不同的名字提出,但后来统一被称为SVRG方法,遵循[28]和[84]的命名。
我们在算法3中对SVRG方法进行了形式化。
利用 (23),我们可以得出梯度估计
g
k
g_k
gk 的方差有界:
E
[
∥
g
k
−
∇
f
(
x
k
)
∥
2
]
≤
E
[
∥
∇
f
i
(
x
k
)
−
∇
f
i
(
x
ˉ
)
∥
2
]
≤
L
max
2
∥
x
k
−
x
ˉ
∥
2
E\left[ \| g_k - \nabla f(x_k) \|^2 \right] \leq E\left[ \| \nabla f_i(x_k) - \nabla f_i(\bar{x}) \|^2 \right] \leq L_{\text{max}}^2 \| x_k - \bar{x} \|^2
E[∥gk−∇f(xk)∥2]≤E[∥∇fi(xk)−∇fi(xˉ)∥2]≤Lmax2∥xk−xˉ∥2
其中第二个不等式使用了每个
f
i
f_i
fi 的
L
i
L_i
Li-平滑性。
值得注意的是,参考点 x ˉ \bar{x} xˉ 越接近当前点 x k x_k xk,梯度估计的方差就越小。
为了让SVRG方法有效,我们需要在频繁更新参考点 x ˉ \bar{x} xˉ(从而需要计算完整梯度)的成本与降低方差的好处之间做出权衡。为此,我们每 t t t 次迭代更新一次参考点,使其接近 x k x_k xk(见算法II-C的第11行)。也就是说,SVRG方法包含两个循环:一个外循环 s s s,其中计算参考梯度 ∇ f ( x ˉ s − 1 ) \nabla f(\bar{x}_{s-1}) ∇f(xˉs−1)(第4行),以及一个内循环,其中固定参考点,并根据随机梯度步骤(22)更新内部迭代 x k x_k xk(第10行)。
与SAG和SAGA不同,SVRG只需要 O ( d ) O(d) O(d) 的内存。SVRG的缺点包括:1) 我们有一个额外的参数 t t t,即内循环的长度,需要调整;2) 每次迭代需要计算两个梯度,并且每次更改参考点时都需要计算完整梯度。
Johnson和Zhang[28]展示了SVRG具有迭代复杂度 O ( ( κ max + n ) log ( 1 / ϵ ) ) O((\kappa_{\text{max}} + n) \log(1/\epsilon)) O((κmax+n)log(1/ϵ)),与SAG和SAGA相似。这是在假设内循环次数 t t t 从集合 { 1 , … , m } \{1, \ldots, m\} {1,…,m} 中均匀抽样的情况下得出的,其中 L max L_{\text{max}} Lmax, μ \mu μ,步长 γ \gamma γ 和 t t t 之间必须满足一定的依赖关系。在实践中,通过使用 γ = O ( 1 / L max ) \gamma = O(1/L_{\text{max}}) γ=O(1/Lmax) 和内循环长度 t = n t = n t=n,SVRG往往表现良好,这正是我们在图1中使用的设置。
现在,有许多原始SVRG方法的变体。例如,有些变体使用
t
t
t 的替代分布[32],有些变体允许形式为
O
(
1
/
L
max
)
O(1/L_{\text{max}})
O(1/Lmax) 的步长[27],[33],[35]。还有一些变体使用
∇
f
(
x
ˉ
)
\nabla f(\bar{x})
∇f(xˉ) 的小批量近似来减少这些完整梯度评估的成本,并增加小批量大小以保持VR属性。还有一些变体在内循环中根据[54]重复更新
g
k
g_k
gk:
[ g_k = \nabla f_{i_k}(x_k) - \nabla f_{i_k}(x_{k-1}) + g_{k-1} \quad (25) ]
这提供了更局部的近似。使用这种连续更新变体 (25) 在最小化非凸函数时显示出独特的优势,正如我们在第四节简要讨论的。最后,注意SVRG可以利用
∇
f
(
x
ˉ
s
)
\nabla f(\bar{x}_s)
∇f(xˉs) 的值来帮助决定何时终止算法。
算法 3 SVRG方法
- 参数:步长 γ > 0 \gamma > 0 γ>0
- 初始化参考点 x ˉ 0 = x 0 ∈ R d \bar{x}_0 = x_0 \in \mathbb{R}^d xˉ0=x0∈Rd
- 进行外循环
s
=
1
,
2
,
…
s = 1, 2, \ldots
s=1,2,…:
a. 计算并存储 ∇ f ( x ˉ s − 1 ) \nabla f(\bar{x}_{s-1}) ∇f(xˉs−1)
b. 设 x 0 = x ˉ s − 1 x_0 = \bar{x}_{s-1} x0=xˉs−1
c. 选择内循环迭代次数 t t t
d. 进行内循环 k = 0 , 1 , … , t − 1 k = 0, 1, \ldots, t - 1 k=0,1,…,t−1:
i. 随机抽取 i k ∈ { 1 , … , n } i_k \in \{1, \ldots, n\} ik∈{1,…,n}
ii. 计算 g k = ∇ f i k ( x k ) − ∇ f i k ( x ˉ s − 1 ) + ∇ f ( x ˉ s − 1 ) g_k = \nabla f_{i_k}(x_k) - \nabla f_{i_k}(\bar{x}_{s-1}) + \nabla f(\bar{x}_{s-1}) gk=∇fik(xk)−∇fik(xˉs−1)+∇f(xˉs−1)
iii. 更新 x k + 1 = x k − γ g k x_{k+1} = x_k - \gamma g_k xk+1=xk−γgk
e. 更新参考点 x ˉ s = x t \bar{x}_s = x_t xˉs=xt
2.4. SDCA及其变体
SAG和SVRG方法的一个不足之处在于,它们的步长依赖于可能在某些问题中未知的
L
max
L_{\text{max}}
Lmax。在SVRG之前,SDCA方法[70]作为最早的VR方法之一,将坐标下降方法的研究扩展到了有限和问题。SDCA及其变体背后的理念是,梯度的坐标提供了一种自然的方差减少梯度估计。具体来说,设
j
∈
{
1
,
…
,
d
}
j \in \{1, \ldots, d\}
j∈{1,…,d},并且定义
∇
j
f
(
x
)
:
=
(
∂
f
(
x
)
∂
x
j
)
e
j
\nabla_j f(x) := \left( \frac{\partial f(x)}{\partial x_j} \right) e_j
∇jf(x):=(∂xj∂f(x))ej 为 (f(x)) 的第
j
j
j 个坐标方向的导数,其中
e
j
∈
R
d
e_j \in \mathbb{R}^d
ej∈Rd 是第
j
j
j 个单位向量。坐标导数的一个关键特性是
∇
j
f
(
x
∗
)
=
0
\nabla_j f(x^*) = 0
∇jf(x∗)=0,这是因为我们知道
∇
f
(
x
∗
)
=
0
\nabla f(x^*) = 0
∇f(x∗)=0。这与每个数据点的导数
∇
f
j
\nabla f_j
∇fj 不同,后者在
x
∗
x^*
x∗ 处可能不为零。因此,我们有:
∥
∇
f
(
x
)
−
∇
j
f
(
x
)
∥
2
→
0
当
x
→
x
∗
(
26
)
\| \nabla f(x) - \nabla_j f(x) \|^2 \rightarrow 0 \quad \text{当} \quad x \rightarrow x^* \quad (26)
∥∇f(x)−∇jf(x)∥2→0当x→x∗(26)
这意味着坐标导数满足了方差减少属性(12)。此外,我们还可以使用
∇
j
f
(
x
)
\nabla_j f(x)
∇jf(x) 来构建
∇
f
(
x
)
\nabla f(x)
∇f(x) 的无偏估计。例如,设
j
j
j 是从集合
{
1
,
…
,
d
}
\{1, \ldots, d\}
{1,…,d} 中均匀随机选取的索引。因此,对于任何
i
∈
{
1
,
…
,
d
}
i \in \{1, \ldots, d\}
i∈{1,…,d},我们有
P
[
j
=
i
]
=
1
d
P[j = i] = \frac{1}{d}
P[j=i]=d1。因此,
d
×
∇
j
f
(
x
)
d \times \nabla_j f(x)
d×∇jf(x) 是
∇
f
(
x
)
\nabla f(x)
∇f(x) 的无偏估计,因为:
E
[
d
∇
j
f
(
x
)
]
=
d
∑
i
=
1
d
P
[
j
=
i
]
∂
f
(
x
)
∂
x
i
e
i
=
∑
i
=
1
d
∂
f
(
x
)
∂
x
i
e
i
=
∇
f
(
x
)
E\left[ d \nabla_j f(x) \right] = d \sum_{i=1}^{d} P[j = i] \frac{\partial f(x)}{\partial x_i} e_i = \sum_{i=1}^{d} \frac{\partial f(x)}{\partial x_i} e_i = \nabla f(x)
E[d∇jf(x)]=di=1∑dP[j=i]∂xi∂f(x)ei=i=1∑d∂xi∂f(x)ei=∇f(x)
因此, ∇ j f ( x ) \nabla_j f(x) ∇jf(x) 具有我们期望的VR估计完整梯度的所有理想属性,而且不需要使用协变量。使用这种坐标梯度的一个缺点是,对于我们的和问题(2),它的计算成本很高。这是因为计算 ∇ j f ( x ) \nabla_j f(x) ∇jf(x) 需要遍历整个数据集,因为 ∇ j f ( x ) = 1 n ∑ i = 1 n ∇ j f i ( x ) \nabla_j f(x) = \frac{1}{n} \sum_{i=1}^{n} \nabla_j f_i(x) ∇jf(x)=n1∑i=1n∇jfi(x)。因此,使用坐标导数似乎与我们的和问题的结构不兼容。然而,我们可以经常将原始问题(2)重写为所谓的对偶公式,其中坐标导数可以利用固有的结构。
例如,L2正则化线性模型(15)的对偶公式为:
v
∗
∈
arg
max
v
∈
R
n
1
n
∑
i
=
1
n
−
ℓ
i
∗
(
−
v
i
)
−
λ
2
∥
1
λ
∑
i
=
1
n
v
i
a
i
∥
2
(
27
)
v^* \in \arg\max_{v \in \mathbb{R}^n} \frac{1}{n} \sum_{i=1}^{n} -\ell_i^*(-v_i) - \frac{\lambda}{2} \left\| \frac{1}{\lambda} \sum_{i=1}^{n} v_i a_i \right\|^2 \quad (27)
v∗∈argv∈Rnmaxn1i=1∑n−ℓi∗(−vi)−2λ
λ1i=1∑nviai
2(27)
其中
ℓ
i
∗
(
v
)
\ell_i^*(v)
ℓi∗(v) 是
ℓ
i
\ell_i
ℓi 的凸共轭。我们可以使用映射
x
=
1
λ
∑
i
=
1
n
v
i
a
i
x = \frac{1}{\lambda} \sum_{i=1}^{n} v_i a_i
x=λ1∑i=1nviai 来恢复原始问题(15)中的
x
x
x 变量。将解
v
∗
v^*
v∗ 代入上述映射的右侧可以得到(15)的解
x
∗
x^*
x∗。
注意,这个对偶问题有 n n n 个实变量 v i ∈ R v_i \in \mathbb{R} vi∈R,每个训练样本对应一个。此外,每个对偶损失函数 ℓ i ∗ \ell_i^* ℓi∗ 仅是 v i v_i vi 的函数。也就是说,损失函数中的第一项在坐标上是可分离的。这种在坐标上的可分离性,加上第二项的简单形式,允许我们有效实现坐标上升方法。实际上,Shalev-Shwartz和Zhang展示了在这个问题上的坐标上升具有与SAG、SAGA和SVRG类似的迭代复杂度 O ( ( κ max + n ) log ( 1 / ϵ ) ) O((\kappa_{\text{max}} + n) \log(1/\epsilon)) O((κmax+n)log(1/ϵ))。
迭代成本和算法结构也非常相似:通过跟踪求和 ∑ i = 1 n v i a i \sum_{i=1}^{n} v_i a_i ∑i=1nviai 来处理(27)中的第二项,每个对偶坐标上升迭代只需要考虑一个训练样本,并且每次迭代的成本与 n n n 无关。此外,我们可以使用一维线搜索有效地计算步长,以最大限度地提高作为 v i v_i vi 函数的对偶目标。这意味着,即使没有 L max L_{\text{max}} Lmax 或相关量的了解,也可以实现VR方法的快速最坏情况运行时间。
3.方差缩小的实践问题
为了实现基本的方差减少(VR)方法并取得合理的性能,必须解决几个实施问题。在本节中,我们将讨论上述未涉及的若干问题。
3.1.SAG/SAGA/SVRG设置步长
在优化算法的领域,特别是随机平均梯度(SAG)、随机平均梯度算法(SAGA)和随机梯度(SVRG)等变分减少方法中,步长的设置是一个关键问题。虽然对于随机对偶坐标上升(SDCA)方法,我们可以使用对偶目标来确定步长,但是SAG、SAGA和SVRG这些原始变量方法的理论依据是步长应为 γ = O ( 1 L max ) \gamma = O\left(\frac{1}{L_{\text{max}}}\right) γ=O(Lmax1) 的形式。然而,在实际应用中,我们往往不知道 L max L_{\text{max}} Lmax 的确切值,而且使用其他步长可能会得到更好的性能。
全梯度下降(full-GD)方法中设置步长的一种经典策略是Armijo线搜索。给定当前点
x
k
x_k
xk 和搜索方向
g
k
g_k
gk,Armijo线搜索在
γ
k
\gamma_k
γk 的线上进行,该线定义为
γ
k
∈
{
γ
:
x
k
+
γ
g
k
}
\gamma_k \in \{\gamma : x_k + \gamma g_k\}
γk∈{γ:xk+γgk},并且要求函数有充分的减少,即:
f
(
x
k
+
γ
k
g
k
)
<
f
(
x
k
)
−
c
γ
k
∥
∇
f
(
x
k
)
∥
2
f(x_k + \gamma_k g_k) < f(x_k) - c \gamma_k \|\nabla f(x_k)\|^2
f(xk+γkgk)<f(xk)−cγk∥∇f(xk)∥2
然而,这种方法需要在多个候选步长
γ
k
\gamma_k
γk 上计算
f
(
x
k
+
γ
k
g
k
)
f(x_k + \gamma_k g_k)
f(xk+γkgk),这在评估
f
(
x
)
f(x)
f(x) 需要遍历整个数据集时成本过高。
为了解决这个问题,可以采用随机变体的方法,寻找满足以下条件的
γ
k
\gamma_k
γk:
f
i
k
(
x
k
+
γ
k
g
k
)
<
f
i
k
(
x
k
)
−
c
γ
k
∥
∇
f
i
k
(
x
k
)
∥
2
f_{ik}(x_k + \gamma_k g_k) < f_{ik}(x_k) - c \gamma_k \|\nabla f_{ik}(x_k)\|^2
fik(xk+γkgk)<fik(xk)−cγk∥∇fik(xk)∥2
这种方法在实践中通常效果良好,尤其是在
∥
∇
f
i
k
(
x
k
)
∥
\|\nabla f_{ik}(x_k)\|
∥∇fik(xk)∥ 不接近零的情况下,尽管目前还没有理论支持这种方法。
另外,Mairal 提出了一种在实践中设置步长的“Bottou技巧”。这种方法通过取数据集的一小部分(例如5%)进行二分搜索,以尝试找到在通过这个样本进行一次遍历时的最优步长。与Armijo线搜索类似,这种方法在实践中通常表现良好,但同样缺乏理论基础。
请注意,上述内容是对原文的重新表述,使用了Markdown格式来表示数学公式和变量。
然而,SDCA方法也有一些缺点。首先,它需要计算凸共轭 ℓ i ∗ \ell_i^* ℓi∗ 而不是简单的梯度。我们没有凸共轭的自动微分等价物,所以这可能会增加实现工作量。最近的工作已经提出了不需要共轭的“无对偶”SDCA方法,而是直接使用梯度。然而,在这些方法中,不再可能跟踪对偶目标以设置步长。其次,尽管SDCA只需要 O ( n + d ) O(n + d) O(n+d) 的内存来解决(15)问题,但对于这个问题类别,SAG/SAGA也只需要 O ( n + d ) O(n + d) O(n+d) 的内存(见第三节)。适用于更一般问题的SDCA变体具有SAG/SAGA的 O ( n d ) O(nd) O(nd) 内存,因为 v i v_i vi 成为具有 d d d 个元素的向量。SDCA的一个最后的微妙缺点是它隐含地假设强凸性常数 μ \mu μ 等于 λ \lambda λ。对于 μ \mu μ 大于 λ \lambda λ 的问题,原始的VR方法通常显著优于SDCA。
3.2. 终止条件的确定
在算法优化领域,我们通常依赖于迭代复杂度的理论结果来预测算法达到特定精度所需的最坏情况下的迭代次数。但是,这些理论界限往往依赖于一些我们无法预知的常数,而在实际应用中,算法往往能在更少的迭代次数内达到预期精度。因此,我们需要设立一些测试标准来决定何时应该结束算法的运行。
在传统的全梯度下降(full-GD)方法中,我们通常根据梯度的范数 ∥ ∇ f ( x k ) ∥ \| \nabla f(x_k) \| ∥∇f(xk)∥ 或者与此相关的其他量来决定何时停止迭代。对于SVRG方法,我们可以采用相同的准则,但使用 ∥ ∇ f ( x ˉ s ) ∥ \| \nabla f(\bar{x}_s) \| ∥∇f(xˉs)∥ 来作为判断依据。对于SAG/SAGA方法,尽管我们没有显式地计算完整的梯度,但量 $ g_{\bar{k}} $ 会逐渐逼近 ∇ f ( x k ) \nabla f(x_k) ∇f(xk),因此,使用 ∥ g k ˉ ∥ \| g_{\bar{k}} \| ∥gkˉ∥ 作为停止条件是一种合理的启发式方法。
在SDCA方法中,通过一些额外的记录工作,我们可以在不增加额外渐近成本的情况下跟踪对偶目标的梯度。此外,一种更为系统的方法是跟踪对偶间隙,虽然这会增加每迭代 O ( n ) O(n) O(n) 的成本,但它能够提供具有对偶间隙证明的终止条件。另外,基于强凸目标的最优性条件,MISO方法采用了一种基于二次下界[41]的原则性方法。
以下是使用Markdown格式表示的数学公式和变量:
- 梯度范数: ∥ ∇ f ( x k ) ∥ \| \nabla f(x_k) \| ∥∇f(xk)∥
- SVRG方法中的梯度范数: ∥ ∇ f ( x ˉ s ) ∥ \| \nabla f(\bar{x}_s) \| ∥∇f(xˉs)∥
- SAG/SAGA方法中逼近梯度的量:$ g_{\bar{k}} $
- 每迭代增加的成本: O ( n ) O(n) O(n)
- MISO方法
- 二次下界
请注意,上述内容是对原文的重新表述,使用了Markdown格式来表示数学公式和变量。
3.3. 减少内存需求
尽管随机梯度变分减少(SVRG)算法消除了早期变分减少方法的内存需求,但在实际应用中,SAG(随机平均梯度下降)和SAGA(带梯度累积的随机平均梯度下降)算法在很多问题上往往比SVRG算法需要更少的迭代次数。这引发了一个思考:是否存在某些问题,使得SAG/SAGA能够在 O ( n d ) O(nd) O(nd) 内存需求以下实现。本节将探讨线性模型类别,这类模型的内存需求可以显著降低。
考虑线性模型,其中每个函数
f
i
(
x
)
f_i(x)
fi(x) 可以表示为
ξ
i
(
a
i
⊤
x
)
\xi_i(\mathbf{a}_i^\top x)
ξi(ai⊤x)。对
x
x
x 求导得到梯度形式:
∇
f
i
(
x
)
=
ξ
′
(
a
i
⊤
x
)
a
i
\nabla f_i(x) = \xi'(\mathbf{a}_i^\top x) \mathbf{a}_i
∇fi(x)=ξ′(ai⊤x)ai
这里,
ξ
′
\xi'
ξ′ 表示
ξ
\xi
ξ 的导数。假设我们可以直接访问特征向量
a
i
\mathbf{a}_i
ai,那么为了实现SAG/SAGA方法,我们只需要存储标量
ξ
(
a
i
⊤
x
)
\xi(\mathbf{a}_i^\top x)
ξ(ai⊤x)。这样,内存需求就从
O
(
n
d
)
O(nd)
O(nd) 减少到了
O
(
n
)
O(n)
O(n)。SVRG算法也可以利用梯度的这种结构:通过存储这
n
n
n 个标量,我们可以将SVRG“内部”迭代中每次所需的梯度评估次数减少到1,对于这一类问题。
还有其他类型的问题,例如概率图模型,它们也提供了降低内存需求的可能性[66]。通过特定的数据结构和算法优化,可以进一步减少算法在运行时所需的内存资源。
以下是使用Markdown格式表示的数学公式和变量:
- 线性模型函数: f i ( x ) = ξ i ( a i ⊤ x ) f_i(x) = \xi_i(\mathbf{a}_i^\top x) fi(x)=ξi(ai⊤x)
- 梯度表达式: ∇ f i ( x ) = ξ ′ ( a i ⊤ x ) a i \nabla f_i(x) = \xi'(\mathbf{a}_i^\top x) \mathbf{a}_i ∇fi(x)=ξ′(ai⊤x)ai
- 特征向量: a i \mathbf{a}_i ai
- 内存需求从 O ( n d ) O(nd) O(nd) 减少到 O ( n ) O(n) O(n)。
3.4. 稀疏梯度的处理
在某些问题中,梯度 ∇ f i ( x ) \nabla f_i(x) ∇fi(x) 可能包含大量零值,例如具有稀疏特征的线性模型。在这种情况下,传统的随机梯度下降(SGD)算法可以高效实现,其计算复杂度与梯度中非零元素的数量成线性关系,这通常远小于问题维度 d d d。然而,在标准的变分减少(VR)方法中,这种优势并没有被利用。幸运的是,存在两种已知的改进方法。
第一种改进方法由Schmidt等人提出,它利用了更新过程的简单性,实现了一种“即时”计算的变体,使得每次迭代的成本与非零元素的数量成正比。以SAG为例(但这种方法适用于所有变体),具体做法是在每次迭代后不存储完整的向量 v i k v_{ik} vik,而是只计算对应于非零元素的 v i k j v_{ik_j} vikj,通过更新自上次该元素非零以来的每个变量 v i k j v_{ik_j} vikj。
第二种改进方法由Leblond等人为SAGA提出,它在更新公式 x k + 1 = x k − γ ( ∇ f i k ( x k ) − ∇ f i k ( x ˉ i k ) + g ˉ k ) x_{k+1} = x_k - \gamma(\nabla f_{ik}(x_k) - \nabla f_{ik}(\bar{x}_{ik}) + \bar{g}_k) xk+1=xk−γ(∇fik(xk)−∇fik(xˉik)+gˉk) 中引入了额外的随机性。这里, ∇ f i k ( x k ) \nabla f_{ik}(x_k) ∇fik(xk) 和 ∇ f i k ( x ˉ i k ) \nabla f_{ik}(\bar{x}_{ik}) ∇fik(xˉik) 是稀疏的,而 g ˉ k \bar{g}_k gˉk 是密集的。在这个方法中,密集项 ( g ˉ k ) j (\bar{g}_k)_j (gˉk)j 的每个分量被替换为 w j ( g ˉ k ) j w_j (\bar{g}_k)_j wj(gˉk)j,其中 w ∈ R d w \in \mathbb{R}^d w∈Rd 是一个随机稀疏向量,其支持集包含在 ∇ f i k ( x k ) \nabla f_{ik}(x_k) ∇fik(xk) 中,并且期望上是一个所有元素都为1的常数向量。这样,更新过程保持了无偏性(尽管现在是稀疏的),并且增加的方差不会影响算法的收敛速率。Leblond等人提供了更多的细节。
以下是使用Markdown格式表示的数学公式和变量:
- 梯度: ∇ f i ( x ) \nabla f_i(x) ∇fi(x)
- SGD更新: x k + 1 = x k − γ ( ∇ f i k ( x k ) − ∇ f i k ( x ˉ i k ) + g ˉ k ) x_{k+1} = x_k - \gamma(\nabla f_{ik}(x_k) - \nabla f_{ik}(\bar{x}_{ik}) + \bar{g}_k) xk+1=xk−γ(∇fik(xk)−∇fik(xˉik)+gˉk)
- 稀疏梯度: ∇ f i k ( x k ) \nabla f_{ik}(x_k) ∇fik(xk) 和 ∇ f i k ( x ˉ i k ) \nabla f_{ik}(\bar{x}_{ik}) ∇fik(xˉik)
- 密集梯度: g ˉ k \bar{g}_k gˉk
- 随机稀疏向量: w w w
- 期望常数向量:所有元素都为1的向量。