SciPyのLeastSqを使用して、実験的なスペクトルを理論的な期待値に適合させています。もちろん、実験値に関連するエラーがあります。これらをLeastSqにフィードするにはどうすればよいですか、それとも別のルーチンが必要ですか?ドキュメントには何も見つかりません。
2 に答える
scipy.optimize.leastsq関数には、重みを組み込むための組み込みの方法がありません。ただし、scipy.optimize.curve_fit関数にはsigma
、各yデータポイントの分散を示すために使用できるパラメーターがあります。
curve_fit
1.0/sigma
重みとして使用します。ここでsigma
、は長さの配列N
(と同じ長さydata
)にすることができます。
したがって、どういうわけか、エラーバーのサイズに基づいて各ydataポイントの分散を推測し、それを使用してを決定する必要がありますsigma
。
たとえば、エラーバーの長さの半分が1標準偏差を表すと宣言した場合、分散(いわゆる)curve_fit
はsigma
標準偏差の2乗になります。
sigma = (length_of_error_bar/2)**2
参照:
私は自分でこれを行っている最中なので、私が行ったことを共有し、おそらくコミュニティからいくつかのコメントを得ることができます。標準偏差を計算した一定の時間間隔で取得されたデータポイントのコレクションがあります。これらの点をsin関数に当てはめたいと思います。Leastsqは、残差、または一連のパラメーターpに基づくデータポイントと近似関数の差を最小化することによってこれを行います。残差を分散または標準偏差の2乗で割ることにより、残差に重みを付けることができます。
次のように:
from scipy.optimize import leastsq
import numpy as np
from matplotlib import pyplot as plt
def sin_func(t, p):
""" Returns the sin function for the parameters:
p[0] := amplitude
p[1] := period/wavelength
p[2] := phase offset
p[3] := amplitude offset
"""
y = p[0]*np.sin(2*np.pi/p[1]*t+p[2])+p[3]
return y
def sin_residuals(p, y, t, std):
err = (y - p[0]*np.sin(2*np.pi/p[1]*t+p[2])-p[3])/std**2
return err
def sin_fit(t, ydata, std, p0):
""" Fits a set of data, ydata, on a domain, t, with individual standard
deviations, std, to a sin curve given the initial parameters, p0, of the form:
p[0] := amplitude
p[1] := period/wavelength
p[2] := phase offset
p[3] := amplitude offset
"""
# optimization #
pbest = leastsq(sin_residuals, p0, args=(ydata, t, std), full_output=1)
p_fit = pbest[0]
# fit to data #
fit = p_fit[0]*np.sin(2*np.pi/p_fit[1]*t+p_fit[2])+p_fit[3]
return p_fit