python - SciPy + Numpy : Finding the slope of a sigmoid curve

Question

我有一些数据遵循 sigmoid 分布，如下图所示:

对数据进行标准化和缩放后，我使用 scipy.optimize.curve_fit 和一些初始参数调整了底部的曲线:

popt, pcov = curve_fit(sigmoid_function, xdata, ydata, p0 = [0.05, 0.05, 0.05])
>>> print popt
[  2.82019932e+02  -1.90996563e-01   5.00000000e-02]

所以popt，根据the documentation ，返回*“参数的最佳值，使 f(xdata, popt) - ydata 的平方误差之和最小化”。我在这里理解，没有用curve_fit计算斜率，因为我不认为这条平缓曲线的斜率是282，都不是负数。

然后我尝试使用scipy.optimize.leastsq，因为文档说它返回“解决方案(或失败调用的最后一次迭代的结果)。 ”，所以我认为坡度会被归还。像这样:

p, cov, infodict, mesg, ier = leastsq(residuals, p_guess, args = (nxdata, nydata), full_output=True)
>>> print p
Param(x0=281.73193626250207, y0=-0.012731420027056234, c=1.0069006606656596, k=0.18836680131910222)

但是，我再次没有得到我所期望的。 curve_fit 和 leastsq 返回几乎相同的值，我想这并不奇怪，因为 curve_fit 正在使用最小二乘法的实现找到曲线。但不会向后倾斜......除非我忽略了一些东西。

那么，如何计算 X = 285 且 Y = 0.5 的点的斜率？

我试图避免手动方法，例如 calculating the derivative比如说，(285.5, 0.55) 和 (284.5, 0.45) 并减去和除以结果等等。我想知道是否有更自动的方法。

谢谢大家!

编辑#1

这是我的“sigmoid_function”，由 curve_fit 和 lesssq 方法使用:

def sigmoid_function(xdata, x0, k, p0): # p0 not used anymore, only its components (x0, k)
    # This function is called by two different methods: curve_fit and leastsq,
    # this last one through function "residuals". I don't know if it makes sense
    # to use a single function for two (somewhat similar) methods, but there 
    # it goes.

    # p0:
    #   + Is the initial parameter for scipy.optimize.curve_fit. 
    #   + For residuals calculation is left empty
    #   + It is initialized to [0.05, 0.05, 0.05]
    # x0:
    #   + Is the convergence parameter in X-axis and also the shift
    #   + It starts with 0.05 and ends up being around ~282 (days in a year)
    # k:
    #   + Set up either by curve_fit or leastsq
    #   + In least squares it is initially fixed at 0.5 and in curve_fit
    #   + to 0.05. Why? Just did this approach in two different ways and 
    #   + it seems it is working. 
    #   + But honestly, I have no clue on what it represents
    # xdata: 
    #   + Positions in X-axis. In this case from 240 to 365

# Finally I changed those parameters as suggested in the answer. 
# Sigmoid curve has 2 degrees of freedom, therefore, the initial 
# guess only needs to be this size. In this case, p0 = [282, 0.5]


    y = np.exp(-k*(xdata-x0)) / (1 + np.exp(-k*(xdata-x0)))
    return y

def residuals(p_guess, xdata, ydata):
    # For the residuals calculation, there is no need of setting up the initial parameters
    # After fixing the initial guess and sigmoid_function header, remove [] 
    # return ydata - sigmoid_function(xdata, p_guess[0], p_guess[1], [])
    return ydata - sigmoid_function(xdata, p_guess[0], p_guess[1], [])

如果我在描述参数时犯了错误或混淆了技术术语，我很抱歉。我对 numpy 很陌生，而且我已经很多年没有学习数学了，所以我正在重新学习。

那么，对于计算此数据集的 X = 285、Y = 0.5(大致为中点)的斜率，您有何建议？谢谢!!

编辑#2

感谢 Oliver W.，我按照他的建议更新了我的代码，并且更好地理解了问题。

最后一个细节我还没有完全明白。显然，curve_fit 返回一个 popt 数组 (x0, k)，其中包含拟合的最佳参数:

• x0 似乎是通过指示曲线的中心点来确定曲线的移动程度
• k 参数是 y = 0.5 时的斜率，也在曲线的中心(我认为!)

为什么如果 sigmoid 函数是一个增长函数，popt 中的导数/斜率是负数？有道理吗？

我使用 sigmoid_derivative 来计算斜率，是的，我获得了与 popt 相同的结果，但带有正号。

# Year 2003, 2005, 2007. Slope in midpoint.
k = [-0.1910, -0.2545, -0.2259] # Values coming from popt
slope = [0.1910, 0.2545, 0.2259] # Values coming from sigmoid_derivative function

我知道这有点夸张，因为我可以同时使用两者。相关数据在那里，但带有负号，但我想知道为什么会发生这种情况。

因此，只有当我需要知道 y = 0.5 以外的其他点的斜率时，才需要按照您的建议计算导数函数。仅对于中点，我可以使用 popt。

感谢您的帮助，它节省了我很多时间。 :-)

Answer 1

您从未使用过传递给 sigmoid 函数的参数 p0。因此，曲线拟合不会有任何好的方法来找到收敛性，因为它可以为此参数取任何值。您应该首先像这样重写 sigmoid 函数:

def sigmoid_function(xdata, x0, k):

    y = np.exp(-k*(xdata-x0)) / (1 + np.exp(-k*(xdata-x0)))
    return y

这意味着您的模型(S形)只有两个自由度。这将在 popt 中返回:

initial_guess = [282, 1]  # (x0, k): at x0, the sigmoid reaches 50%, k is slope related
popt, pcov = curve_fit(sigmoid_function, xdata, ydata, p0=initial_guess)

现在，popt 将是一个元组(或 2 个值的数组)，即最好的 x0 和 k。

要获得该函数在任意点的斜率，老实说，我只会象征性地计算导数，因为 sigmoid 并不是一个很难的函数。你最终会得到:

def sigmoid_derivative(x, x0, k):
    f = np.exp(-k*(x-x0))
    return -k / f

如果您将曲线拟合的结果存储在 popt 中，则可以轻松地将其传递给此函数:

print(sigmoid_derivative(285, *popt))

这将为您返回 x=285 处的导数。但是，因为您专门要求中点，所以当 x==x0 和 y==.5 时，您会看到(从 sigmoid_derivative)那里的导数只是 -k，可以从您已经获得的 curve_fit 输出中立即观察到。在您显示的输出中，该值约为 0.19。