Как вычислить производную, используя Numpy?

Как рассчитать производную функции, например

у = х2+1

с помощью numpy?

Допустим, я хочу значение производной при х = 5...

Источник

9 ответов

Решение

У вас есть четыре варианта

  1. Вы можете использовать конечные различия
  2. Вы можете использовать автоматические производные
  3. Вы можете использовать символическое дифференцирование
  4. Вы можете вычислить производные вручную.

Конечные различия не требуют внешних инструментов, но подвержены числовым ошибкам и, если вы находитесь в многомерной ситуации, могут занять некоторое время.

Символическая дифференциация идеальна, если ваша задача достаточно проста. Символические методы становятся достаточно надежными в наши дни. SymPy - отличный проект для этого, который хорошо интегрируется с NumPy. Посмотрите на функции autowrap или lambdify или посмотрите пост Дженсена в блоге о подобном вопросе.

Автоматические производные очень крутые, не подвержены числовым ошибкам, но требуют некоторых дополнительных библиотек (для этого есть несколько хороших вариантов). Это самый надежный, но и самый сложный / сложный в настройке выбор. Если ты в порядке, ограничиваясь 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.]
Источник

Самый простой способ, который я могу придумать, - это использовать функцию градиента 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) вектор размера.

Источник

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)
Источник

Предполагая, что вы хотите использовать numpyВы можете численно вычислить производную функции в любой точке, используя определение Rigorous:

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
Источник

Я брошу другой метод на кучу...

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)
Источник

Вы можете использовать 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
Источник

Для расчета градиентов сообщество машинного обучения использует 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))

Он также может вычислять градиенты сложных функций, например, многомерных функций.

Источник

В зависимости от уровня точности, который вам требуется, вы можете решить это самостоятельно, используя простое доказательство дифференциации:

>>> (((5 + 0.1) ** 2 + 1) - ((5) ** 2 + 1)) / 0.1
10.09999999999998
>>> (((5 + 0.01) ** 2 + 1) - ((5) ** 2 + 1)) / 0.01
10.009999999999764
>>> (((5 + 0.0000000001) ** 2 + 1) - ((5) ** 2 + 1)) / 0.0000000001
10.00000082740371

на самом деле мы не можем взять предел градиента, но это довольно забавно. Вы должны остерегаться, хотя, потому что

>>> (((5+0.0000000000000001)**2+1)-((5)**2+1))/0.0000000000000001
0.0
Источник

Чтобы вычислить производную числовой функции, используйте эту схему конечных разностей второго порядка, как показано на:https://youtu.be/5QnToSn_oxk?t=1804 .

      dx = 0.01
x = np.arange(-4, 4+dx, dx)
y = np.sin(x)
n = np.size(x)

yp = np.zeros(n)
yp[0] = (-3*y[0] + 4*y[1] - y[2]) / (2*dx)
yp[n-1] = (3 * y[n-1] - 4*y[n-2] + y[n-3]) / (2*dx)
for j in range(1,n-1):
    yp[j] = (y[j+1] - y[j-1]) / (2*dx)

Или, если вы хотите использовать более высокий порядок, используйте:https://youtu.be/5QnToSn_oxk?t=1374

Все это из лекций Натана Куца по курсу «Начало научных вычислений».

Источник
Другие вопросы по тегам