Как указать ошибку Y-переменной при подгонке с lmfit? - PullRequest
0 голосов
/ 23 февраля 2019

Я почти новичок в Python и пытаюсь собрать данные из колледжа, используя lmfit.Переменная Y имеет переменную ошибку 3%.Как добавить эту ошибку в процесс подбора?Я изменяю от подгонки кривой scipy, и в scipy это было действительно легко сделать, просто создав массив со значениями ошибок и указав ошибку при подгонке, добавив текст "sigma = [yourarray]". Это мой текущий код:

from lmfit import Minimizer, Parameters, report_fit
import matplotlib.pyplot as plt
w1, V1, phi1, scal1 = np.loadtxt("./FiltroPasaBajo_1.txt", delimiter = "\t", unpack = True)

t = w1
eV= V1*0.03 + 0.01

def funcion(parametros, x, y):
    R = parametros['R'].value
    C = parametros['C'].value

    modelo = 4/((1+(x**2)*(R**2)*(C**2))**1/2)
    return modelo - y
parametros = Parameters()
parametros.add('R', value = 1000, min= 900, max = 1100)
parametros.add('C', value = 1E-6, min = 1E-7, max = 1E-5)

fit = Minimizer(funcion, parametros, fcn_args=(t,V1))
resultado = fit.minimize()


final = V1 + resultado.residual
report_fit(resultado)

try:
    plt.plot(t, V1, 'k+')
    plt.plot(t, final, 'r')
    plt.show()
except ImportError:
    pass


V1 - значения, которые я измерил, а eV - массив ошибок.t является координатой х.Спасибо за ваше время

Ответы [ 2 ]

0 голосов
/ 24 февраля 2019

Функция minimize() минимизирует массив в смысле наименьших квадратов, настраивая параметры переменной, чтобы минимизировать (resid**2).sum() для массива resid, возвращаемого вашей целевой функцией.Он действительно ничего не знает о неопределенности в ваших данных или даже о ваших данных.Чтобы использовать неопределенности в вашей аппроксимации, вы должны передать свой массив eV точно так же, как вы передаете t и V1, а затем использовать это при вычислении минимизируемого массива.

Обычно требуется минимизировать Sum[ (data-model)^2/epsilon^2 ], где epsilon - это неопределенность данных (ваш eV), поэтому остаточный массив следует изменить с data-model на (data-model)/epsilon.Для вашей подгонки, вы хотели бы

def funcion(parametros, x, y, eps):
    R = parametros['R'].value
    C = parametros['C'].value

    modelo = 4/((1+(x**2)*(R**2)*(C**2))**1/2)
    return (modelo - y)/eps

, а затем использовать это с

fit = Minimizer(funcion, parametros, fcn_args=(t, V1, eV))
resultado = fit.minimize()
...

Если вы используете интерфейс lmfit.Model (разработан для подгонки кривой), то вы можете передатьв массиве weights, который умножается на data -model, и поэтому будет 1.0 / eV для представления весов для неопределенностей (как указано выше с minimize).Использование интерфейса lmfit.Model и предоставление неопределенностей будет выглядеть следующим образом:

from lmfit import Model
# model function, to model the data
def func(t, r, c):
    return 4/((1+(t**2)*(r**2)*(c**2))**1/2)

model  = Model(func)
parametros = model.make_params(r=1000, c=1.e-6)

parametros['r'].set(min=900, max=1100)
parametros['c'].set(min=1.e-7, max=1.e-5)

resultado = model.fit(V1, parametros, t=t, weights=1.0/eV)

print(resultado.fit_report())

plt.errorbar(t, V1, eV, 'k+', label='data')
plt.plot(t, resultado.best_fit, 'r', label='fit')
plt.legend()
plt.show()

надеюсь, что это поможет ....

0 голосов
/ 23 февраля 2019

Я думаю, что вы не можете предоставить сигму в fit.minimize () напрямую.

Однако я вижу, что fit.minimize () использует метод scipy leastsq (по умолчанию), который является тем же методом, который используется scipy's curve_fit.

Если вы посмотрите на Scipy's curve_fit source , он следит за сигмой (для 1-го случая).

 transform = 1.0 / sigma
 jac = _wrap_jac(jac, xdata, transform)
 res = leastsq(func, p0, Dfun=jac, full_output=1, **kwargs)

Поскольку fit.minimize () позволяет вам передавать kwargs (Dfun) для leastsq, вы можете передать jac так, как это делается в scipy curve_fit.

...