如何根据 Python 中数据集的最佳拟合确定微分方程的未知参数?

问题描述

我正在尝试使用 python 将不同的微分方程拟合到给定的数据集.为此，我分别使用了 scipy 包和 solve_ivp 函数.只要我对微分方程中包含的参数(b = 0.005)有一个粗略的估计，这对我来说很好用，例如:

I am trying to fit different differential equations to a given data set with python. For this reason, I use the scipy package, respectively the solve_ivp function.This works fine for me, as long as I have a rough estimate of the parameters (b= 0.005) included in the differential equations, e.g:

import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
import numpy as np


def f(x, y, b):
    dydx= [-b[0] * y[0]]
    return dydx

xspan= np.linspace(1, 500, 25)
yinit= [5]
b= [0.005]


sol= solve_ivp(lambda x, y: f(x, y, b),
               [xspan[0], xspan[-1]], yinit, t_eval= xspan)

print(sol)
print("\n")

print(sol.t)
print(sol.y)

plt.plot(sol.t, sol.y[0], "b--")

然而，我想要实现的是，参数 b(或更多参数)是自动"确定的.基于求解的微分方程对给定数据集(x 和 y)的最佳拟合.有没有办法做到这一点，例如将这个例子与 scipy 的 curve_fit 函数结合起来，这看起来如何?

However, what I like to achieve is, that the parameter b (or more parameters) is/are determined "automatically" based on the best fit of the solved differential equation to a given data set (x and y). Is there a way this can be done, for example by combining this example with the curve_fit function of scipy and how would this look?

先谢谢你！

推荐答案

是的，你的想法应该可行，应该很容易插在一起.你想打电话

Yes, what you think about should work, it should be easy to plug together. You want to call

popt, pcov = scipy.optimize.curve_fit(curve, xdata, ydata, p0=[b0])
b = popt[0]

您现在必须定义一个函数 curve(x,*p)，该函数根据唯一的参数 b 将任何点列表转换为值列表.

where you now have to define a function curve(x,*p) that transforms any list of point into a list of values according to the only parameter b.

def curve(x,b):
    res = solve_ivp(odefun, [1,500], [5], t_eval=x, args = [b])
    return res.y[0]

根据需要为容错添加可选参数.

Add optional arguments for error tolerances as necessary.

为了使这更现实，还将初始点设为参数.然后，期望列表的位置和单个参数的位置也变得更加明显.为了获得适当的拟合任务，向测试数据添加一些随机噪声.也让归零不要那么快，这样最终的剧情还是有点意思的.

To make this more realistic, make also the initial point a parameter. Then it also becomes more obvious where a list is expected and where single arguments. To get a proper fitting task add some random noise to the test data. Also make the fall to zero not so fast, so that the final plot still looks somewhat interesting.

from scipy.integrate import solve_ivp
from scipy.optimize import curve_fit

xmin,xmax = 1,500

def f(t, y, b):
    dydt= -b * y
    return dydt

def curve(t, b, y0):
    sol= solve_ivp(lambda t, y: f(t, y, b),
               [xmin, xmax], [y0], t_eval= t)
    return sol.y[0]


xdata = np.linspace(xmin, xmax, 25)
ydata = np.exp(-0.02*xdata)+0.02*np.random.randn(*xdata.shape)

y0 = 5
b= 0.005
p0 = [b,y0]

popt, pcov = curve_fit(curve, xdata, ydata, p0=p0)
b, y0 = popt

print(f"b={b}, y0 = {y0}")

返回

b=0.019975693539459473, y0 = 0.9757709108115179

现在根据拟合曲线绘制测试数据

Now plot the test data against the fitted curve

这篇关于如何根据 Python 中数据集的最佳拟合确定微分方程的未知参数?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！