JavaScript is required

分布式处理与计算：数值优化方法

牛顿迭代法与收敛性理论基础

笔记#分布式处理与计算#数值优化#牛顿法#收敛性

一、前介

大 O 符号

定义:对于两个有相同定义域的函数 ff 与 gg,如果存在常数 C>0C > 0 以及 x0x_0,使得当 x≥x0x \ge x_0 时 ∣f(x)∣≤C ∣g(x)∣|f(x)| \le C\,|g(x)| ，记 f(x)=O(g(x))f(x) = O(g(x))

性质:f=O(1)f = O(1) 表示 ff 是有界的;f=O(g)f = O(g) 表示 ff 的增长不超过 gg 的量级。

证明思路:假设存在 C1,C2C_1, C_2 使得 f1=O(g1)f_1 = O(g_1)、f2=O(g2)f_2 = O(g_2),则对和式有 ∣f1+f2∣≤C1∣g1∣+C2∣g2∣≤(C1+C2)max⁡(∣g1∣,∣g2∣)|f_1+f_2| \le C_1|g_1| + C_2|g_2| \le (C_1+C_2)\max(|g_1|,|g_2|),令 C=C1+C2C = C_1+C_2 即得 f1+f2=O(max⁡(g1,g2))f_1+f_2 = O(\max(g_1,g_2)),得证。

一些简单关系:O(g1)+O(g2)=O(max⁡(g1,g2))O(g_1) + O(g_2) = O(\max(g_1, g_2)),O(g1)⋅O(g2)=O(g1g2)O(g_1)\cdot O(g_2) = O(g_1 g_2)。

小 o 符号

定义:如果对每一个常数 ϵ>0\epsilon > 0 都存在 x0x_0,使得当 x≥x0x \ge x_0 时，∣f(x)∣≤ϵ ∣g(x)∣|f(x)| \le \epsilon\,|g(x)| ， 记 f(x)=o(g(x))f(x) = o(g(x))

等价于 lim⁡x→∞f(x)g(x)=0\lim_{x\to\infty}\frac{f(x)}{g(x)} = 0,即 ff 的量级严格低于 gg

性质:如果 f=o(g)f = o(g),则也有 f=O(g)f = O(g)(小 o 比大 O 更强)。

一些简单关系:o(g1)+o(g2)=o(max⁡(g1,g2))o(g_1) + o(g_2) = o(\max(g_1,g_2)),o(g1)⋅o(g2)=o(g1g2)o(g_1)\cdot o(g_2) = o(g_1 g_2)。

Weierstrass 定理

Weierstrass 定理:取值于闭区间 [a,b][a,b] 上的连续函数 ff,其最大值、最小值一定存在,并且在区间中取得。

平衡点(驻点):使得梯度 ∇f(x)=0\nabla f(x) = 0 的点称为平衡点。但 ∇f(x)=0\nabla f(x)=0 并不能保证是全局最优点,它可能只是局部最优点或鞍点。

凸函数

凸函数:如果对任意 x1,x2x_1, x_2 及 λ∈[0,1]\lambda \in [0,1],有

f(λx1+(1−λ)x2)≤λf(x1)+(1−λ)f(x2)f\big(\lambda x_1 + (1-\lambda)x_2\big) \le \lambda f(x_1) + (1-\lambda)f(x_2)

那么 ff 被称为凸函数。

严格凸函数:如果对任意 x1≠x2x_1 \ne x_2 及 λ∈(0,1)\lambda \in (0,1),把上式的 ≤\le 换成严格的 <<,那么 ff 被称为严格凸函数。凸函数的重要意义在于:局部最优即全局最优。

中值定理

假设连续函数 ff 取值在闭区间 [a,b][a,b] 上,并且在 (a,b)(a,b) 上可导,那么一定可以找到一点 ξ∈(a,b)\xi \in (a,b) 使得

f′(ξ)=f(b)−f(a)b−af'(\xi) = \frac{f(b) - f(a)}{b - a}

泰勒展开

假设 ff 在点 x0x_0 附近的区间中 nn 阶可导,那么对任意 xx,都可以把 ff 在 x0x_0 处泰勒展开:

f(x)=f(x0)+f′(x0)(x−x0)+f′′(x0)2!(x−x0)2+⋯+f(n)(x0)n!(x−x0)n+o((x−x0)n)f(x) = f(x_0) + f'(x_0)(x-x_0) + \frac{f''(x_0)}{2!}(x-x_0)^2 + \cdots + \frac{f^{(n)}(x_0)}{n!}(x-x_0)^n + o\big((x-x_0)^n\big)

算法的收敛速率

如果存在常数 qq 使得

lim⁡k→∞∥xk+1−x∗∥∥xk−x∗∥ q=C(0<C<∞)\lim_{k\to\infty}\frac{\|x_{k+1} - x^*\|}{\|x_k - x^*\|^{\,q}} = C \quad (0 < C < \infty)

