Python | 从零实现10种优化算法并比较

🚀SGD

算法原理

随机梯度下降（Stochastic Gradient Descent，SGD）是一种基于梯度的优化算法，用于在机器学习和深度学习中最小化目标函数。以下是其基本原理和参数更新数学公式的简要说明：

• 基本原理： SGD旨在通过迭代的方式，沿着目标函数梯度的反方向更新模型参数，以逐步找到使目标函数最小化的最优参数值。

• 在每次迭代中，随机选取一个训练样本，计算该样本对应的损失函数关于模型参数的梯度，然后根据这个梯度来更新模型参数。通过不断地重复这个过程，模型参数逐渐向最优解靠近，使得损失函数的值不断减小。

• 参数更新数学公式：设 $\theta$ 为模型参数，$J(\theta )$为目标函数，$\eta$为学习率，在第$t$次迭代时，随机选取一个样本$(x_i,y_i)$，计算该样本损失函数$L(y_i,f(x_i;\theta ))$关于参数$\theta$的梯度$\nabla L(y_i,f(x_i;\theta ))$，则参数更新公式为：
  ${\theta }_{t+1}={\theta }_t-\eta \nabla L(y_i,f(x_i;{\theta }_t))$
  其中，${\theta }_t$是第$t$次迭代时的参数值，${\theta }_{t+1}$是第$t+1$次迭代时的参数值。学习率$\eta$是一个超参数，用于控制参数更新的步长，通常需要根据具体问题进行调整。

代码实现

本文这里实际实现的是GD算法，因为样例优化问题比较简单：$f(x,y)=3x^2+4y^2$

class SGD:
    def __init__(self, lr=0.01):
        self.lr = lr

    def update(self, w, grad_w):
        return w - self.lr * grad_w

可视效果

优化的损失函数为：$f(x,y)=3x^2+4y^2$

🚀 Adagrad

算法原理

• 基本原理：Adagrad根据每个参数的历史梯度平方和来调整学习率，对于频繁出现的参数，其学习率会逐渐减小，而对于不常出现的参数，学习率则相对较大，从而实现了自适应学习率的调整。
• 参数更新数学公式：
  $G_{t,ii}={\sum }_{k=1}^t{g}_{k,i}^2$
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{G_{t,ii}+ϵ}}{g}_t$
  其中，$g_{k,i}$是第$k$步参数$i$的梯度，$G_{t,ii}$是参数$i$到第$t$步为止的梯度平方和，$\eta$是初始学习率，$ϵ$是一个小常数，用于防止分母为零。

代码实现

import numpy as np

class Adagrad:
    def __init__(self, eps=1e-8):
        self.eps = eps
        self.G = None

    def update(self, w, grad_w):
        if self.G is None:
            self.G = np.zeros_like(grad_w)

        self.G += grad_w ** 2
        dx = grad_w / (np.sqrt(self.G) + self.eps)
        return w - dx

可视效果

优化的损失函数为：$f(x,y)=3x^2+4y^2$

🚀 Adadelta

算法原理

• 基本原理：Adadelta是Adagrad的一种改进算法，旨在解决Adagrad学习率单调递减的问题。它通过引入一个随时间衰减的平均梯度平方项和一个自适应的学习率来调整参数更新步长，使得在训练过程中学习率能够根据问题的复杂度自适应调整。
• 参数更新数学公式：
  $E[g^2{]}_t=\rho E[g^2{]}_{t-1}+(1-\rho )g_t^2$
  $\Delta {\theta }_t=-\frac{\sqrt{E[\Delta {\theta }^2{]}_{t-1}+ϵ}}{\sqrt{E[g^2{]}_t+ϵ}}{g}_t$
  ${\theta }_{t+1}={\theta }_t+\Delta {\theta }_t$
  其中，$g_t$是第$t$步的梯度，$E[g^2{]}_t$是梯度平方的指数移动平均，$\rho$是衰减率，$\Delta {\theta }_t$是参数更新量，$E[\Delta {\theta }^2{]}_{t-1}$是参数更新量平方的指数移动平均，$ϵ$是一个小常数，用于防止分母为零。

