如何使用Numpy计算导数？

Question

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

y = x ² +1

用numpy？

比方说,我希望x = 5时导数的值...

Answer 1

你有四个选择

1. 您可以使用有限差异
2. 您可以使用自动衍生品
3. 您可以使用符号差异化
4. 您可以手动计算衍生物.

有限差异不需要外部工具,但容易出现数值误差,如果您处于多变量情况,可能需要一段时间.

如果您的问题足够简单,符号差异是理想的选择.如今,符号方法变得相当强大.SymPy是一个很好的项目,与NumPy很好地集成.查看autowrap或lambdify函数或查看Jensen关于类似问题的博文.

自动衍生品非常酷,不容易出现数字错误,但确实需要一些额外的库(谷歌为此,有一些不错的选择).这是最强大但也是最复杂/最难设置的选择.如果你很好地限制自己的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.]

Answer 2

我能想到的最直接的方法是使用numpy的渐变函数:

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 3

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函数值数组的abscissae数组,你可以计算导数的近似值

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

Answer 4

假设您想使用numpy,您可以使用严格的定义在任何点数字计算函数的导数 :

def d_fun(x):
    h = 1e-5 #in theory h is an infinitesimal
    return (fun(x+h)-fun(x))/h

您还可以使用对称导数获得更好的结果:

def d_fun(x):
    h = 1e-5
    return (fun(x+h)-fun(x-h))/(2*h)

使用您的示例,完整代码应如下所示:

def fun(x):
    return x**2 + 1

def d_fun(x):
    h = 1e-5
    return (fun(x+h)-fun(x-h))/(2*h)

现在,您可以在数字上找到衍生物x=5:

In [1]: d_fun(5)
Out[1]: 9.999999999621423

Answer 5

您可以使用scipy，这非常简单：

scipy.misc.derivative(func, x0, dx=1.0, n=1, args=(), order=3)

求函数在一点处的 n 阶导数。

在你的情况下：

from scipy.misc import derivative

def f(x):
    return x**2 + 1

derivative(f, 5, dx=1e-6)
# 10.00000000139778

Answer 6

我会在堆上扔另一种方法......

scipy.interpolate许多插值样条都能够提供衍生物.因此,使用线性样条(k=1),样条的导数(使用该derivative()方法)应等于前向差异.我不完全确定,但我相信使用三次样条导数类似于居中差分导数,因为它使用前后的值来构造三次样条.

from scipy.interpolate import InterpolatedUnivariateSpline

# Get a function that evaluates the linear spline at any x
f = InterpolatedUnivariateSpline(x, y, k=1)

# Get a function that evaluates the derivative of the linear spline at any x
dfdx = f.derivative()

# Evaluate the derivative dydx at each x location...
dydx = dfdx(x)

Answer 7

为了计算梯度，机器学习社区使用 Autograd：

"高效计算 numpy 代码的导数。 "

安装：

pip install autograd

下面是一个例子：

import autograd.numpy as np
from autograd import grad

def fct(x):
    y = x**2+1
    return y

grad_fct = grad(fct)
print(grad_fct(1.0))

它还可以计算复杂函数的梯度，例如多元函数。