那么称 xkx_k 收敛到 x∗x^* 的速率为 qq 阶。当 q=1q=1 时称线性收敛速率,q=2q=2 时称二次收敛速率(牛顿法在良好条件下即二次收敛)。

梯度与 Hessian 矩阵

梯度向量(一阶偏导组成的列向量):

∇f(x)=(∂f∂x1,…,∂f∂xp)⊤\nabla f(x) = \left(\frac{\partial f}{\partial x_1}, \dots, \frac{\partial f}{\partial x_p}\right)^{\top}

Hessian 矩阵(二阶偏导组成的对称方阵):

∇2f(x)=[∂2f∂xi∂xj]p×p\nabla^2 f(x) = \left[\frac{\partial^2 f}{\partial x_i \partial x_j}\right]_{p\times p}

二、牛顿迭代法(Newton's method)

对于优化问题 min⁡xf(x)\min_x f(x),这里要求 ff 为二阶可导的凸函数。在给定当前估计 xkx_k 时,把 ff 在 xkx_k 处做二阶泰勒展开并对其极小化,可得迭代式:

xk+1=xk−[∇2f(xk)]−1∇f(xk)x_{k+1} = x_k - \big[\nabla^2 f(x_k)\big]^{-1}\nabla f(x_k)

证明:对二阶泰勒近似 f(x)≈f(xk)+∇f(xk)⊤(x−xk)+12(x−xk)⊤∇2f(xk)(x−xk)f(x)\approx f(x_k)+\nabla f(x_k)^{\top}(x-x_k)+\frac{1}{2}(x-x_k)^{\top}\nabla^2 f(x_k)(x-x_k) 求导并令导数为 0,得 ∇f(xk)+∇2f(xk)(x−xk)=0\nabla f(x_k) + \nabla^2 f(x_k)(x - x_k) = 0,解出 x=xk−[∇2f(xk)]−1∇f(xk)x = x_k - [\nabla^2 f(x_k)]^{-1}\nabla f(x_k)。

NR 算法的一些问题

  • Hessian 矩阵求逆的计算量在某些情况下过大,并且 Hessian 可能不可逆;
  • 在 NR 算法中,函数值 f(xk)f(x_k) 并不保证随每次迭代递减(步长可能太大冲过头);
  • NR 算法对初值的依赖程度较高,初值不好可能不收敛。

提出背景:为了解决 NR 算法在某些情况下迭代步长太大、导致不收敛的问题,用回溯线搜索筛选出既足够长、又能保证 ff 递减的步长。

算法步骤:

  1. 给定参数 α∈(0,0.5)\alpha \in (0, 0.5)、β∈(0,1)\beta \in (0,1),令步长 t=1t = 1;
  2. 判断下面的充分下降条件是否成立:

f(x+t Δx)≤f(x)+α t ∇f(x)⊤Δxf(x + t\,\Delta x) \le f(x) + \alpha\,t\,\nabla f(x)^{\top}\Delta x

若为真则进行第 3 步(接受这个 tt),否则进行第 4 步;

  1. 接受步长,更新 x←x+tΔxx \leftarrow x + t\Delta x;
  2. 缩小步长 t←βtt \leftarrow \beta t,并重新回到第 2 步

解释:充分下降条件成立时,说明 xx 与 x+tΔxx+t\Delta x 之间连线的斜率,比 xx 点处斜率的 α\alpha 倍还要平缓——若不满足就说明步子迈太大,需要缩小步长。

收敛标准

  • 绝对收敛标准:∥xk+1−xk∥<ϵ\|x_{k+1} - x_k\| < \epsilon;
  • 相对收敛标准:∥xk+1−xk∥∥xk∥<ϵ\dfrac{\|x_{k+1} - x_k\|}{\|x_k\|} < \epsilon;
  • 调整的相对收敛标准:∥xk+1−xk∥∥xk∥+1<ϵ\dfrac{\|x_{k+1} - x_k\|}{\|x_k\| + 1} < \epsilon(避免 ∥xk∥\|x_k\| 太小时分母失稳)。

例子:逻辑回归

假设有一组样本 (xi,yi)(x_i, y_i),考虑逻辑回归 pi=11+e−xi⊤βp_i = \frac{1}{1+e^{-x_i^{\top}\beta}}。似然函数与对数似然分别为

L(β)=∏ipiyi(1−pi)1−yi,ℓ(β)=∑i[yilog⁡pi+(1−yi)log⁡(1−pi)]L(\beta) = \prod_i p_i^{y_i}(1-p_i)^{1-y_i}, \qquad \ell(\beta) = \sum_i\big[y_i\log p_i + (1-y_i)\log(1-p_i)\big]

对 β\beta 分别求一阶导(梯度)、二阶导(Hessian),得到

∇ℓ(β)=X⊤(y−p),∇2ℓ(β)=−X⊤WX\nabla\ell(\beta) = X^{\top}(y - p), \qquad \nabla^2\ell(\beta) = -X^{\top}W X

其中 W=diag(pi(1−pi))W = \text{diag}\big(p_i(1-p_i)\big)。因而 β\beta 可以通过下面的迭代(牛顿法 / IRLS)来估计

βk+1=βk+(X⊤WX)−1X⊤(y−p)\beta_{k+1} = \beta_k + (X^{\top}W X)^{-1}X^{\top}(y - p)

每次迭代后计算变化量并与给定阈值比较,判断是否收敛。

三、交替方向乘子法(ADMM)

考虑带等式约束的问题 min⁡x,zf(x)+g(z)\min_{x,z} f(x) + g(z) s.t. Ax+Bz=cAx + Bz = c。定义增广拉格朗日函数:

Lρ(x,z,λ)=f(x)+g(z)+λ⊤(Ax+Bz−c)+ρ2∥Ax+Bz−c∥2L_\rho(x,z,\lambda) = f(x) + g(z) + \lambda^{\top}(Ax + Bz - c) + \frac{\rho}{2}\|Ax + Bz - c\|^2

增广拉格朗日函数可以理解为在普通拉格朗日函数基础上加了一个二次惩罚项,所以该方法是拉格朗日函数法与罚函数法的结合。

标准 ADMM 的优化迭代步骤(交替更新):

xk+1=arg⁡min⁡xLρ(x,zk,λk)x_{k+1} = \arg\min_x L_\rho(x, z_k, \lambda_k)zk+1=arg⁡min⁡zLρ(xk+1,z,λk)z_{k+1} = \arg\min_z L_\rho(x_{k+1}, z, \lambda_k)λk+1=λk+ρ (Axk+1+Bzk+1−c)\lambda_{k+1} = \lambda_k + \rho\,(Ax_{k+1} + Bz_{k+1} - c)

例子:Lasso

Lasso 是带 L1L_1 正则的回归:min⁡β12∥y−Xβ∥2+λ∥β∥1\min_\beta \frac{1}{2}\|y - X\beta\|^2 + \lambda\|\beta\|_1。引入辅助变量 zz,写成约束形式 min⁡12∥y−Xβ∥2+λ∥z∥1\min \frac{1}{2}\|y-X\beta\|^2 + \lambda\|z\|_1 s.t. β−z=0\beta - z = 0。增广拉格朗日函数为

Lρ=12∥y−Xβ∥2+λ∥z∥1+u⊤(β−z)+ρ2∥β−z∥2L_\rho = \frac{1}{2}\|y - X\beta\|^2 + \lambda\|z\|_1 + u^{\top}(\beta - z) + \frac{\rho}{2}\|\beta - z\|^2

对 β\beta 求导并令其为 0,得到 β\beta 更新是一个岭回归式解:

βk+1=(X⊤X+ρI)−1(X⊤y+ρzk−uk)\beta_{k+1} = (X^{\top}X + \rho I)^{-1}\big(X^{\top}y + \rho z_k - u_k\big)

而 zz 的更新公式是软阈值(soft-thresholding)算子:

zk+1=Sλ/ρ(βk+1+uk/ρ),Sκ(a)=sign(a)max⁡(∣a∣−κ, 0)z_{k+1} = S_{\lambda/\rho}\big(\beta_{k+1} + u_k/\rho\big), \qquad S_\kappa(a) = \text{sign}(a)\max(|a| - \kappa,\ 0)

例子:分位数回归

分位数的定义与损失:分位数回归希望找一条直线拟合样本的 τ\tau 分位数。它用的是检验损失(pinball loss):

ρτ(u)=u(τ−1[u<0])={τ u,u≥0(τ−1) u,u<0\rho_\tau(u) = u\big(\tau - \mathbb{1}[u < 0]\big) = \begin{cases} \tau\,u, & u \ge 0 \\ (\tau-1)\,u, & u < 0 \end{cases}

可以理解为:给小于分位数的点(残差为负)一个权重,给大于分位数的点另一个权重(给多数点小权重、少数点大权重)。

几何意义:希望回归曲线之下能够包含 τ\tau 比例的数据点。

优化:对样本 (xi,yi)(x_i, y_i),考虑 min⁡β∑iρτ(yi−xi⊤β)\min_\beta \sum_i \rho_\tau(y_i - x_i^{\top}\beta)。同样引入辅助变量写成约束形式,构造增广拉格朗日函数,按 ADMM 交替更新 β\beta、辅助变量、对偶变量;对 β\beta 子问题求偏导并令其为 0,得到 β\beta 的闭式更新(形如 β=(X⊤X)−1X⊤(⋅)\beta = (X^{\top}X)^{-1}X^{\top}(\cdot)),辅助变量更新则用与 pinball 损失对应的阈值算子。

四、梯度下降家族

梯度下降算法(GD)

给定一组样本 (xi,yi)(x_i, y_i),考虑 min⁡β1n∑iℓ(xi,yi;β)\min_\beta \frac{1}{n}\sum_i \ell(x_i, y_i; \beta)。若只考虑线性模型 y^i=xi⊤β\hat y_i = x_i^{\top}\beta,目标就是均方误差。对当前估计 βk\beta_k,梯度下降的迭代步骤为

βk+1=βk−t ∇f(βk),∇f(βk)=1n∑i=1n∇ℓ(xi,yi;βk)\beta_{k+1} = \beta_k - t\,\nabla f(\beta_k), \qquad \nabla f(\beta_k) = \frac{1}{n}\sum_{i=1}^{n}\nabla\ell(x_i, y_i; \beta_k)

其中 tt 为步长。

步长 tt 的选择:通常用 Backtracking line search 筛选出足够长、又能保证目标函数递减的步长。给定当前估计,若充分下降条件 f(β−t∇f)≤f(β)−αt∥∇f∥2f(\beta - t\nabla f) \le f(\beta) - \alpha t\|\nabla f\|^2 不成立,就令 t←βtt \leftarrow \beta t 缩小步长再试。

随机梯度下降法(SGD)

给定当前估计值 βk\beta_k,随机梯度下降每步随机抽取一个样本 (xi,yi)(x_i, y_i),只用它算梯度来更新:

βk+1=βk−t ∇ℓ(xi,yi;βk)\beta_{k+1} = \beta_k - t\,\nabla\ell(x_i, y_i; \beta_k)

直到事先给定的收敛准则被满足后停止。

GD 与 SGD 对比

  • 梯度下降法:每步用全体样本求平均梯度;
  • 随机梯度下降法:每步随机抽取一个样本估计梯度。

当样本数量巨大时,GD 每步都要算全体样本的平均,计算量极大。SGD 不能保证每一步都让目标函数下降,相比 GD 往往需要更多迭代次数,但每次迭代的计算量小得多。折中办法是小批量(mini-batch):每次随机抽取一小组样本估计梯度,兼顾稳定与速度。

BFGS 与 L-BFGS

BFGS 算法适用条件:当 Hessian 矩阵维度超大时,求逆非常困难、计算复杂度极高。BFGS 通过迭代逼近 Hessian 的逆,即在牛顿迭代过程中,用一个矩阵 BkB_k 来代替 [∇2f(xk)]−1[\nabla^2 f(x_k)]^{-1},每步用梯度差和位置差对 BkB_k 做秩更新,避免直接求逆。

L-BFGS 算法:由于 BFGS 每次都要存储 p×pp\times p 的近似矩阵,内存开销很大,数据维度高时机器可能吃不消。L-BFGS 只保存最近 mm 次迭代的信息(向量对)用于隐式重构更新方向,大幅节省内存,适合超高维问题。

五、总结

  • 大 O / 小 o:OO 是"不超过量级"(∣f∣≤C∣g∣|f|\le C|g|);oo 是"严格更低量级"(f/g→0f/g\to0);o⇒Oo \Rightarrow O。
  • 驻点 ≠ 最优:∇f=0\nabla f=0 可能是局部最优或鞍点;凸函数才能保证局部即全局。
  • 收敛速率:q=1q=1 线性、q=2q=2 二次;牛顿法二次收敛。
  • 牛顿法:xk+1=xk−[∇2f]−1∇fx_{k+1}=x_k-[\nabla^2 f]^{-1}\nabla f;问题是 Hessian 求逆贵/可能不可逆、值不一定递减、依赖初值。
  • 回溯线搜索:t=1t=1 起,不满足充分下降条件 f(x+tΔx)≤f(x)+αt∇f⊤Δxf(x+t\Delta x)\le f(x)+\alpha t\nabla f^\top\Delta x 就 t←βtt\leftarrow\beta t 缩小。
  • 收敛标准:绝对 / 相对 / 调整相对(分母 ∥xk∥+1\|x_k\|+1)。
  • ADMM:增广拉格朗日 = 拉格朗日 + 二次罚项;交替更新 x、z、对偶 λ。Lasso 的 z 更新是软阈值算子;分位数回归用 pinball 损失。
  • GD vs SGD 全样本、稳但慢;SGD 单样本、快但抖、迭代多;mini-batch 折中。
  • BFGS / L-BFGS:用矩阵迭代逼近 Hessian 逆,避免直接求逆;L-BFGS 只存最近 m 次信息省内存。