代码实现

import numpy as np

class Adadelta:
    def __init__(self, rho=0.95, eps=1e-6):
        self.rho = rho
        self.eps = eps
        self.Eg2 = None
        self.Edx2 = None
        self.delta_x = None

    def update(self, w, grad_w):
        if self.Eg2 is None:
            self.Eg2 = np.zeros_like(grad_w)
            self.Edx2 = np.zeros_like(w)
            self.delta_x = np.zeros_like(w)

        self.Eg2 = self.rho * self.Eg2 + (1 - self.rho) * grad_w ** 2
        dx = np.sqrt((self.Edx2 + self.eps) / (self.Eg2 + self.eps)) * grad_w
        self.Edx2 = self.rho * self.Edx2 + (1 - self.rho) * dx ** 2
        self.delta_x = -dx
        return w + self.delta_x

可视效果

🚀 RMSProp

算法原理

• 基本原理：RMSProp通过对梯度平方的指数移动平均来调整学习率，使得学习率能够根据参数的重要性和梯度的变化情况自适应调整，在一定程度上解决了Adagrad中学习率单调递减的问题，提高了训练的稳定性和收敛速度。
• 参数更新数学公式：
  $E[g^2{]}_t=\rho E[g^2{]}_{t-1}+(1-\rho )g_t^2$
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{E[g^2{]}_t+ϵ}}{g}_t$
  其中，各参数含义与Adadelta中的类似。

代码实现

import numpy as np

class RMSProp:
    def __init__(self, lr=0.01, rho=0.9, eps=1e-8):
        self.lr = lr
        self.rho = rho
        self.eps = eps
        self.Eg2 = None

    def update(self, w, grad_w):
        if self.Eg2 is None:
            self.Eg2 = np.zeros_like(grad_w)

        self.Eg2 = self.rho * self.Eg2 + (1 - self.rho) * grad_w ** 2
        dx = self.lr * grad_w / (np.sqrt(self.Eg2) + self.eps)
        return w - dx

可视效果

🚀 Adam

算法原理

• 基本原理：Adam结合了Adagrad和RMSProp的优点，通过计算梯度的一阶矩估计（均值）和二阶矩估计（方差）来调整学习率，能够在训练初期利用较大的学习率快速收敛，在后期根据梯度的变化情况自适应地调整学习率，从而提高训练的稳定性和收敛速度。
• 参数更新数学公式：
  $m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t$
  $v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2$
  ${\hat{m}}_t=\frac{m_t}{1-{\beta }_1^t}$
  ${\hat{v}}_t=\frac{v_t}{1-{\beta }_2^t}$
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{\hat{v}}_t}+ϵ}{\hat{m}}_t$
  其中，$m_t$是一阶矩估计，$v_t$是二阶矩估计，${\beta }_1$和${\beta }_2$是指数衰减率，${\hat{m}}_t$和${\hat{v}}_t$是偏差修正后的一阶矩估计和二阶矩估计，$\eta$是学习率，$ϵ$是一个小常数，用于防止分母为零。

代码实现

import numpy as np

class Adam:
    def __init__(self, lr=0.01, beta1=0.9, beta2=0.999, eps=1e-8):
        self.lr = lr
        self.beta1 = beta1
        self.beta2 = beta2
        self.eps = eps
        self.m = None
        self.v = None
        self.t = 0

    def update(self, w, grad_w):
        self.t += 1
        if self.m is None:
            self.m = np.zeros_like(grad_w)
            self.v = np.zeros_like(grad_w)

        self.m = self.beta1 * self.m + (1 - self.beta1) * grad_w
        self.v = self.beta2 * self.v + (1 - self.beta2) * grad_w ** 2

        m_hat = self.m / (1 - self.beta1 ** self.t)
        v_hat = self.v / (1 - self.beta2 ** self.t)

        dx = self.lr * m_hat / (np.sqrt(v_hat) + self.eps)
        return w - dx

可视效果

🚀 AdamW

算法原理

