カーブフィッティング手法 scipy.optimize.curve_fit の使い方を理解する


Pythonを使ってカーブフィッティング（曲線近似）する方法として、 scipy.optimize.curve_fit を使う方法がありますが、使い方が少し理解しにくいと思ったので整理してみました。

# 用いる実験値

<a href="https://qiita.com/maskot1977/items/f9e84fb6e62f7015d6c4">Numpy.polyfit を使ったカーブフィッティング</a> で用いたのと同じデータを用いて、比較してみましょう。


```python
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 型に変換しておきましょう。


```python
import numpy as np
x_observed = np.array(x_observed)
y_observed = np.array(y_observed)
```

# １次式近似

まず、１次式に近似してみましょう。関数 func1 を定義します。a と b が、求めたい値になります。


```python
def func1(X, a, b): # １次式近似
    Y = a + b * X
    return Y
```

関数 func1 の使用例はこちらです。


```python
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 を使って近似してみましょう。


```python
from scipy.optimize import curve_fit  
popt, pcov = curve_fit(func1,x_observed,y_observed) # poptは最適推定値、pcovは共分散
popt
```




    array([125.77023172,  -0.1605313 ])



ここで得られた popt が最適推定値を格納しています。<a href="https://qiita.com/maskot1977/items/f9e84fb6e62f7015d6c4">Numpy.polyfit を使ったカーブフィッティング</a> で得られた推定値と比較してみましょう。

# ２次式近似

次に、２次式に近似してみましょう。関数 func2 を定義します。a と b と c が、求めたい値になります。


```python
def func2(X, a, b, c): # 2次式近似
    Y = a + b * X + c * X ** 2
    return Y
```

関数 func2 の使用例はこちらです。


```python
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 を使って近似してみましょう。


```python
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 が最適推定値を格納しています。<a href="https://qiita.com/maskot1977/items/f9e84fb6e62f7015d6c4">Numpy.polyfit を使ったカーブフィッティング</a> で得られた推定値と比較してみましょう。

次のようにすれば、最適推定値を用いた近似曲線が描けます。


```python
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 が、求めたい値になります。


```python
def func3(X, a, b, c, d): # 3次式近似
    Y = a + b * X + c * X ** 2 + d * X ** 3
    return Y
```

さて、scipy.optimize.curve_fit を使って近似してみましょう。


```python
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 が最適推定値を格納しています。<a href="https://qiita.com/maskot1977/items/f9e84fb6e62f7015d6c4">Numpy.polyfit を使ったカーブフィッティング</a> で得られた推定値と比較してみましょう。

次のようにすれば、最適推定値を用いた近似曲線が描けます。


```python
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])



# 何次式になっても使える汎用的な関数を作りたい

今まで１次式、２次式、３次式と別々の関数を作りましたが、これだと面倒すぎるので、何次式になっても使える汎用的な関数を作りましょう。何次式に近似するかはパラメータの数で自動的に判定するようにします。


```python
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 による計算結果と同じであることを確認しましょう。


```python
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])




```python
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= で初期値を必要数だけ入力する必要が生じます。


```python
from scipy.optimize import curve_fit  
popt, pcov = curve_fit(func,x_observed,y_observed, p0=[1, 1]) 
popt
```




    array([125.77023172,  -0.1605313 ])




```python
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])




```python
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])




```python
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 を使ったカーブフィッティングとの比較

上記の計算結果を、<a href="https://qiita.com/maskot1977/items/f9e84fb6e62f7015d6c4">Numpy.polyfit を使ったカーブフィッティング</a> で得られた推定値と比較してみましょう。

Numpy.polyfit は、多項式近似するだけなら、便利で使いやすいですが、多項式近似しかできません。もっと他の関数で近似したい場合は、scipy.optimize.curve_fit の使い方を理解するのが良いでしょう。


# 少ない観測値を補間してから、正規分布の線形和で近似する

<a href="https://qiita.com/maskot1977/items/d16a385e23ed07e28037">少ない観測値を補間してから、正規分布の線形和で近似する</a> では、scipy.optimize.curve_fitを使って、正規分布（ガウス関数）の線形和で近似する方法について記述していますが、ここで説明した scipy.optimize.curve_fit の使い方を踏まえれば、理解が容易になるのではないかと期待します。