Существует раздражающее отсутствие качественных пакетов Python для интеграции SDE.У вас есть какая-то необычная установка, в которой ваша случайная переменная имеет явную зависимость от вашей подынтегральной функции (по крайней мере, кажется, что y зависит от p).Из-за этого, вероятно, будет трудно найти существующую реализацию, отвечающую вашим потребностям.
К счастью, самый простой способ интеграции SDE, Euler-Maruyama , очень прост в реализации., как в функции eulmar ниже:
from matplotlib import pyplot as plt
import numpy as np
def eulmar(func, randfunc, x0, tinit, tfinal, dt):
times = np.arange(tinit, tfinal + dt, dt)
x = np.zeros(times.size, dtype=float)
x[0] = x0
for i,t in enumerate(times[1:]):
x[i+1] = x[i] + func(x[i], t) + randfunc(x[i], t)
return times, x
Затем вы можете использовать eulmar для интеграции вашего SDE следующим образом:
def func(x, t):
return 1 - 2*x
def randfunc(x, t):
return np.random.rand()
times,x = eulmar(func, randfunc, 0, 0, 500, 5)
plt.plot(times, x)
Вы должны будете предоставить свои randfunc.Как и выше, это должна быть функция, которая принимает x и t в качестве аргументов и возвращает одну выборку из вашей случайной величины y.Если у вас возникли проблемы с поиском способа генерирования выборок y, поскольку вы знаете PDF, вы всегда можете использовать выборку отклонения (хотя она, как правило, довольно неэффективна).
Примечания
- Это не особенно эффективная реализация Эйлера-Маруямы.Например, случайные выборки обычно генерируются одновременно (например,
np.random.rand(500)).Однако вы не можете предварительно сгенерировать случайные выборки, поскольку y зависит от p.