• 基本原理：AdamW在Adam的基础上，对权重衰减进行了更合理的处理，将权重衰减从优化器的参数更新中分离出来，直接作用于参数本身，而不是像传统的方法那样将权重衰减与学习率耦合在一起，从而在一定程度上提高了模型的泛化能力和训练效果。
• 参数更新数学公式：与Adam类似，只是在更新参数时多了一项权重衰减项：
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{\hat{v}}_t}+ϵ}{\hat{m}}_t-\lambda {\theta }_t$
  其中，$\lambda$是权重衰减系数，其他参数含义与Adam中相同。

代码实现

import numpy as np

class AdamW:
    def __init__(self, lr=0.001, beta1=0.9, beta2=0.999, eps=1e-8, weight_decay=0.01):
        self.lr = lr
        self.beta1 = beta1
        self.beta2 = beta2
        self.eps = eps
        self.weight_decay = weight_decay
        self.m = None
        self.v = None
        self.t = 0

    def update(self, w, grad_w):
        self.t += 1
        if self.m is None:
            self.m = np.zeros_like(grad_w)
            self.v = np.zeros_like(grad_w)

        self.m = self.beta1 * self.m + (1 - self.beta1) * grad_w
        self.v = self.beta2 * self.v + (1 - self.beta2) * grad_w ** 2

        m_hat = self.m / (1 - self.beta1 ** self.t)
        v_hat = self.v / (1 - self.beta2 ** self.t)

        # 加入权重衰减项
        update = (self.lr * m_hat / (np.sqrt(v_hat) + self.eps)) + self.weight_decay * w
        return w - update

可视效果

🚀 SparseAdam

算法原理

• 基本原理：SparseAdam是Adam的一种变体，专门用于处理稀疏数据。它在计算梯度的一阶矩估计和二阶矩估计时，只对非零梯度进行更新，从而提高了计算效率，同时也能在一定程度上保证模型的训练效果。
• 参数更新数学公式：与Adam类似，但在计算$m_t$和$v_t$时只考虑非零梯度：
  $m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t$
  $v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2$
  ${\hat{m}}_t=\frac{m_t}{1-{\beta }_1^t}$
  ${\hat{v}}_t=\frac{v_t}{1-{\beta }_2^t}$
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{\hat{v}}_t}+ϵ}{\hat{m}}_t$
  其中，$g_t$为稀疏梯度，其他参数含义与Adam中相同。

代码实现

import numpy as np

class SparseAdam:
    def __init__(self, lr=0.01, beta1=0.9, beta2=0.999, eps=1e-8):
        self.lr = lr
        self.beta1 = beta1
        self.beta2 = beta2
        self.eps = eps
        self.m = None
        self.v = None
        self.t = 0

    def update(self, w, grad_w):
        self.t += 1
        if self.m is None:
            self.m = {}
            self.v = {}

        for index in np.ndindex(grad_w.shape):
            if index not in self.m:
                self.m[index] = np.zeros(1)
                self.v[index] = np.zeros(1)

            g = grad_w[index]
            self.m[index] = self.beta1 * self.m[index] + (1 - self.beta1) * g
            self.v[index] = self.beta2 * self.v[index] + (1 - self.beta2) * g ** 2

            m_hat = self.m[index] / (1 - self.beta1 ** self.t)
            v_hat = self.v[index] / (1 - self.beta2 ** self.t)

            dx = self.lr * m_hat / (np.sqrt(v_hat) + self.eps)
            w[index] -= dx

        return w

可视效果

🚀 Adamax

算法原理

• 基本原理：Adamax是Adam的一种变体，它对Adam中的二阶矩估计进行了修改，使用了无穷范数来代替原来的$L_2$范数，在一些情况下可能具有更好的性能和稳定性。
• 参数更新数学公式：
  $m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t$
  $u_t=\max ({\beta }_2{u}_{t-1},|g_t|)$
  ${\hat{m}}_t=\frac{m_t}{1-{\beta }_1^t}$
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{u_t}{\hat{m}}_t$
  其中，$m_t$、${\beta }_1$、${\hat{m}}_t$的含义与Adam中相同，$u_t$是无穷范数的指数移动平均，$|g_t|$是第$t$步梯度的绝对值，$\eta$是学习率。

