pythonでフィッティングをする

pythonでfittingをする方法。

例えば、

というをパラメータとする関数でデータ点を

が最小になるようにfittingしたいとする（最小二乗法）。
scipy.optimizeのcurve_fitを使うのが楽（scipy.optimizeにはleastsqという関数もあり、こちらでも同じことができるが、curve_fitの方が分かりやすい）。

import numpy as np
import scipy.optimize
import matplotlib.pylab as plt

# data which you want to fit
xdata = np.array([0.0,1.0,2.0,3.0,4.0,5.0,6.0])
ydata = np.array([0.1,0.9,2.2,2.8,4.2,5.9,7.4])
# initial guess for the parameters
parameter_initial = np.array([0.0, 0.0, 0.0]) #a, b, c
# function to fit
def func(x, a, b, c):
    return a + b*x + c*x*x

paramater_optimal, covariance = scipy.optimize.curve_fit(func, xdata, ydata, p0=parameter_initial)
print "paramater =", paramater_optimal

y = func(xdata,paramater_optimal[0],paramater_optimal[1],paramater_optimal[2])
plt.plot(xdata, ydata, 'o')
plt.plot(xdata, y, '-')
plt.show()

xdata, ydataはfittingしたいデータ点（上の例では手で作ってるが、普通は実験などで得られる）。
parameter_initialはパラメータの初期値（なくてもいい）。
funcがfittingする関数で、第一引数がxで、第二引数以降がパラメータ。返り値はy。
paramater_optimalがfittingで最適化されたパラメータで、covarianceは共分散。

結果は

paramater = [ 0.14761905  0.70357143  0.08452381]

もし、scipy.optimize.leastsqを使うなら、以下のようにすればいい。

xdata = np.array([0.0,1.0,2.0,3.0,4.0,5.0,6.0])
ydata = np.array([0.1,0.9,2.2,2.8,4.2,5.9,7.4])
parameter_initial = np.array([0.0, 0.0, 0.0])
def func(x, a, b, c):
    return a + b*x + c*x*x
def residual_func(parameter, x, y):
    residual = y - func(x, parameter[0], parameter[1], parameter[2])
    return residual
result = scipy.optimize.leastsq(residual_func, parameter_initial, args=(xdata,ydata))
paramater_optimal = result[0]
y = func(xdata,paramater_optimal[0],paramater_optimal[1],paramater_optimal[2])
plt.plot(xdata, ydata, 'o')
plt.plot(xdata, y, '-')
plt.show()

leastsqの場合は初期値は絶対必要。

Ref:
http://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html
http://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.leastsq.html
http://python4mpia.github.io/fitting_data/least-squares-fitting.html