使用 AutoDiff 库进行自动微分

• 下载源代码 - 56.79 KB

引言

在某些情况下，我们需要使用数学工具和方程来解决现实世界的问题。不幸的是，在许多情况下，我们无法使用高中时学到的简单笔算技术（是的，高中时他们骗了我们！）来解决问题，而需要依赖数值方法。例如，下面的简单问题无法用笔算解决

e^(-x) + x - 2 = 0

许多数值方法依赖于用户计算函数梯度（偏导数）的能力。这是一项繁琐的任务，通常可以由计算机自动完成，我们将演示 AutoDiff 库如何实现这一目标。它允许开发人员专注于他们的业务逻辑——对他们试图解决的问题进行建模，并节省手动开发梯度所浪费的时间。

背景

自动微分（AD）是一项非常古老的技术，可能指代截然不同的事物。但它们都有一个共同点——所有方法都提供函数的某种表示，并允许在任何给定点计算梯度。有关各种 AD 技术的更多信息，请参阅 此处。

不幸的是，许多人对此并不了解，并付出了巨大的努力来高效地计算导数，以便开发基于梯度的算法。一个非常著名的例子是神经网络的反向传播算法。在使用 AD 时，绝对无需开发学习算法。任何 AD 库都可以有效地为您计算梯度，您只需要使用其结果来更新权重。

本文介绍的 `AutoDiff` 库采用了运算符重载方法。表达式（或项）使用对象来表示。运算符重载允许从现有项创建新项。一旦定义了项，该库就可以为变量的给定赋值高效地计算其值或梯度。

Using the Code

快速入门

在您 下载 了该库之后，您可以创建一个 C# 项目，添加引用并开始尝试！首先，让我们提供一个简单的“hello world”介绍。这是一个定义函数并计算其值和梯度的简单程序。

using System;
using AutoDiff;

