curve_fit 没有优化参数之一答案

【问题标题】：curve_fit not optimizing one of the parameterscurve_fit 没有优化参数之一
【发布时间】：2018-06-17 05:58:53
【问题描述】：

我需要拟合scipy.optimize.curve_fit 一些看起来像图中点的数据。我使用一个函数y(x)（参见下面的定义），它为x<x0 提供一个常数y(x)=c，否则为多项式（例如第二条倾斜线y1 = mx+q）。

我对参数(x0, c, m, q)给出了一个合理的初步猜测，如图所示。拟合结果表明，除第一个参数x0 外，所有参数均已优化。

为什么会这样？是我如何定义函数testfit(x, *p)，其中x0 (=p[0]) 出现在另一个函数中？

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

# generate some data:
x = np.linspace(0,100,1000)
y1 = np.repeat(0, 500)
y2 = x[500:] - 50
y = np.concatenate((y1,y2))
y = y + np.random.randn(len(y))


def testfit(x, *p):
    ''' piecewise function used to fit 
        it's a constant (=p[1]) for x < p[0]
        or a polynomial for x > p[0]     
    '''
    x = x.astype(float)
    y = np.piecewise(x, [x < p[0], x >= p[0]], [p[1], lambda x: np.poly1d(p[2:])(x)])
    return y

# initial guess, one horizontal and one tilted line:
p0_guess = (30, 5, 0.3, -10)

popt, pcov = curve_fit(testfit, x, y, p0=p0_guess)

print('params guessed  : '+str(p0_guess))
print('params from fit : '+str(popt))

plt.plot(x,y, '.')
plt.plot(x, testfit(x, *p0_guess), label='initial guess')
plt.plot(x, testfit(x, *popt), label='final fit')
plt.legend()

输出

params guessed  : (30, 5, 0.3, -10) 
params from fit : [ 30. 0.04970411   0.80106256 -34.17194401] 

OptimizeWarning: Covariance of the parameters could not be estimated category=OptimizeWarning)

【问题讨论】：

我相信您面临的问题本质上与this question 相似。您有两种选择：1. 使用不同的拟合算法。 2. 尽量避免曲线中的垂直边缘，方法是 a 给它一个小斜率，或者 b 让多项式与常数平滑连接值）。
正如@kazemakase 建议的那样，curve_fit 不会很好地处理像您的 p[0] 这样的离散变量，因为它无法为该变量创建偏导数。

标签： python scipy curve-fitting

【解决方案1】：

按照 kazemakase 的建议，我通过在我用来拟合的两个函数（一条水平线后跟一个多项式）之间的平滑过渡解决了这个问题。诀窍是将一个函数乘以sigmoid(x)，将另一个乘以1-sigmoid(x)（其中sigmoid(x) 定义如下）。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

x = np.linspace(0,100,1000)
y1 = np.repeat(0, 500)
y2 = x[500:] - 50
y = np.concatenate((y1,y2))
y = y + np.random.randn(len(y))

def testfit(x, *p):
    ''' function to fit the indentation curve 
    p = [x0,c, poly1d_coeffs ]'''
    x = x.astype(float)
    y = p[1]*(1-sigmoid(x-p[0],k=1)) + np.poly1d(p[2:])(x) * sigmoid(x-p[0],k=1)
    return y

def sigmoid(x, k=1):
    return 1/(1+np.exp(-k*x))

p0_guess = (30, 5, 0.3, -10 )
popt, pcov = curve_fit(testfit, x, y, p0=p0_guess)
print('params guessed  : '+str(p0_guess))
print('params from fit : '+str(popt))


plt.figure(1)
plt.clf()
plt.plot(x,y, 'y.')
plt.plot(x, testfit(x, *p0_guess), label='initial guess')
plt.plot(x, testfit(x, *popt), 'k', label='final fit')
plt.legend()

【讨论】：

【解决方案2】：

我遇到了类似的问题。我最终使用 np.gradient 和卷积来平滑曲线，然后绘制它。比如：

def mov_avg(n, data):
    return np.convolve(data, np.ones((n,))/n, mode='valid')

如果你想要更直接的方法，你可以试试这个：

def find_change(data):
    def test_flag(pos):
        grad = np.gradient(data) - np.gradient(data).mean()
        return (grad[:pos]<0).sum() + (grad[pos:]>0).sum()
    return np.vectorize(test_flag)(np.arange(len(data)-1)).argmax()

def find_gradient(pos, data):
    return np.gradient(data[:pos]).mean(), np.gradient(data[pos:]).mean()

pos=find_change(x2)
print(pos, find_gradient(pos, data))

第一个函数通过比较点梯度与平均梯度来计算梯度变化的点，并找到梯度“大部分为正”的点。

希望对你有帮助

【讨论】：