使用python修复拐点估计

问题描述

我正在尝试使用python查找曲线上的拐点.曲线的数据在此处: https://www.dropbox. com/s/rig8frgewde8i5n/fitted.txt?dl = 0 .请注意，曲线已拟合到原始数据.原始数据可在此处获取: https://www.dropbox.com/s/1lskykdi1ia1lu7/ww.txt?dl = 0

I am trying to find the inflexion points on a curve using python. The data for the curve is here: https://www.dropbox.com/s/rig8frgewde8i5n/fitted.txt?dl=0. Please note that the curve has been fitted to the raw data. Raw data is available here: https://www.dropbox.com/s/1lskykdi1ia1lu7/ww.txt?dl=0

import numpy as np
# Read in array from text file
arr = np.loadtxt(path_to_file)

inflexion_point_1 = np.diff(arr).argmin()
inflexion_point_2 = np.diff(arr).argmax()

这两个拐点在附图中显示为红线.但是，它们的位置似乎不正确.第一个拐点应靠近黑色箭头指示的区域.我该如何解决?

These 2 inflexion points are shows as red lines in the attached plot. However, their locations do not seem to be right. The first inflexion point should be close to the area indicated by the black arrow. How do I fix this?

此外，这是一个差异图:

Also, here is a plot of the differential:

plt.axvline(np.gradient(arr[:365]).argmax())

如您所见，代码的行为与编码相同，即找到数组的np.diff的argmax.但是，我想找到一个接近第110天的位置，即距argmax大约一半的位置.

As you can see, the code is behaving as coded i.e. it finds the argmax of np.diff of the array. However, I want to find a position closer to day 110 or so, i.e. about half-way to argmax.

-编辑-

这是另一个显示原始数据和拟合曲线(使用二次函数)的图.

Also, here is another plot showing thee raw data and the fitted curve (using a quadratic function).

推荐答案

是否有理由不直接在渐变上使用单变量样条线?

Any reason not to use uni-variate spline directly on the gradient?

from scipy.interpolate import UnivariateSpline

#raw data
data = np.genfromtxt('ww.txt')

plt.plot(np.gradient(data), '+')

spl = UnivariateSpline(np.arange(len(data)), np.gradient(data), k=5)
spl.set_smoothing_factor(1000)
plt.plot(spl(np.arange(len(data))), label='Smooth Fct 1e3')
spl.set_smoothing_factor(10000)
plt.plot(spl(np.arange(len(data))), label='Smooth Fct 1e4')
plt.legend(loc='lower left')

max_idx = np.argmax(spl(np.arange(len(data))))
plt.vlines(max_idx, -5, 9, linewidth=5, alpha=0.3)

我们也可以数值求解最大值:

Also we can solve for the maximum numerically:

In [122]:

import scipy.optimize as so
F = lambda x: -spl(x)
so.fmin(F, 102)
Optimization terminated successfully.
         Current function value: -3.339112
         Iterations: 20
         Function evaluations: 40
Out[122]:
array([ 124.91303558])

这篇关于使用python修复拐点估计的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！