在给定曲线中查找相关点

Question

考虑以下 Python 代码，该代码绘制曲线并对其进行分析以找到一些点：

%matplotlib inline
import numpy as np
from numpy.polynomial.polynomial import Polynomial
from scipy.interpolate import UnivariateSpline
from scipy.signal import savgol_filter
import scipy.stats
import scipy.optimize
import matplotlib.pyplot as plt

# Create curve X axis.
x = np.arange(10000)
y = np.zeros_like(x, dtype=float)

# Create a straight line for the first 1000 points.
y[:1000] = (10.0,) * 1000
# The remaining points, create an exponential.
y[1000:] = 20 * np.exp(-x[1000:]/1200.0) 
# Add noise.
y += np.random.normal(0, 0.5, x.size)

# Create polynomial.
poly = Polynomial.fit(x, y, 15)
# Create Y values for the polynomial.
py = poly(x)

# Calculate first and second derivatives.
dxdy = np.gradient(py)
dx2dy2 = np.gradient(dxdy)

# Plot original curve and fit on top.
plt.plot(x, y, '', color='black')
plt.plot(x, py, '', color='red')

# Plot first order and second order derivatives.
plt.plot(x, dxdy * 100, '', color='blue')
plt.plot(x, dx2dy2, '', color='green')

我已经在图片上方标记了我要计算的点：

我想我可以用一阶导数来计算第一个，也就是过渡的开始。

我不太确定如何计算第二个，即过渡结束或曲线何时变平。我尝试使用曲线中最后 100 个点的平均值进行计算，并找到曲线中低于该平均值的第一个值，但是，它似乎不太可靠。

编辑1：

在研究一阶导数时，我想出了以下可能的解决方案，即找到峰值左右两侧的符号变化，我用一阶导数和信号的图像进行说明，并在下面拟合：

Answer 1

我更喜欢处理过滤（平滑）而不是插值数据。

我找到的第一点：

• 找到平滑数据的最大值
• 找到第一个值为最大值 90% 的点
• 返回寻找导数 >= 0 的第一个点

我找到的第二点

• 找到平滑数据的最小值
• 找到第一个值minimum + (maximum - minimum) * 1%的点

import numpy as np
import matplotlib.pyplot as plt
import scipy.ndimage

# Create curve X axis.
x = np.arange(10000)
y = np.zeros_like(x, dtype=float)

# Create a straight line for the first 1000 points.
y[:1000] = (10.0,) * 1000
# The remaining points, create an exponential.
y[1000:] = 20 * np.exp(-x[1000:]/1200.0)
# Add noise.
y += np.random.normal(0, 0.5, x.size)

y_smoothed = scipy.ndimage.uniform_filter1d(y, 250, mode='reflect')
diff_y_smoothed = np.diff(y_smoothed)
minimum = np.min(y_smoothed)
maximum = np.max(y_smoothed)
almost_first_point_idx = np.where(y_smoothed < 0.9 * maximum)[0][0]
first_point_idx = almost_first_point_idx - np.where(diff_y_smoothed[:almost_first_point_idx][::-1] >= 0)[0][0]
second_point_idx = np.where(y_smoothed < 0.01 * (maximum - minimum) + minimum)[0][0]

# Plot original curve and fit on top.
plt.plot(x, y, '', color='black')
plt.plot(x, y_smoothed, '', color='red')
plt.axvline(x[first_point_idx])
plt.axvline(x[second_point_idx])
plt.show()

您可能不会将这些参数调整为运行窗口大小（硬编码 250）和 1% 用于查找第二点。