Pythonを使ってカーブフィッティング（曲線近似）する方法として、 scipy.optimize.curve_fit を使う方法がありますが、使い方が少し理解しにくいと思ったので整理してみました。
用いる実験値
Numpy.polyfit を使ったカーブフィッティング で用いたのと同じデータを用いて、比較してみましょう。
x_observed = [9, 28, 38, 58, 88, 98, 108, 118, 128, 138, 148, 158, 168, 178, 188, 198, 208, 218, 228, 238, 278, 288, 298]
y_observed = [51, 80, 112, 294, 286, 110, 59, 70, 56, 70, 104, 59, 59, 72, 87, 99, 64, 60, 74, 151, 157, 57, 83]
Python の List 型から、Numpy の Array 型に変換しておきましょう。
import numpy as np
x_observed = np.array(x_observed)
y_observed = np.array(y_observed)
１次式近似
まず、１次式に近似してみましょう。関数 func1 を定義します。a と b が、求めたい値になります。
def func1(X, a, b): # １次式近似
Y = a + b * X
return Y
関数 func1 の使用例はこちらです。
func1(x_observed, 2, 2)
array([ 20, 58, 78, 118, 178, 198, 218, 238, 258, 278, 298, 318, 338,
358, 378, 398, 418, 438, 458, 478, 558, 578, 598])
さて、scipy.optimize.curve_fit を使って近似してみましょう。
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func1,x_observed,y_observed) # poptは最適推定値、pcovは共分散
popt
array([125.77023172, -0.1605313 ])
ここで得られた popt が最適推定値を格納しています。Numpy.polyfit を使ったカーブフィッティング で得られた推定値と比較してみましょう。
２次式近似
次に、２次式に近似してみましょう。関数 func2 を定義します。a と b と c が、求めたい値になります。
def func2(X, a, b, c): # 2次式近似
Y = a + b * X + c * X ** 2
return Y
関数 func2 の使用例はこちらです。
func2(x_observed, 2, 2, 2)
array([ 182, 1626, 2966, 6846, 15666, 19406, 23546, 28086,
33026, 38366, 44106, 50246, 56786, 63726, 71066, 78806,
86946, 95486, 104426, 113766, 155126, 166466, 178206])
さて、scipy.optimize.curve_fit を使って近似してみましょう。
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func2,x_observed,y_observed) # poptは最適推定値、pcovは共分散
popt
array([ 1.39787684e+02, -4.07809516e-01, 7.97183226e-04])
ここで得られた popt が最適推定値を格納しています。Numpy.polyfit を使ったカーブフィッティング で得られた推定値と比較してみましょう。
次のようにすれば、最適推定値を用いた近似曲線が描けます。
func2(x_observed, 1.39787684e+02, -4.07809516e-01, 7.97183226e-04)
array([136.1819702 , 128.9940092 , 125.44205497, 118.81645644,
110.07383349, 107.47849913, 105.04260142, 102.76614035,
100.64911593, 98.69152815, 96.89337701, 95.25466253,
93.77538468, 92.45554348, 91.29513893, 90.29417102,
89.45263976, 88.77054514, 88.24788717, 87.88466585,
88.02614699, 88.46010889, 89.05350743])
３次式近似
同様に、３次式に近似してみましょう。関数 func3 を定義します。a と b と c と d が、求めたい値になります。
def func3(X, a, b, c, d): # 3次式近似
Y = a + b * X + c * X ** 2 + d * X ** 3
return Y
さて、scipy.optimize.curve_fit を使って近似してみましょう。
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func3,x_observed,y_observed) # poptは最適推定値、pcovは共分散
popt
array([ 7.84214107e+01, 1.88213104e+00, -1.74165777e-02, 3.89638087e-05])
ここで得られた popt が最適推定値を格納しています。Numpy.polyfit を使ったカーブフィッティング で得られた推定値と比較してみましょう。
次のようにすれば、最適推定値を用いた近似曲線が描けます。
func3(x_observed, 7.84214107e+01, 1.88213104e+00, -1.74165777e-02, 3.89638087e-05)
array([ 93.97825188, 118.32181643, 126.93087413, 136.59795028,
135.72770915, 132.27386543, 127.62777811, 122.02323006,
115.69400413, 108.87388316, 101.79665001, 94.69608754,
87.80597859, 81.36010602, 75.59225267, 70.73620141,
67.02573508, 64.69463654, 63.97668864, 65.10567422,
92.76660851, 106.63700433, 123.75703075])
何次式になっても使える汎用的な関数を作りたい
今まで１次式、２次式、３次式と別々の関数を作りましたが、これだと面倒すぎるので、何次式になっても使える汎用的な関数を作りましょう。何次式に近似するかはパラメータの数で自動的に判定するようにします。
import numpy as np
def func(X, *params):
Y = np.zeros_like(X)
for i, param in enumerate(params):
Y = Y + np.array(param * X ** i)
return Y
次の２つのコードを実行して、上記の func2 による計算結果と、func3 による計算結果と同じであることを確認しましょう。
func(x_observed, 1.39787684e+02, -4.07809516e-01, 7.97183226e-04)
array([136.1819702 , 128.9940092 , 125.44205497, 118.81645644,
110.07383349, 107.47849913, 105.04260142, 102.76614035,
100.64911593, 98.69152815, 96.89337701, 95.25466253,
93.77538468, 92.45554348, 91.29513893, 90.29417102,
89.45263976, 88.77054514, 88.24788717, 87.88466585,
88.02614699, 88.46010889, 89.05350743])
func(x_observed, 7.84214107e+01, 1.88213104e+00, -1.74165777e-02, 3.89638087e-05)
array([ 93.97825188, 118.32181643, 126.93087413, 136.59795028,
135.72770915, 132.27386543, 127.62777811, 122.02323006,
115.69400413, 108.87388316, 101.79665001, 94.69608754,
87.80597859, 81.36010602, 75.59225267, 70.73620141,
67.02573508, 64.69463654, 63.97668864, 65.10567422,
92.76660851, 106.63700433, 123.75703075])
何次式になっても使える汎用的な関数による多項式近似
ここまでくれば、何次式にでも多項式近似できます。ただし、パラメータの数を指定するために p0= で初期値を必要数だけ入力する必要が生じます。
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func,x_observed,y_observed, p0=[1, 1])
popt
array([125.77023172, -0.1605313 ])
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func,x_observed,y_observed, p0=[1, 1, 1])
popt
array([ 1.39787684e+02, -4.07809516e-01, 7.97183226e-04])
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func,x_observed,y_observed, p0=[1, 1, 1, 1])
popt
array([ 7.84214107e+01, 1.88213104e+00, -1.74165777e-02, 3.89638087e-05])
from scipy.optimize import curve_fit
popt, pcov = curve_fit(func,x_observed,y_observed, p0=[1, 1, 1, 1, 1])
popt
array([-7.54495613e+01, 1.13179580e+01, -1.50591743e-01, 7.02970546e-04,
-1.07313869e-06])
Numpy.polyfit を使ったカーブフィッティングとの比較
上記の計算結果を、Numpy.polyfit を使ったカーブフィッティング で得られた推定値と比較してみましょう。
Numpy.polyfit は、多項式近似するだけなら、便利で使いやすいですが、多項式近似しかできません。もっと他の関数で近似したい場合は、scipy.optimize.curve_fit の使い方を理解するのが良いでしょう。
少ない観測値を補間してから、正規分布の線形和で近似する
少ない観測値を補間してから、正規分布の線形和で近似する では、scipy.optimize.curve_fitを使って、正規分布（ガウス関数）の線形和で近似する方法について記述していますが、ここで説明した scipy.optimize.curve_fit の使い方を踏まえれば、理解が容易になるのではないかと期待します。