Python NumPyのlegfit()でルジャンドル級数をデータに最小二乗フィットする方法

データに対してルジャンドル級数の最小二乗フィットを求めるには、Python NumPy の legendre.legfit() メソッドを使用します。このメソッドは、低次から高次の順に並べられたルジャンドル係数を返します。また、y が 2 次元配列の場合、y の k 列目のデータに対応する係数は、結果の k 列目に格納されます。

主なパラメータ

x： M 個のサンプル点 (x[i], y[i]) の x 座標を指定します。

y： サンプル点の y 座標です。同じ x 座標を共有する複数組のサンプル点を、1 回の呼び出しで独立にフィットさせることも可能です。その場合は、各列に 1 つのデータセットを含む 2 次元配列を y として渡します。

deg： フィッティング多項式の次数です。deg に単一の整数を指定すると、deg 次までのすべての項がフィットに含まれます。

rcond： フィットの相対条件数です。最大特異値と比較して rcond より小さい特異値は無視されます。デフォルト値は len(x)*eps で、eps はプラットフォームの float 型の相対精度（多くの場合 約 2e-16）です。

full： 戻り値の内容を決めるスイッチです。False（デフォルト）の場合は係数のみを返し、True の場合は特異値分解による診断情報も一緒に返します。

w： 重みです。None 以外の値を指定した場合、重み w[i] は x[i] における残差 y[i] − ŷ[i] に適用されます。理想的には、w[i]*y[i] の誤差がすべて同じ分散になるよう重みを選択します。逆分散重み付けを行う場合は w[i] = 1/sigma(y[i]) を指定します。デフォルト値は None です。

手順

まず、必要なライブラリをインポートします。

import numpy as np
from numpy.polynomial import legendre as L

x 座標を作成します。

x = np.linspace(-1,1,51)

x 座標を表示します。

print("X Co-ordinate...\n",x)

y 座標を作成して表示します。

y = x**3 - x + np.random.randn(len(x))
print("\nY Co-ordinate...\n",y)

ルジャンドル級数の最小二乗フィットを求めるために、legendre.legfit() メソッドを使用します。full=True を指定することで、診断情報も取得できます。

c, stats = L.legfit(x,y,3,full=True)
print("\nResult...\n",c)
print("\nResult...\n",stats)

完全なコード例

import numpy as np
from numpy.polynomial import legendre as L

# x座標の作成
x = np.linspace(-1,1,51)

# x座標の表示
print("X Co-ordinate...\n",x)

# y座標の作成
y = x**3 - x + np.random.randn(len(x))
print("\nY Co-ordinate...\n",y)

# legfit()メソッドで最小二乗フィットを実行
c, stats = L.legfit(x,y,3,full=True)
print("\nResult...\n",c)
print("\nResult...\n",stats)

実行結果

X Co-ordinate...
  [-1.   -0.96 -0.92 -0.88 -0.84 -0.8  -0.76 -0.72 -0.68 -0.64 -0.6  -0.56
   -0.52 -0.48 -0.44 -0.4  -0.36 -0.32 -0.28 -0.24 -0.2  -0.16 -0.12 -0.08
   -0.04  0.    0.04  0.08  0.12  0.16  0.2   0.24  0.28  0.32  0.36  0.4
    0.44  0.48  0.52  0.56  0.6   0.64  0.68  0.72  0.76  0.8   0.84  0.88
    0.92  0.96  1. ]

   Y Co-ordinate...
  [-5.28795520e-02 -7.61252904e-03  7.35194215e-02 -1.33072588e-01
   -1.21785636e+00  7.75679385e-02  6.55168668e-01  1.42872448e+00
    8.42326214e-01  2.49667989e+00  9.58942508e-01 -2.67332869e-01
   -7.85575928e-01  1.93333045e+00  7.32492468e-01  5.23576961e-01
   -1.91529521e+00 -1.41434385e+00  4.44787373e-01  3.81831261e-01
    3.74128321e-01  1.20562789e+00  1.44870029e+00  1.01091575e-03
    8.94334713e-01  1.22342199e+00  9.52055370e-01 -7.29520012e-01
   -2.42648820e-01 -9.78434555e-02  1.27468237e-01  9.39489448e-01
    1.08795136e+00  2.31230197e+00  1.93107556e-02 -6.13335407e-01
    1.93170835e-01 -8.77958854e-01 -3.59868085e-01  4.31331759e-01
    7.24929856e-01 -2.22736540e-01 -1.29623093e+00  4.13226024e-01
    7.82155644e-01 -1.56618537e-01  1.25043737e+00  6.32386988e-01
   -2.75716271e-01  8.80669895e-02 -3.20225560e-01]

Result...
 [ 0.29249467 -0.10521942 -0.24847572 0.2010877 ]

Result...
 [array([39.35467561]), 4, array([1.0425003 , 1.02126704, 0.97827074, 0.95561139]), 1.1324274851176597e-14]