scipy curve_fit无法适合tophat函数 [英] scipy curve_fit cannot fit a tophat function
问题描述
我正在尝试将高礼帽功能适合某些数据，即f（x）对于整个实线是恒定的，除了一个有限长度的段等于另一个常数。我的参数是tophat函数的两个常量，中点和宽度，我正在尝试使用scipy.optimize.curve_fit来获取所有这些常量。不幸的是，curve_fit无法获得帽子的宽度。无论我做什么，它都会拒绝测试除我开始使用的那个宽度以外的任何宽度值，并且会非常糟糕地拟合其余数据。以下代码段说明了该问题：
I am trying to fit a top hat function to some data, ie. f(x) is constant for the entire real line, except for one segment of finite length which is equal to another constant. My parameters are the two constants of the tophat function, the midpoint, and the width and I'm trying to use scipy.optimize.curve_fit to get all of these. Unfortunately, curve_fit is having trouble obtaining the width of the hat. No matter what I do, it refuses to test any value of the width other than the one I start with, and fits the rest of the data very badly. The following code snippet illustrates the problem:
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
def tophat(x, base_level, hat_level, hat_mid, hat_width):
ret=[]
for xx in x:
if hat_mid-hat_width/2. < xx < hat_mid+hat_width/2.:
ret.append(hat_level)
else:
ret.append(base_level)
return np.array(ret)
x = np.arange(-10., 10., 0.01)
y = tophat(x, 1.0, 5.0, 0.0, 1.0)+np.random.rand(len(x))*0.2-0.1
guesses = [ [1.0, 5.0, 0.0, 1.0],
[1.0, 5.0, 0.0, 0.1],
[1.0, 5.0, 0.0, 2.0] ]
plt.plot(x,y)
for guess in guesses:
popt, pcov = curve_fit( tophat, x, y, p0=guess )
print popt
plt.plot( x, tophat(x, popt[0], popt[1], popt[2], popt[3]) )
plt.show()
为什么curve_fit如此正确如此可怕，我该如何解决呢？
Why is curve_fit so extremely terrible at getting this right, and how can I fix it?
推荐答案
首先， tophat 的定义可以使用 numpy.where 而不是循环：
First, the definition of tophat could use numpy.where instead of a loop:
def tophat(x, base_level, hat_level, hat_mid, hat_width):
return np.where((hat_mid-hat_width/2. < x) & (x < hat_mid+hat_width/2.), hat_level, base_level)
第二，棘手的不连续目标函数抵制 curve_fit 调用的优化算法。 Nelder-Mead方法通常对于粗糙函数更可取，但是 curve_fit 似乎无法使用它。因此，我建立了一个目标函数（只是偏差的绝对值之和），并将其最小化：
Second, the tricky discontinuous objective function resists the optimization algorithms that curve_fit calls. The Nelder-Mead method is usually preferable for rough functions, but it looks like curve_fit cannot use it. So I set up an objective function (just the sum of absolute values of deviations) and minimize that:
def objective(params, x, y):
return np.sum(np.abs(tophat(x, *params) - y))
plt.plot(x,y)
for guess in guesses:
res = minimize(objective, guess, args=(x, y), method='Nelder-Mead')
print(res.x)
plt.plot(x, tophat(x, *(res.x)))
结果是更好的做法是，从宽度为2的宽帽子开始，使其缩小到正确的大小（请参阅三个猜测的最后一个）。
The results are better, in that starting with a too-wide hat of width 2 makes it shrink down to the correct size (see the last of three guesses).
[9.96041297e-01 5.00035502e+00 2.39462103e-04 9.99759984e-01]
[ 1.00115808e+00 4.94088711e+00 -2.21340843e-05 1.04924153e-01]
[9.95947108e-01 4.99871040e+00 1.26575116e-03 9.97908018e-01]
不幸的是，何时最初的猜测是太窄的帽子，优化器仍然卡住了。
Unfortunately, when the starting guess is a too-narrow hat, the optimizer is still stuck.
您可以尝试其他优化方法/目标函数组合但我还没有找到能使帽子可靠展开的帽子。
You can try other optimization method / objective function combinations but I haven't found one that makes the hat reliably expand.
要尝试的一件事是不要使用太接近真实水平的参数。这有时可能会受伤。
One thing to try is not to use the parameters that are too close to the true levels; this sometimes might hurt. With
guesses = [ [1.0, 1.0, 0.0, 1.0],
[1.0, 1.0, 0.0, 0.1],
[1.0, 1.0, 0.0, 2.0] ]
我曾经设法
[ 1.00131181 4.99156649 -0.01109271 0.96822019]
[ 1.00137925 4.97879423 -0.05091561 1.096166 ]
[ 1.00130568 4.98679988 -0.01133717 0.99339777]
这三个宽度都正确。但是，这只是几次尝试中的一些（优化过程的初始化中存在一些随机性）。其他一些具有相同初始点的尝试均失败了；该过程不够鲁棒。
which is correct for all three widths. However, this was only on some of several tries (there is some randomness in the initialization of the optimizing procedure). Some other attempts with the same initial points failed; the process is not robust enough.
这篇关于scipy curve_fit无法适合tophat函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
