实现二阶导数的自动微分：遍历计算图的算法？

Question

我正在尝试为一个Python统计包实现自动微分（这个问题的形式跟优化问题的形式很像）。

我通过重载操作符和工厂函数来生成计算图，比如用sum()、exp()等操作。我已经实现了反向累积的方法来计算梯度的自动微分。不过，我发现实现二阶导数（也就是Hessian）要困难得多。我知道怎么计算单独的二阶偏导数，但在遍历计算图和进行累积时遇到了麻烦。有没有人知道好的文章，能提供关于二阶导数的自动微分算法，或者有没有开源库可以让我学习一下？

Answer 1

在三维空间中，通常用来近似海森矩阵的方法是BFGS。

L-BFGS方法和BFGS很相似。

你可以在这里找到L-BFGS的源代码（这个代码会计算海森矩阵，作为解决常微分方程的中间结果），支持多种编程语言（比如C#、C++、VBA等），不过没有Python版本。我觉得把它翻译成Python不太容易。

如果你打算把算法从其他语言翻译过来，特别要注意数值误差，并进行敏感性分析（你需要计算海森矩阵的逆矩阵）。

Answer 2

首先，你需要决定是想计算一个稀疏的Hessian矩阵，还是更接近于一个完全密集的Hessian矩阵。

如果你想要稀疏的Hessian矩阵，现在有两种比较有效的方法。第一种是聪明地利用计算图，通过一次反向遍历计算图，你可以使用边推算法来计算Hessian矩阵：

http://www.tandfonline.com/doi/full/10.1080/10556788.2011.580098

另一种方法是尝试图着色技术，将你的Hessian矩阵压缩成更少列的矩阵，然后使用反向累积来计算每一列。

http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.66.2603

如果你想要的是一个密集的Hessian矩阵（在实际中不太常见），那么你可能更适合一次计算Hessian的一列，使用反向累积的方法（可以搜索BRUCE CHRISTIANSON和反向累积）。