class Program
{
	static void Main(string[] args)
	{
		// we will use a function of two variables
		Variable x = new Variable();
		Variable y = new Variable();

		// Define our function: func(x, y) = (x + y) * exp(x - y)
		Term func = (x + y) * TermBuilder.Exp(x - y);

		// define the ordering of variables, and a point where 
		// the function will be evaluated/differentiated
		Variable[] vars = { x, y };
		double[] point = { 1, -2 };

		// calculate the value and the gradient at the point (x, y) = (1, -2)
		double eval = func.Evaluate(vars, point);
		double[] gradient = func.Differentiate(vars, point);

		// write the results
		Console.WriteLine("f(1, -2) = " + eval);
		Console.WriteLine("Gradient of f at (1, -2) = 
			({0}, {1})", gradient[0], gradient[1]);
	}
}

输出如下

f(1, -2) = -20.0855369231877
Gradient of f at (1, -2) = (0, 40.1710738463753)

以下是一些需要注意的关键点

• 变量是 `Variable` 类的对象。
• 更广泛地说，我们使用 `Term` 基类的对象来构建我们的函数。变量也是项。
• 我们使用运算符重载和 `TermBuilder` 类从现有项创建新项。请参见上面代码片段中 `func` 的定义。
• 使用扩展方法 `Evaluate` 和 `Differentiate` 分别计算给定项定义的函数的*值*和*梯度*。

AutoDiff 升级版

`Evaluate` 和 `Differentiate` 函数类似于运行解释代码。它们有效，但速度较慢。我们可以通过“编译”我们的函数来加快计算速度，并使用编译后的版本来高效地评估或区分函数。当我们需要多次计算同一函数的梯度时，编译会很有益。而这正是大多数数值算法中发生的情况。这是另一个执行相同操作的程序，但首先会编译函数。

using System;
using AutoDiff;

class Program
{
	static void Main(string[] args)
	{
		// we will use a function of two variables
		Variable x = new Variable();
		Variable y = new Variable();

		// Define our function: func(x, y) = (x + y) * exp(x - y)
		Term func = (x + y) * TermBuilder.Exp(x - y);

		// create a compiled version of func.
		ICompiledTerm compiledFunc = func.Compile(x, y);

		// calculate value and gradient at (x, y) = (1, -2)
		Tuple<double[], double> gradientAndValue =
			compiledFunc.Differentiate(1, -2);

		double eval = gradientAndValue.Item2;
		double[] gradient = gradientAndValue.Item1;

		// write the results
		Console.WriteLine("f(1, -2) = " + eval);
		Console.WriteLine("Gradient of f at 
			(1, -2) = ({0}, {1})", gradient[0], gradient[1]);
	}
}

普通项和编译项之间有一个小区别。它们的 `Differentiate` 方法返回一个包含值和梯度的元组。原因有两个

• 计算编译形式的梯度也会评估函数。那么为什么不返回结果呢？
• 大多数数值算法在给定点都需要函数的值和梯度。因此，这只是节省了一次额外的 `Evaluate` 调用。

方程求解示例

让我们创建一个程序来求解上面给出的方程

e^(-x) + x - 2 = 0

我们将使用众所周知的 牛顿-拉夫逊 方程求解方法。我们将称之为 NR 方法。首先，让我们创建一个小型方法，给定一个函数 *f*，求解方程 *f*( *x* ) = 0。我们的函数将以编译后的项形式给出。NR 方法需要一个初始值 `x0`，该值足够接近真实解，还需要一个停止条件。在下面的示例中，我们将使用迭代次数作为停止条件。

static double NewtonRaphson(ICompiledTerm function, double x0, int iterationsCount = 10)
{
	double x = x0;
	for(int i = 0; i < iterationsCount; ++i)
	{
		var valueAndGradient = function.Differentiate(x);

		// take the function's value - f(x)
		var fx = valueAndGradient.Item2;

		// we are dealing with single-variable functions. 
		// Therefore the first gradient element is the derivative
		var dfx = valueAndGradient.Item1[0];

		x = x - fx / dfx;
	}

	return x;
}

现在我们可以使用此方法来求解我们的方程。首先，让我们绘制 `f(x) = e^(-x) + x - 2` 的图

我们可以看到我们的方程有两个解——一个在 `x = -1` 附近，另一个在 `x = 2` 附近。因此，我们可以编写一个小型程序来求解我们的方程。

static void Main(string[] args)
{
	// define and compile our function
	var x = new Variable();
	var func = TermBuilder.Exp(-x) + x - 2;
	var compiledFunc = func.Compile(x);

	// compute the two solutions
	var solution1 = NewtonRaphson(compiledFunc, -1);
	var solution2 = NewtonRaphson(compiledFunc, 2);

	Console.WriteLine("X1 = {0}, X2 = {1}", solution1, solution2);
}

这是输出。您可以使用您喜欢的计算器来验证解的正确性

X1 = -1.14619322062058, X2 = 1.84140566043696

更多信息

该库最初是为了在数值优化算法中使用而编写的。大多数此类算法都依赖于用户能够提供函数值及其在任何给定点的梯度。算法种类繁多（BFGS、L-BFGS、共轭梯度、梯度下降），但它们都需要梯度。它已成功用于我的论文研究项目，并取得了优异的结果。

它被创建出来是为了让开发人员能够专注于对目标函数和约束进行建模，以解决他们的现实世界问题，从而使他们免于执行繁琐的数学任务。将 AD 与优化库结合使用，开发人员可以快速轻松地尝试许多不同的函数。该库将此能力带到了 .NET 平台。

您可以下载该库并在您的研究项目或实际应用中使用它，欢迎您提供反馈。我们需要反馈！

一些算法细节

该库使用反向累积算法。该算法允许以线性时间计算函数梯度，其复杂度与函数的大小成正比。更准确地说，我们可以用一个 DAG（有向无环图）来定义我们的函数，它描述了它的计算方式。这类似于编译器的语法树。终端节点是变量和常量，数学运算符是内部节点。它非常类似于用于 LINQ 的 .NET 表达式树。该算法以线性时间计算梯度，其复杂度与此图的大小（节点数+边数）成正比。

令人惊讶的是，线性时间梯度计算实际上比数值近似梯度*更快*。当我们进行近似时，我们需要为每个变量评估一次函数（通过稍微改变其值，然后评估函数）。这意味着梯度的近似计算是二次时间的——每次评估需要线性时间，而我们有 *n* 次这样的评估。您将获得免费午餐——既能获得实际计算梯度的精度，又能获得反向累积算法的速度。

不幸的是，反向累积不适合单变量函数。您之前看到的 NR 示例并不是该库的良好用途，而前向模式 AD 库是更好的选择。然而，NR 示例对于本文来说足够简单，因此被选中。该库的真正性能和易用性优势在于多变量优化算法。

历史

• 2011年4月8日：初始版本