Стандартная ошибка в нелинейной регрессии - PullRequest
7 голосов
/ 19 августа 2011

Я проводил некоторые физические симуляции Монте-Карло с Python и не смог определить стандартную ошибку для коэффициентов нелинейного приближения методом наименьших квадратов.

Изначально я использовал SciPy scipy.stats.linregress для своей модели, так как думал, что это будет линейная модель, но заметил, что на самом деле это какая-то степенная функция. Затем я использовал NumPy polyfit со степенями свободы, равными 2, но я никак не могу определить стандартную ошибку коэффициентов.

Я знаю, что gnuplot может определить ошибки для меня, но мне нужно выполнить подборки для более чем 30 различных случаев. Мне было интересно, знает ли кто-нибудь, что Python читает стандартную ошибку из gnuplot, или есть какая-то другая библиотека, которую я могу использовать?

Ответы [ 2 ]

10 голосов
/ 20 августа 2011

Наконец-то нашел ответ на этот давно заданный вопрос!Я надеюсь, что это может, по крайней мере, сэкономить кому-то несколько часов безнадежных исследований на эту тему.Scipy имеет специальную функцию кривой_fit в секции оптимизации.Он использует метод наименьших квадратов для определения коэффициентов и, что лучше всего, он дает ковариационную матрицу.Ковариационная матрица содержит дисперсию каждого коэффициента.Точнее, диагональ матрицы - это дисперсия, и путем квадратного укоренения значений можно определить стандартную ошибку каждого коэффициента!Scipy не располагает достаточным количеством документации для этого, поэтому вот пример кода для лучшего понимания:

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


def func(x,a,b,c):
    return a*x**2 + b*x + c #Refer [1]

x = np.linspace(0,4,50)
y = func(x,2.6,2,3) + 4*np.random.normal(size=len(x)) #Refer [2] 


coeff, var_matrix = curve_fit(func,x,y)
variance = np.diagonal(var_matrix) #Refer [3]

SE = np.sqrt(variance) #Refer [4]

#======Making a dictionary to print results========
results = {'a':[coeff[0],SE[0]],'b':[coeff[1],SE[1]],'c':[coeff[2],SE[2]]}

print "Coeff\tValue\t\tError"
for v,c in results.iteritems():
    print v,"\t",c[0],"\t",c[1]
#========End Results Printing=================

y2 = func(x,coeff[0],coeff[1],coeff[2]) #Saves the y values for the fitted model

plot.plot(x,y)
plot.plot(x,y2)

plot.show()
  1. То, что возвращает эта функция, очень важно, потому что оно определяет, что будет использоваться для подгонки к модели
  2. Использование функции для создания произвольных данных + некоторый шум
  3. Сохраняет диагональ ковариационной матрицы в 1D-матрице, которая является просто нормальным массивом
  4. Квадрат, корень отклонения, чтобы получить стандартошибка (SE)
2 голосов
/ 20 августа 2011

похоже, что gnuplot использует levenberg - marquardt и доступна реализация Python - вы можете получить оценки ошибок из атрибута mpfit.covar (случайноВы должны беспокоиться о том, что «ошибки» имеют в виду оценки ошибок - можно ли корректировать другие параметры, например, для компенсации?)

...