scipy curve_fit повышает "OptimizeWarning: ковариация параметров не может быть оценена" - PullRequest
0 голосов
/ 16 мая 2018

Я пытаюсь приспособить эту функцию к некоторым данным:

enter image description here

Но когда я использую свой код

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

def f(x, start, end):
    res = np.empty_like(x)
    res[x < start] =-1
    res[x > end] = 1
    linear = np.all([[start <= x], [x <= end]], axis=0)[0]
    res[linear] = np.linspace(-1., 1., num=np.sum(linear))
    return res

if __name__ == '__main__':

    xdata = np.linspace(0., 1000., 1000)
    ydata = -np.ones(1000)
    ydata[500:1000] = 1.
    ydata = ydata + np.random.normal(0., 0.25, len(ydata))

    popt, pcov = curve_fit(f, xdata, ydata, p0=[495., 505.])
    print(popt, pcov)
    plt.figure()
    plt.plot(xdata, f(xdata, *popt), 'r-', label='fit')
    plt.plot(xdata, ydata, 'b-', label='data')
    plt.show()

Я получаю сообщение об ошибке

OptimizeWarning: Ковариация параметров не может быть оценена

Вывод:

enter image description here

В этом примере начало и конец должны быть ближе к 500, но они не меняются по сравнению с моим первоначальным предположением.

1 Ответ

0 голосов
/ 16 мая 2018

Предупреждение (не ошибка)

OptimizeWarning: Covariance of the parameters could not be estimated

означает, что подгонка не может определить неопределенности (дисперсию) параметров подгонки.

Основная проблема заключается в том, что функция модели f обрабатывает параметры start и end как дискретные значения - они используются как целочисленные местоположения для изменения функциональной формы.curve_fit от scipy (и все другие процедуры оптимизации в scipy.optimize) предполагают, что параметры являются непрерывными переменными, а не дискретными.

Процедура подбора попытается сделать небольшие шаги (обычно в точности станка) по параметрам, чтобы получить числовую производную от остатка по переменным (якобиан).При значениях, используемых в качестве дискретных переменных, эти производные будут равны нулю, и процедура подбора не будет знать, как изменить значения для улучшения подгонки.

Похоже, вы пытаетесь подогнать пошаговую функцию к некоторым данным,Позвольте мне порекомендовать попробовать lmfit (https://lmfit.github.io/lmfit-py), который обеспечивает высокоуровневый интерфейс для подгонки кривой и имеет много встроенных моделей. Например, он включает StepModel, который должен иметь возможность моделировать вашdata.

Для небольшого изменения ваших данных (так, чтобы он имел конечный шаг), следующий скрипт с lmfit может соответствовать таким данным:

#!/usr/bin/python
import numpy as np
from lmfit.models import StepModel, LinearModel
import matplotlib.pyplot as plt

np.random.seed(0)
xdata = np.linspace(0., 1000., 1000)
ydata = -np.ones(1000)
ydata[500:1000] = 1.
# note that a linear step is added here:
ydata[490:510] = -1 + np.arange(20)/10.0
ydata = ydata + np.random.normal(size=len(xdata), scale=0.1)

# model data as Step + Line
step_mod = StepModel(form='linear', prefix='step_')
line_mod = LinearModel(prefix='line_')

model = step_mod + line_mod

# make named parameters, giving initial values:
pars = model.make_params(line_intercept=ydata.min(),
                         line_slope=0,
                         step_center=xdata.mean(),
                         step_amplitude=ydata.std(),
                         step_sigma=2.0)

# fit data to this model with these parameters
out = model.fit(ydata, pars, x=xdata)

# print results
print(out.fit_report())

# plot data and best-fit
plt.plot(xdata, ydata, 'b')
plt.plot(xdata, out.best_fit, 'r-')
plt.show()

, который печатаетотчет

[[Model]]
    (Model(step, prefix='step_', form='linear') + Model(linear, prefix='line_'))
[[Fit Statistics]]
    # fitting method   = leastsq
    # function evals   = 49
    # data points      = 1000
    # variables        = 5
    chi-square         = 9.72660131
    reduced chi-square = 0.00977548
    Akaike info crit   = -4622.89074
    Bayesian info crit = -4598.35197
[[Variables]]
    step_sigma:      20.6227793 +/- 0.77214167 (3.74%) (init = 2)
    step_center:     490.167878 +/- 0.44804412 (0.09%) (init = 500)
    step_amplitude:  1.98946656 +/- 0.01304854 (0.66%) (init = 0.996283)
    line_intercept: -1.00628058 +/- 0.00706005 (0.70%) (init = -1.277259)
    line_slope:      1.3947e-05 +/- 2.2340e-05 (160.18%) (init = 0)
[[Correlations]] (unreported correlations are < 0.100)
    C(step_amplitude, line_slope)     = -0.875
    C(step_sigma, step_center)        = -0.863
    C(line_intercept, line_slope)     = -0.774
    C(step_amplitude, line_intercept) =  0.461
    C(step_sigma, step_amplitude)     =  0.170
    C(step_sigma, line_slope)         = -0.147
    C(step_center, step_amplitude)    = -0.146
    C(step_center, line_slope)        =  0.127

и создает график enter image description here

Lmfit имеет множество дополнительных функций. Например, если вы хотите установить границы для некоторых изЗначения параметров или некоторые из них могут быть разными: вы можете сделать следующее:

# make named parameters, giving initial values:
pars = model.make_params(line_intercept=ydata.min(),
                         line_slope=0,
                         step_center=xdata.mean(),
                         step_amplitude=ydata.std(),
                         step_sigma=2.0)

# now set max and min values for step amplitude"
pars['step_amplitude'].min = 0
pars['step_amplitude'].max = 100

# fix the offset of the line to be -1.0
pars['line_offset'].value = -1.0
pars['line_offset'].vary = False

# then run fit with these parameters
out = model.fit(ydata, pars, x=xdata)

Если вы знаете, что модель должна быть Step+Constant, а постоянная должна быть фиксированной, вы также можете изменить модель так, чтобы

from lmfit.models import ConstantModel
# model data as Step + Constant
step_mod = StepModel(form='linear', prefix='step_')
const_mod = ConstantModel(prefix='const_')

model = step_mod + const_mod

pars = model.make_params(const_c=-1,
                         step_center=xdata.mean(),
                         step_amplitude=ydata.std(),
                         step_sigma=2.0)
pars['const_c'].vary = False
...