自动微积分

发布于 2024-06-11 00:15:11 字数 7357 浏览 17 评论 0

几种微分求解的方法

手动求解法(Manual Differentiation)
数值微分法(Numerical Differentiation)
符号微分法(Symbolic Differentiation)
自动微分法(Automatic Differentiation)

自动微分

在数学和计算代数领域，automatic differentiation (AD) 又称为 algorithmic differentiation 或者 computational differentiation。AD 是一个可以对程序代码表示的数学函数进行自动微分的技术。AD 利用链式法则来达到自动求解的目录，AD 有两种主要的方法：

代码转换（source-code transformation）（R. Giering and T. Kaminski. 1998）：
- 利用一个代码转换编译器，这个编译器会分析源代码，然后产生一个和源代码对应的伴随模式(adjoint model) 程序，编译时的代码生成（如用 flex-bison 做词法、语法分析）；
- 优点是静态生成效率高（原始算法的 3~4 倍) ;
- 一次生成，多次使用，缺点是学习门槛较高（编译原理…）；
- 很多比较好的工具非免费；
- 对现代编程语言特性的限制（如 C++类、模板等）；
运算符重载(operator overloading)
- 应用比较广泛，很多编程语言特性可以很好的工作；
- 优点是简单直接，缺点是动态生成成本较高（代表性的工具效率是原始算法的 10~35 倍）。
- 较多免费开源 C++ 工具 (e.g. ADOL-C, CppAD, Sacado)；

AD 这两种实现方式：运算符重载与代码生成，两种方式的原理都一样， 链式法则 。AD 相关工具，请到这个 http://www.autodiff.org/ 页面。自动微分（AD）是计算导数的最优方法，比符号计算、有限微分更快更精确，AD 已经广泛应用在优化领域，包括人工神经网络的训练算法 back-propagation（BP）等。

AD 基本原理

链式法则

AD 基本原理是 链式法则 。链式法则又分为正向和反向。如下例子

$$y=f(g(h(x)))=f(g(h(w_0)))=f(g(w_1))=f(w_2)=w_3$$

那么链式法则就表示为:

$$\frac{dy}{dx}=\frac{dy}{dw_2}\frac{dw_2}{dw_1}\frac{dw_1}{dx}$$

正向链式法则：则从链里往外算（即，先算$\frac{dw_1}{dx}$，再算$\frac{dw_2}{dw_1}$，最后算$\frac{dy}{dw_2}$）。

反向链式法则：则从链外往外里（即，先算$\frac{dy}{dw_2}$，再算$\frac{dw_2}{dw_1}$，最后算$\frac{dw_1}{dx}$）。

简洁的表示为：

正向链式法则：计算递归式$$\frac{dw_i}{dx}=\frac{dw_i}{dw_{i-1}}\frac{dw_{i-1}}{dx}，且 w_3=y$$

反向链式法则：计算递归式$$\frac{dy}{dw_i}=\frac{dy}{dw_{i+1}}\frac{dw_{i+1}}{dw_i}，且 w_0=y$$

通常，正向跟反向链式法则都是通过计算图来表示和计算，这样更加的方便。

正向链式法则计算图

将函数转化为一个 DAG（有向无环图），就能很容易的求解每一步的值。

例 1：

假设有计算式子：$y=x_0*x_1+x_1$ (1)，我们要求解$\frac{\partial y}{\partial x_1}$。

首先，其计算图表示为：

计算图 1.png

把（1）式展开为链式计算序列为：

c=x0*x1
d=c+x1
y=d

那么其对应的计算图则为：

计算图 2.png

这里假设$x_0,x_1$的初值分别为：1, 2。 (1) 式展开计算子序列为：

x0=1
x1=2
c=x0*x1
d=c+x1
y=d

那么其对应计算图为：

计算图 3.png

然后对该计算图求解偏导数：

计算图 4.png

由链式法则得到：$$\frac{\partial y}{\partial x_1}=\frac{\partial y}{\partial d}\frac{\partial d}{\partial c}\frac{\partial c}{\partial x_1}+\frac{\partial y}{\partial d}\frac{\partial d}{\partial x_1}$$

结合计算图中对应的值可知：$\frac{\partial y}{\partial x_1}=1\times1\times1+1\times1=2$

例 2：

对于下列函数：$$f(x_1,x_2)=ln(x_1)+x_1*x_2-sin(x_2)$$

转化为计算图：

计算图 5.png

那么求每一步的导数值就可以表示为：

计算图 6.png

上表，左半部分是从左往右计算图每个节点的求值结果，右半部分是每个节点对于$x_1$的求导结果，比如$\dot{v_1}=\frac{dv}{dx_1}$，注意到每一步的求导都会利用到上一步的求导结果。

对于自动微分的正向模式，如果函数输入输出为：$R \rightarrow R^m $

那么正向模式只需要计算一次上表右侧过程即可，非常高效。但是对于输入输出映射为：$R^n\rightarrow R^m$，这样一个有$n$个输入的函数，对于函数的梯度求解则需要处理$n$遍上述过程。而且再实际算法模型中，通常输入输出是极度不成比例的，也就是$n>>m$，那么利用正向模式进行自动微分的效率就太低了，因此有了反向模式的出现。

反向链式法则计算图

自动微分的反向模式其实就是一种通用的 BackPropagation ( 反向传播算法 )，即 backpropagation 是自动微分反向模式的一种特殊形式。

反向模式从最终结果开始求导，利用最终输出对每一个节点进行求导，其过程如下计算图所示：

计算图 7.png

其具体计算过程如下表所示：

计算图 8.png

上表左边和之前的正向模式一致，用于求解函数值，右边则是反向模式的计算过程，须从下向上看，也就是一开始先计算输出$y$对于节点$v_5$的导数，用$\bar{v_5}=\frac{dy}{dv_5}$，这样的记号可以强调我们对当前计算结果进行缓存，以便用于后续计算，而不必重复计算。再由链式法则我们可以计算输出节点对于图种每个节点的导数。

比如对节点$v_3$：

$$\frac{y}{v_3}=\frac{dy}{dv_5}\frac{dv_5}{dv_3}$$

用$\bar{v_i}=\frac{dy}{dv_i}$记法，则有：

$$\frac{y}{v_3}=\bar{v_5}\frac{dv_5}{dv_3}$$

比如对于节点$v_0$：

$$\frac{y}{v_0}=\frac{dy}{dv_2}\frac{dv_2}{dv_0}+\frac{dy}{dv_3}\frac{dv_3}{dv_0}$$

用$\bar{v_i}=\frac{dy}{dv_i}$记法，则有：

$$\frac{y}{v_0}=\bar{v_2}\frac{dv_2}{dv_0}+\bar{v_3}\frac{dv_3}{dv_0}$$

和 backpropagation 算法一样，我们必须记住前向时当前节点发出的边，然后在反向传播时，可以搜集所有受到当前节点影响的节点。

如上的计算过程，对于像神经网络这种模型，输入通常是上万到上百万维，而输出损失函数是 1 维的模型，则只需要一遍反向模式计算过程，便可以求出输出对于各个输入的导数，从而轻松求取梯度用于后续优化更新。