代码实现

import numpy as np

class Adamax:
    def __init__(self, lr=0.01, beta1=0.9, beta2=0.999, eps=1e-8):
        self.lr = lr
        self.beta1 = beta1
        self.beta2 = beta2
        self.eps = eps
        self.m = None
        self.u = None
        self.t = 0

    def update(self, w, grad_w):
        self.t += 1
        if self.m is None:
            self.m = np.zeros_like(grad_w)
            self.u = np.zeros_like(grad_w)

        self.m = self.beta1 * self.m + (1 - self.beta1) * grad_w
        self.u = np.maximum(self.beta2 * self.u, np.abs(grad_w))

        m_hat = self.m / (1 - self.beta1 ** self.t)
        dx = self.lr * m_hat / (self.u + self.eps)
        return w - dx

可视效果

🚀 NAdam

算法原理

• 基本原理：NAdam是在Adam的基础上结合了Nesterov加速梯度（NAG）的思想，通过对未来梯度的近似估计来调整参数更新的方向，使得参数更新能够更快速地向最优解收敛，进一步提高了训练的效率和收敛速度。
• 参数更新数学公式：
  $m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t$
  ${\hat{m}}_t=\frac{m_t}{1-{\beta }_1^t}$
  $g_t^{\prime }=\frac{{\beta }_1{m}_{t-1}}{1-{\beta }_1^t}+(1-{\beta }_1)g_t$
  $v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2$
  ${\hat{v}}_t=\frac{v_t}{1-{\beta }_2^t}$
  ${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{\hat{v}}_t}+ϵ}{g}_t^{\prime }$
  其中，各参数含义与Adam中类似。

代码实现

import numpy as np

class NAdam:
    def __init__(self, lr=0.01, beta1=0.9, beta2=0.999, eps=1e-8):
        self.lr = lr
        self.beta1 = beta1
        self.beta2 = beta2
        self.eps = eps
        self.m = None
        self.v = None
        self.t = 0

    def update(self, w, grad_w):
        self.t += 1
        if self.m is None:
            self.m = np.zeros_like(grad_w)
            self.v = np.zeros_like(grad_w)

        self.m = self.beta1 * self.m + (1 - self.beta1) * grad_w
        self.v = self.beta2 * self.v + (1 - self.beta2) * grad_w ** 2

        m_hat = self.beta1 * self.m / (1 - self.beta1 ** self.t) + (1 - self.beta1) * grad_w / (1 - self.beta1 ** self.t)
        v_hat = self.v / (1 - self.beta2 ** self.t)

        dx = self.lr * m_hat / (np.sqrt(v_hat) + self.eps)
        return w - dx

可视效果

🚀 RProp

算法原理

