用积分函数拟合数据

Question

使用Python中的curve_fit函数从scipy.optimize库拟合数据时，首先需要定义拟合函数（例如二次多项式）如下：

1. def f(x, a, b): return a*x**2+b*x
2. 然后进行拟合 popt, pcov = curve_fit(f,x,y)

但问题是，如果函数包含积分或离散求和怎么办，例如：

实验数据仍然给出x和f(x)，因此我想一旦能够在Python中定义f(x)，点2应该是类似的。顺便说一下，这里假设g(t)具有已知形式，并且包含拟合参数，即像多项式示例中给出的a和b等参数。非常感谢您的帮助。问题实际上应该是通用的，本文中使用的函数只是随机示例。

- user4587874

---

@cfh 哦，我明白了，这是真的，但如果它没有任何封闭形式的解决方案，那么数值积分究竟包含什么？它不是假设所有参数都应该知道吗？ - user4587874

@cfh 当然可以，但在那个例子中，我不需要进行数值积分，因此在拟合之前我不需要知道a和b。但是，根据您的建议“首先评估积分”，我必须在拟合之前知道a和b（以便进行数值积分），而我不知道... - user4587874

不。当 curve_fit 函数调用你的 f 时，它将始终为 a 和 b 提供特定值。您可以使用这些值来评估多项式、计算积分或进行任何其他操作。 - cfh

@cfh 哦，现在我明白你的意思了，没错！你介意（或者如果有时间的话）给个简单的例子吗？比如一个数值积分相关的例子。 - user4587874

这真的很大程度上取决于您尝试集成的函数。一些用于执行此操作的函数包含在scipy.integrate中。 - cfh

Answer 1

这里是一个拟合基于积分定义的曲线的例子。该曲线是从0到Pi上对于t的sin(t*w)/t+p积分。我们的x数据点对应于w，并且我们调整p参数以使得数据拟合。

import math, numpy, scipy.optimize, scipy.integrate

def integrand(t, args):
    w, p = args
    return math.sin(t * w)/t + p

def curve(w, p):
    res = scipy.integrate.quad(integrand, 0.0, math.pi, [w, p])
    return res[0]

vcurve = numpy.vectorize(curve, excluded=set([1]))

truexdata = numpy.asarray([0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0])
trueydata = vcurve(truexdata, 1.0)

xdata = truexdata + 0.1 * numpy.random.randn(8)
ydata = trueydata + 0.1 * numpy.random.randn(8)

popt, pcov = scipy.optimize.curve_fit(vcurve,
                                      xdata, ydata,
                                      p0=[2.0])
print popt

这将会输出一个非常接近1.0的结果，这也是我们在创建“trueydata”时所使用的“p”的值。

请注意，我们在曲线函数上使用“numpy.vectorize”来生成一个向量化版本，以便与“scipy.optimize.curve_fit”兼容。

- Jay Kominek

---

这很有趣、有用，而且我认为它足够通用，但是我们能不能做得更好，知道这是积分的一部分呢？在你的解决方案中，“curve”可以是任何东西，并且并不受到它是积分的事实的好处。也许问题应该改为“使用非向量化函数拟合数据”？ - dashesy

Answer 2

有时你可能很幸运，能够对积分进行解析评估。在以下示例中，将从0到无穷大积分 h(t)=exp(-(t-x)**2/2) 和一个二次多项式g(t) 的乘积。使用 Sympy 来评估积分并生成可用于 curve_fit() 的函数：

import sympy as sy
sy.init_printing()  # LaTeX-like pretty printing of IPython


t, x = sy.symbols("t, x", real=True)

h = sy.exp(-(t-x)**2/2)

a0, a1, a2 = sy.symbols('a:3', real=True)  # unknown coefficients
g = a0 + a1*t + a2*t**2

gh = (g*h).simplify()  # the intgrand
G = sy.integrate(gh, (t, 0, sy.oo)).simplify()  # integrate from 0 to infinty

# Generate numeric function to be usable by curve_fit()
G_opt = sy.lambdify((x, t, a0, a1, a2), G)

print(G_opt(1, 2, 3, 4, 5))  # example usage

需要注意的是，通常该问题往往不够明确，因为积分在解的足够大的邻域内不一定收敛（curve_fit() 假设解在该邻域内）。

- Dietrich