演習：Ridge回帰では，ハイパーパラメータは0.0と0.1の2通りで固定してRidge回帰を行なった． この演習では，交差検証（cross-validation: CV）によりハイパーパラメータチューニングを行う．
下記のようにライブラリのインポートが必要
#%matplotlib inline # jupyter notebookでは必要．JupyterLabでは不要
import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import PolynomialFeatures, StandardScaler
from sklearn.linear_model import Ridge
from sklearn.metrics import root_mean_squared_error
from sklearn.model_selection import train_test_split
from sklearn.model_selection import KFold
演習1：k-fold CVとLOO
演習：Ridge回帰では訓練データとテストデータを分けて準備したが， 普通はデータセットが一つあって，それを訓練データとテストデータに分割して回帰等を行う． データの分け方は色々あるがここでは以下のサイトの上から2つ目の図のようなデータ分割を考える．
https://scikit-learn.org/stable/modules/cross_validation.html
K-分割交差検証（k-fold CV）とleave one out（LOO）がどういったものか調べよ．
また，下記サイトのExample（3.1.2.1.1. K-fol）を試してsklearn.model_selection.KFoldの使い方を学習せよ．
https://scikit-learn.org/stable/modules/cross_validation.html#k-fold
インデックスが返ってくると思うので，どのようにインデックスを利用すれば良いか理解すること． 下記のようにnumpyのarrayにした方がわかりやすい．
hogex = np.array(["a", "b", "c", "d"])
hogey = np.array([1, 2, 3, 4])
kf = KFold(n_splits=2)
i = 1
for train_indx, test_indx in kf.split(hogex):
print('loop: ', i)
print('train: ', hogex[train_indx], hogey[train_indx])
print('test: ', hogex[test_indx], hogey[test_indx])
i += 1
演習2：データの分割
まず演習：Ridge回帰で行なった時と同様にデータをいくつか準備する．
データ数はなんでも良いが，すぐに変えられるようにここではndataという変数を用いている．
データ数は自分で何度か変更してみて調節すること．
ここではデータを20点生成し，訓練データ16点とテストデータ4点に分ける．
ndata = 20
xa = ???
yt = ???
x = ??? # データの数はndata
tt = ???
err = ??? # データの数はndata
t = tt + err # sin関数にノイズを加える
plt.plot(xa, yt) # 元になるsin関数をプロット
plt.plot(x, t, 'o') # データ点をプロット
sklearn.model_selection.train_test_splitを用いて，このデータを訓練データとテストデータに分割せよ．
ここではテストデータの割合は20%程度とする（10%とすることが多いかもしれない）．
まずは下記ページのExamplesを試してみると良い．
sklearn.model_selection.train_test_split
sklearn.model_selection.train_test_split(*arrays, test_size=None, train_size=None, random_state=None, shuffle=True, stratify=None)
元のデータ数を確認．
x.shape
データを分割
x_train, x_test, t_train, t_test = ???
分割後のデータ数を確認．
print(x_train.shape)
print(x_test.shape)
print(t_train.shape)
print(t_test.shape)
ここで分割したテストデータは最後の方に利用する．
演習3：CVによるハイパーパラメータチューニング
k分割交差検証を用いて，演習：Ridge回帰で行なった時と同様にRidge回帰を行い， ハイパーパラメータチューニングを実施せよ．
これまで通り，計画行列を計算する関数を準備しておく．
def design_matrix(x, m):
pf = PolynomialFeatures(degree=m)
dmat = pf.fit_transform(x.reshape(-1, 1))
return dmat
ハイパーパラメータ alphaはnp.logspaceを使って色々変化させるとよい．
n_alphas = 100
alphas = np.logspace(-10, 0, n_alphas)
alphas # jupyterで表示用
上記alphasでループさせながら，k分割交差検証でRMSEをk回計算して，その平均を求める．
分割数も各自で調節すること．
以下で出てくるテストとか_testは上で分けておいたテストデータではなく，
交差検証の検証用データのことなので注意．
M = 9 # 多項式の次数
nsp = 4 # CVでの分割数
rmse_train = [] # 各alphaの値における訓練データに対するRMSEの平均
rmse_test = [] # 各alphaの値におけるテストデータに対するRMSEの平均
# ---------- loop for alpha
for a in alphas:
kf = KFold(n_splits=nsp)
krmse_train = [] # CV中の訓練データに対するRMSEを格納するためのリスト
krmse_test = [] # CV中のテストデータに対するRMSEを格納するためのリスト
# ---------- k-loop
for ktrain_indx, ktest_indx in kf.split(x_train):
# ------ k分割の訓練データと検証データ
x_ktrain, x_ktest = ???, ??? # x_trainを分割
t_ktrain, t_ktest = ???, ??? # t_trainを分割
# ------ 計画行列の作成と標準化
phi = design_matrix(???, M) # 訓練の計画行列作成
phi = ??? # 標準化前に定数項を外す
ss = StandardScaler()
phi_std = ??? # 標準化
# -- テストデータの計画行列と標準化
phi_test = design_matrix(???, M) # テストデータの計画行列作成
phi_test = ??? # 定数項を外す
phi_test_std = ??? # テストデータはfitはしない。transformだけ
# ------ フィット
ridge = Ridge(alpha=a, fit_intercept=True)
ridge.fit(???, ???) # fit
pre_train = ridge.predict(???) # 訓練データに対する予測値
pre_test = ridge.predict(???) # テストデータに対する予測値
# ------ RMSEの計算
krmse_train.append(root_mean_squared_error(???)) # 訓練
krmse_test.append(root_mean_squared_error(???)) # テスト
# ---------- k-fold CVの平均値の計算
krmse_train = np.array(krmse_train) # list --> array
krmse_test = np.array(krmse_test) # list --> array
rmse_train.append(krmse_train.mean()) # 各alphaに対してCVの平均値を格納
rmse_test.append(krmse_test.mean())
訓練とテスト（CV中の）データに関して，横軸をalpha，縦軸をRMSEにしたグラフをプロットする．
plt.plot(alphas, rmse_train, label='Train')
plt.plot(alphas, rmse_test, label='Test')
plt.xscale('log') # 横軸はlogスケールにする
plt.legend()
plt.xlabel('alpha')
plt.ylabel('RMSE')
plt.title('CV error')
CV中のテスト（検証）データに対するRMSEが最小となるalphaの値を求める．
a_opt = alphas[np.argmin(rmse_test)]
a_opt
演習4：最適化したハイパーパラメータを用いたRidge回帰
ハイパーパラメータチューニングができたので，求めたalphaの値と訓練データを全て用いて改めてRidge回帰を行う． Ridge回帰で得られたモデルと最初にデータ分割しておいたテストデータを比較して，作成したモデルの汎化性能を検証せよ．
phi = ??? # 訓練データに対する計画行列作成
phi = ??? # 定数項を外す
phi_test = ??? # テストデータの計画行列作成
phi_test = ??? # 定数項を外す
ss = StandardScaler()
phi_std = ss.??? # 標準化
phi_test_std = ss.??? # テストデータの標準化。fitしない
ridge = Ridge(alpha=a_opt, fit_intercept=True)
ridge.fit(???, ???) # 訓練データをfit
y_train = ridge.predict(???) # 訓練データに対する予測値
y_test = ridge.predict(???) # テストデータに対する予測値
比較のため，最小二乗法（$\alpha = 0$）の予測モデルも作成．ridge0という変数名にしておく．
ridge0 = Ridge(alpha=0, fit_intercept=True)
ridge0.fit(???, ???) # fit
予測モデルができたら，グラフをプロットして確認しておく． データ点ではなく，0から1の連続した（ように見える）xaに対する予測値を計算する．
phi_xa = ??? # xaの計画行列
phi_xa = ??? # 定数項を外す
phi_xa_std = ss.??? # fitはしない！同じスケーリングパラメーターを使う
ya = ridge.predict(???) # xaに対する予測値
y0 = ridge0.predict(???) # alpha=0の場合の予測値
plt.plot(xa, yt) # 元になるsin関数をプロット
plt.plot(x_train, t_train, 'o') # 訓練データをプロット
plt.plot(x_test, t_test, 'v') # テストデータをプロット
plt.plot(xa, ya) # 予測モデルをプロット
plt.ylim(-1.5, 1.5)
訓練データおよびテストデータ両方に対して，横軸が真値，縦軸が予測値のグラフを作成する．
plt.plot(t_train, y_train, 'o', label="training")
plt.plot(t_test, y_test, 'v', label="test")
plt.plot([-1.5, 1.5], [-1.5, 1.5]) # y = xの基準線
plt.ylim(-1.5, 1.5)
plt.xlabel('True')
plt.ylabel('Predicted')
plt.legend()
plt.title('Ridge')
訓練データ，テストデータに対する予測モデルとの誤差をRMSEを用いて計算する．
print('training: RMSE = ', mean_squared_error(???, ???, ???))
print('test: RMSE = ', mean_squared_error(???, ???, ???))