• 基本原理：RProp基于梯度的符号来调整参数更新的步长，而不依赖于梯度的大小。如果当前梯度与上一次梯度的符号相同，则增加步长；如果符号相反，则减小步长。这种方法可以在一定程度上避免梯度消失或爆炸问题，提高训练的稳定性和收敛速度。
• 参数更新数学公式：
  $\Delta {\theta }_t=\{\begin{array}{ll}-{\eta }^{+}\Delta {\theta }_{t-1},& g_t{g}_{t-1}>0\\ -{\eta }^{-}\Delta {\theta }_{t-1},& g_t{g}_{t-1}<0\\ \Delta {\theta }_{t-1},& g_t{g}_{t-1}=0\end{array}$
  ${\theta }_{t+1}={\theta }_t+\Delta {\theta }_t$
  其中，${\eta }^{+}$和${\eta }^{-}$分别是步长增加和减小的因子，$g_t$是第$t$步的梯度，$\Delta {\theta }_t$是参数更新量。

代码实现

import numpy as np

class RProp:
    def __init__(self, lr_plus=1.2, lr_minus=0.5, delta_max=50, delta_min=1e-6, initial_delta=0.1):
        self.lr_plus = lr_plus
        self.lr_minus = lr_minus
        self.delta_max = delta_max
        self.delta_min = delta_min
        self.initial_delta = initial_delta
        self.delta = None
        self.sign_grad_prev = None

    def update(self, w, grad_w):
        if self.delta is None:
            self.delta = np.ones_like(grad_w) * self.initial_delta
            self.sign_grad_prev = np.sign(grad_w)

        sign_grad = np.sign(grad_w)
        self.delta *= (sign_grad == self.sign_grad_prev) * self.lr_plus + (sign_grad!= self.sign_grad_prev) * self.lr_minus
        self.delta = np.minimum(self.delta_max, np.maximum(self.delta_min, self.delta))
        self.sign_grad_prev = sign_grad

        dx = -self.delta * sign_grad
        return w + dx

可视效果

🚀 代码结构

github项目链接

preject
	|--optim
	|--result
	|--main.py

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from optim import Adadelta, Adagrad, Adam, Adamax, AdamW, NAdam, RMSProp, RProp, SGD, SparseAdam

# 定义二维损失函数
def loss_function(x, y):
    return 3 * x ** 2 + 4 * y ** 2

# 计算损失函数关于x和y的偏导数（针对单个样本的情况）
def gradient(x, y):
    grad_x = 6 * x
    grad_y = 8 * y
    return np.array([grad_x, grad_y])

def train(optimizer_name, num_iterations=500):
    # 初始参数
    params = np.array([4.5, 4.5])
    x, y = params[0], params[1]
    # 创建优化器实例
    if optimizer_name == "Adadelta":
        optimizer = Adadelta()
    elif optimizer_name == "Adagrad":
        optimizer = Adagrad()
    elif optimizer_name == "Adam":
        optimizer = Adam()
    elif optimizer_name == "Adamax":
        optimizer = Adamax()
    elif optimizer_name == "AdamW":
        optimizer = AdamW()
    elif optimizer_name == "NAdam":
        optimizer = NAdam()
    elif optimizer_name == "RMSProp":
        optimizer = RMSProp()
    elif optimizer_name == "RProp":
        optimizer = RProp()
    elif optimizer_name == "SGD":
        optimizer = SGD()
    elif optimizer_name == "SparseAdam":
        optimizer = SparseAdam()
    else:
        raise Exception("Not Support This Optimizer!")

    # 用于存储每次迭代的x、y和损失值，方便可视化
    x_history = []
    y_history = []
    loss_history = []

    for _ in range(num_iterations):
        grads = gradient(x, y)
        params = np.array([x, y])
        updated_params = optimizer.update(params, grads)
        x, y = updated_params[0], updated_params[1]
        loss = loss_function(x, y)
        x_history.append(x)
        y_history.append(y)
        loss_history.append(loss)

    # 创建3D图形对象
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')

    # 生成用于绘制损失函数曲面的网格数据
    x_range = np.linspace(-5, 5, 100)
    y_range = np.linspace(-5, 5, 100)
    X, Y = np.meshgrid(x_range, y_range)
    Z = loss_function(X, Y)

    # 绘制损失函数曲面
    ax.plot_surface(X, Y, Z, alpha=0.5, cmap='viridis')

    # 绘制参数x、y变化以及对应损失函数值的散点图
    ax.plot(x_history, y_history, loss_history, c='r', marker='.', label=f'{optimizer_name} Optimization Path')

    # 设置坐标轴标签
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.set_zlabel('Loss')

    # 设置图形标题
    ax.set_title(f'Parameter and Loss Change with {optimizer_name}')

    # 添加图例
    ax.legend()

    # 显示图形
    plt.savefig(f'result\\{optimizer_name}.png')
    plt.close()

    # 可视化损失下降过程
    plt.plot(loss_history)
    plt.xlabel('Iteration')
    plt.ylabel('Loss')
    plt.title(f'Loss Decrease with {optimizer_name}')
    # 显示图形
    plt.savefig(f'result\\{optimizer_name}_loss.png')
    plt.close()

if __name__ == "__main__":
    names = ["Adadelta", "Adagrad", "Adam", "Adamax", "AdamW", "NAdam", "RMSProp", "RProp", "SGD", "SparseAdam"]
    for name in names:
        train(optimizer_name=name)