如何使用Numpy计算导数？

Question

例如，如何计算函数的导数

y = x²+1

使用numpy？

比方说，我想要导数在x=5时的值。。。

Answer 1

我能想到的最直接的方法是使用numpy's gradient function：

x = numpy.linspace(0,10,1000)
dx = x[1]-x[0]
y = x**2 + 1
dydx = numpy.gradient(y, dx)

这样，dydx将使用中心差分计算，并且具有与y相同的长度，而numpy.diff使用前向差分并返回（n-1）大小向量。

Answer 2

NumPy不提供计算导数的一般功能。但它可以处理多项式的简单特例：

>>> p = numpy.poly1d([1, 0, 1])
>>> print p
   2
1 x + 1
>>> q = p.deriv()
>>> print q
2 x
>>> q(5)
10

如果你想用数值计算导数，你可以在绝大多数应用中使用中心差分商。对于单点的导数，公式如下

x = 5.0
eps = numpy.sqrt(numpy.finfo(float).eps) * (1.0 + x)
print (p(x + eps) - p(x - eps)) / (2.0 * eps * x)

如果有一个横坐标数组x和相应的函数值数组y，则可以计算

numpy.diff(y) / numpy.diff(x)

Answer 3

你有四个选择

1. Finite Differences
2. Automatic Derivatives
3. Symbolic Differentiation
4. 手工计算导数。

有限差分不需要外部工具，但很容易出现数值误差，如果处于多元情况下，可能需要一段时间。

如果你的问题足够简单的话，符号微分是理想的。如今，符号方法变得相当强大。SymPy是一个很好的项目，可以很好地与NumPy集成。看看autowrap或lambdify函数，或者查看Jensen's blogpost about a similar question。

自动衍生工具非常酷，不容易出现数字错误，但确实需要一些额外的库（google提供了一些不错的选择）。这是最稳健但也是最复杂/最难设置的选择。如果您能很好地限制自己使用numpy语法，那么Theano可能是一个不错的选择。

下面是一个使用SymPy的示例

In [1]: from sympy import *
In [2]: import numpy as np
In [3]: x = Symbol('x')
In [4]: y = x**2 + 1
In [5]: yprime = y.diff(x)
In [6]: yprime
Out[6]: 2⋅x

In [7]: f = lambdify(x, yprime, 'numpy')
In [8]: f(np.ones(5))
Out[8]: [ 2.  2.  2.  2.  2.]