パラメータ数と誤差
モデルのパラメータの数を増やしていくと誤差が増大していくとBrunton & Kutz (2019) に説明されている。これは、データに対してモデルの適合度を上げるためにパラメータ数を増やす、すなわちモデルを複雑にしてデータに対する適合度を上げるという考え方に反することである。Brunton & Kutz (2019)で扱われている誤差は、真のモデルからの差であり、モデルを複雑にしてデータとの適合度を上げるときの誤差は、データとの差であることに注意する必要がある。このことを確認する簡単なスクリプトをBrunton & Kutz (2019)のCode 4.14を基にして書いてみた（リスト１）。
真のモデルは2次関数
である。
データは、真のモデルにガウスノイズを加えて作成されている。
データに対して、p次の多項式モデルを適用して適用モデルと誤差との関係が調べられた。
モデルの係数は最小2乗法により推定された。すなわち、
と置いたとき、パラメータの推定値は
で与えられている（cf. 岡本（2009）、p. 26）。
Brunton & Kutz (2019)で論じられている誤差は、適用モデルの真値からの差であり、次式で与えられる。
これに対して、データと適用モデルとの差は次式で与えられる。
２つの誤差
と
の比較は、リスト１を実行して得られる図１に示される通りである。
図１
真値からの差としての誤差はBrunton & Kutz (2019)のFigure 4.15と同じ傾向が確認できる。モデルが複雑になる（多項式モデルの次数が上がる）と真値からの誤差は増加する。しかし、データとの差としての誤差は、パラメータ数が増えてモデルの複雑度が高くなると減少している。
パラメータ数を増やして、モデルが複雑なデータの挙動に対応できるようにすると、モデルは真値から乖離していくことを見るために、リスト２のスクリプトを用意した。シミュレーションの状況設定は上のリスト１の場合と同じであり、真のモデルは２次の多項式である。適用モデルのグラフが多いと重なり合って見難くなるので、適用モデルとして２次、３次、15次の多項式を取り上げ、それぞれ５回のシミュレーションを行い、グラフを描いた。リスト２のスクリプトを実行すると、図２のグラフが表示される。
図２
左側のグラフの左下の部分を示したものが、右側のグラフである。赤のグラフが、真値のグラフであり、データ値のグラフは黒の点線で描かれている。データ値（点線）は真値の周りにかなり変動しているが、2次の適用モデルは真値の近くで安定している。これに対して、高次のモデルである15次の多項式モデルは、データに合わせて変動（オーバーフィッティング）している。モデルの次数が上がるとデータに合わせた変動が容易になるので、オーバーフィッティングのために真値からは遠ざかっていくことが分かる。単純にモデルを複雑化してデータに対する適合度を求めるのではなく、単純さを重視すること（オッカムのレーザー）も留意されるべきであろう。オーバーフィッティングは、一般化の妨げになる虞れがある。
import numpy
as np
import matplotlib.pyplot
as plt
M = 20
n = 1000
L = 4
n_pnts = 100
x = np.linspace(0,
L, n_pnts)
En = np.empty((n,
M))
Err = np.empty((n,
M))
for jj in
range(M):
A = np.empty((n_pnts, jj+1))
for j in range(jj+1):
A[:,j] = x ** j
f = x**2
# True value
for j in range(n):
# True value +
Noise
fn = f + 0.1 * np.random.normal(size=n_pnts)
b = np.linalg.pinv(A) @ fn
fna = A @ b
# Error from the
true values
En[j, jj] = np.linalg.norm(f
- fna) / np.linalg.norm(f)
# Error fron the given data
Err[j, jj] = np.linalg.norm(fn - fna) / np.linalg.norm(fn)
plt.figure(figsize=(12,5))
plt.subplot(1,2,1)
Mean_En = np.mean(En, axis=0)
plt.plot(range(M), Mean_En,
'^', c='r', lw=0, label='from True values')
Mean_Err = np.mean(Err,
axis=0)
plt.plot(range(M), Mean_Err,
'v', c = 'g', lw=0, label='from Data')
plt.boxplot(En,
positions=np.arange(0,M))
plt.boxplot(Err, positions=np.arange(0,M))
plt.xticks(np.arange(0,
M))
plt.xlabel('Polynomial Degree',fontsize=14)
plt.ylabel('Error', fontsize=14)
plt.title('Degree >= 0', fontsize=16)
plt.legend(fontsize=14)
plt.subplot(1,2,2)
plt.plot(np.arange(M)[2:],
Mean_En[2:], '^-', c='r', lw=1,
label='from True values')
plt.plot(np.arange(M)[2:],
Mean_Err[2:], 'v-', c = 'g', lw=1,
label='from Data')
plt.xticks(np.arange(2,
M))
plt.xlabel('Polynomial Degree', fontsize=14)
plt.ylabel('Error', fontsize=14)
plt.boxplot(En[:,2:],
positions=np.arange(2, M))
plt.boxplot(Err[:,2:], positions=np.arange(2, M))
plt.title('Degree >= 2', fontsize=16)
plt.legend()
plt.tight_layout()
plt.show()
import numpy
as np
import matplotlib.pyplot
as plt
M = [2, 3, 15]
# パラメータ数
n = 5
# シミュレーションの回数
n_pnts = 100
# サンプルサイズ
x = np.linspace(0,
4, n_pnts)
# 変域
f = x**2
# 真のモデル
fn = np.empty((n,
n_pnts)) # （シミュレーション、サンプリング点）
fna = np.empty((3,
n, n_pnts)) # （パラメータ数条件、シミュレーション、サンプリング点）
def c_v(jj):
if jj
== 2:
return 'g'
elif
jj == 3:
return 'b'
else:
return 'm'
for j in range(n):
# 真値 ＋ ノイズ
fn[j]
= f + 0.1 * np.random.normal(size=n_pnts)
for i,
jj in enumerate(M):
A = np.empty((n_pnts, jj+1))
for jv in range(jj+1):
A[:,jv] = x ** jv
b = np.linalg.pinv(A) @ fn[j]
fna[i, j] = A @ b
# モデルによる予測値
plt.figure(figsize=(12,5))
plt.subplot(1,2,1)
plt.plot(x, f, c='r')
for j in range(n):
plt.plot(x,
fn[j], c='k', ls=':', lw=2)
for i,
jj in enumerate(M):
plt.plot(x, fna[i, j], c=c_v(jj),
ls='-', alpha=0.5)
plt.xlabel('x', fontsize=16)
plt.xticks(fontsize=14)
plt.ylabel('f(x)', fontsize=16)
plt.yticks(fontsize=14,
rotation=45, va='center')
plt.plot([], c='r', label='True Value')
plt.plot([], c='k', ls=':', lw=2, label='Data')
for jj in M:
plt.plot([],
c=c_v(jj), ls='-',
alpha=0.5, label=f'p={jj}')
plt.legend(fontsize=16)
plt.subplot(1,2,2)
plt.plot(x[:20], f[:20], c='r')
for j in range(n):
plt.plot(x[:20],
fn[j][:20], c='k', ls=':', lw=2)
for i,
jj in enumerate(M):
plt.plot(x[:20], fna[i, j, :20], c=c_v(jj), ls='-', lw=1)
plt.xlabel('x', fontsize=16)
plt.xticks(fontsize=14)
plt.ylabel('f(x)', fontsize=16)
plt.yticks(fontsize=14,
rotation=45, va='center')
plt.plot([], c='r', label='True Value')
plt.plot([], c='k', ls=':', lw=2, label='Data')
for jj in M:
plt.plot([],
c=c_v(jj), ls='-', lw=1, label=f'p={jj}')
plt.legend(fontsize=14)
plt.tight_layout()
plt.show()
S. L. Brunton and
J. N. Kutz (2019) Data-driven
science and engineering: Machine learning, dynamical systems, and control.
Cambridge University Press.
岡本安晴（2009）統計学を学ぶための数学入門［下］．培風館
岡本安晴（2019）いまさら聞けないPythonでデータ分析．丸善出版