'20/12/13更新:PyCaretのverup2.2に伴いコードや図を編集、差し替えしました。
本記事では、機械学習でお馴染みのボストン住宅価格のデータセットを用いて、多数ある回帰モデルからベスト回帰モデルを抽出してオブジェクトファイルに保存する雛形コードを載せました。データ分析する例題のcsvファイルは、次のリンク先を例に作成できます。Python scikit-learn付属のボストン市の住宅価格データ(Boston house prices dataset)をcsvファイル化する - PythonとVBAで世の中を便利にする
本コードを実行すると、はじめに下図のように複数の回帰モデルを実行して各々5つの指標の結果を出力します。そして、指定した指標で並び替えによるランキング付けを行い、トップ5の回帰モデルをオブジェクトに格納します。

次に、そのトップ5を元に下図のようにたった一行で各々の回帰モデルのハイパーパラメータを自動調整します。

調整されたハイパーパラメータを下図のようにテキストファイルに保存して確認できるようにもしています。参考までに、回帰モデルの実使用はオブジェクト(.pklファイル)を扱うため、このテキストファイルを直接操作することはありません。

最終的に、下図のようにautoml()でベストな回帰モデルを自動で選定します。そして、predict_model()で回帰モデルの確認評価をすると、今回のベストな回帰モデルは、CatBoost Regressorということがわかります。

そして、interpret_model(best)と僅か一行のコードで、下図のようにシャープレイ値を可視化します。これにより、目的変数に対する感度の高い説明変数を視覚的に確認することが出来ます。横軸は目的変数に対する影響度(シャープレイ値)、縦軸は各説明変数です。赤色から青色は説明変数(Feature)が大きい、小さい場合に対応します。例えば、RMは部屋数を表し、その赤点が右側に多くプロットされています。このことより、部屋数が多い程、住宅価格が高くなる傾向であると理解することが出来ます。

そして、上記結果より感度の高い因子を元に、下図のように行列散布図で傾向を視覚的に確認することもできます。

下図は今回のベスト回帰モデルCatBoostの適合度を視覚的に確認したグラフです。横軸は生データで縦軸は回帰モデルによる予測結果です。比例y=x線上に点が乗ることが適合度が高いことを意味します。実機データを分析する場合は過剰適合は避けます。理由はばらつきによる飛び値があることが普通のためです。一方、シミュレーション結果をデータ分析する場合は過剰適合するくらいが良いです。理由はシミュレーション結果は理屈による数値計算結果であるため、飛び値などはあり得ないためです。但し、シミュレータのバグや数値計算誤差などによって、そういうのが発生するかもしれないことは頭に入れておき、行列散布図等で計算結果を確認することは必要です。

参考までに、下図は線形重回帰分析による適合度を可視化したグラフです。傾向は再現していますが、絶対値の当てはまりがよくないことが見た目にもわかります。

本コードを実行すると下図のように様々なデータを作成します。分析するデータのサマリー、行列散布図、基本統計量。他には、PyCaretによるベスト5の回帰モデル、ベストに選定された回帰モデルの特徴量の感度を可視化した図等です。PyCaretによるautoml()で選定されたベスト回帰モデルは.pklファイルで保存します。回帰モデルの利用は、そのオブジェクトファイルを読み出して使用します。

最後に、PyCaretが優れているのは実行が簡単なだけでなく、分析結果も簡単に行えます。plot_mode(best, plot='ほげほげ')といった僅か一行のコードで、次のように相関図、残渣図、特徴量のランキング図等を作成することが出来ます。



■本プログラム
Jupyter Labのコードを.pyファイルへ変換して一部編集したものです。Anacona PromptなどのCUIで実行できます。
import pandas as pd
file_path = './boston_dataset.csv'
df = pd.read_csv(file_path,
sep=',',
skiprows=0,
header=0,
encoding='utf-8')
print(df)
print(df.isna().sum())
import datetime
now = datetime.datetime.now()
now = now.strftime("%y%m%d")
df_describe = df.describe()
df_describe.to_csv(now + '_describe.csv',
sep=',',
index=True,
encoding='utf-8')
print(df_describe)
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import warnings
warnings.resetwarnings()
warnings.simplefilter('ignore', UserWarning)
fig = plt.figure(figsize=(8, 8), dpi=120)
ax = fig.add_subplot(1, 1, 1)
plt.subplots_adjust(left=0.1, right=0.95, bottom=0.1, top=0.95)
df.hist(ax=ax)
plt.tight_layout()
plt.savefig(now + '_hist.png')
plt.close()
from pycaret.regression import *
set_data = setup(data = df,
normalize = False,
train_size = 0.75,
session_id = 1,
target = 'PRICE',
silent=True,
)
top5 = compare_models(n_select = 5,
sort = 'MAE',
verbose = False)
with open(now + '_top5.txt', 'a') as f:
for i, my_model in enumerate(top5, start=1):
f.write(str(i) + '\n')
f.write(str(my_model) + '\n')
tuned_top5 = [tune_model(i, verbose = False) for i in top5]
with open(now + '_tuned_top5.txt', 'a') as f:
for i, my_model in enumerate(tuned_top5, start=1):
f.write(str(i) + '\n')
f.write(str(my_model) + '\n')
best = automl(optimize = 'MAE')
predict_model(best)
plot_model(best, plot = 'error')
plt.savefig(now + '_error.png')
plt.close()
plot_model(best, plot = 'residuals')
plt.savefig(now + '_residuals.png')
plt.close()
plot_model(best, plot = 'feature')
plt.savefig(now + '_feature.png')
plt.close()
plot_model(best, plot = 'manifold')
plt.savefig(now + '_manifold.png')
plt.close()
plot_model(best, plot = 'parameter')
plt.savefig(now + '_parameter.png')
plt.close()
try:
interpret_model(best)
plt.savefig(now + '_interpret_model.png')
plt.close()
#interpret_model(best, plot = 'reason', observation = 0)
except:
pass
import seaborn as sns
df2 = df.loc[:, ['LSTAT', 'RM', 'PTRATIO', 'AGE', 'PRICE']]
sns.set_context('talk')
mypairplot = sns.pairplot(df2,
kind="reg",
markers=".",
diag_kind="kde",
diag_kws=dict(shade=True),
)
plt.savefig(now + '_pairplot.png')
plt.close()
save_model(best, model_name = 'best_model_' + now)
load_reg_model = load_model('best_model_' + now)
print(load_reg_model)
predict_df = predict_model(load_reg_model, data=df)
print(predict_df)
from pycaret.utils import check_metric
metric_list = ['MAE', 'MSE', 'RMSE', 'R2', 'RMSLE']
for metric in metric_list:
print(metric, check_metric(predict_df.PRICE, predict_df.Label, metric))
predict_df.to_csv(now + '_predict.csv', index=False, encoding='utf-8')
import numpy as np
from scipy import optimize
import math
def approximation_expression(x, a):
return a * x
def plot_reg(DF, target_Yname, pred_Yname):
Y = DF[target_Yname].values.tolist()
pred_Y = DF[pred_Yname].values.tolist()
popt, pcov = optimize.curve_fit(approximation_expression, Y, pred_Y)
ax = plt.figure(num=0, dpi=120).gca()
ax.set_title("pred vs real ", fontsize=14)
ax.set_xlabel(target_Yname, fontsize=14)
ax.set_ylabel("Pred\n" + target_Yname, fontsize=14)
rp = ax.scatter(x = target_Yname,
y = pred_Y,
data = DF,
facecolors="none",
edgecolors='black')
y_min = DF[target_Yname].min()
y_max = DF[target_Yname].max()
y_pred_min = DF[pred_Yname].min()
y_pred_max = DF[pred_Yname].max()
x_min = min(y_min, y_pred_min)
x_max = max(y_max, y_pred_max)
x_range = x_max - x_min
print("x_range = x_max - x_min = " + str(x_range))
if x_max > 1:
print("range A")
min_lim = 0
if x_range <= 10:
max_lim = math.floor(x_max + 1)
else:
max_lim = math.floor(x_max + 10)
elif x_max >= 0:
print("range B")
min_lim = -0.1
max_lim = 1.1
elif x_min >= -1:
print("range C")
min_lim = -1.1
max_lim = 0.1
else:
print("range D")
max_lim = 0.1
if x_range <= 100:
min_lim = math.floor(x_min - 1)
else:
min_lim = math.floor(x_min - 10)
print(min_lim, max_lim)
rp.axes.set_xlim(min_lim, max_lim)
rp.axes.set_ylim(min_lim, max_lim)
x_approximation = np.linspace(min_lim, max_lim, 10)
y_approximation = popt[0] * x_approximation
line_approximation = ax.plot(x_approximation, y_approximation,
linestyle = 'dashed', linewidth = 3, color='r')
rp.axes.set_aspect('equal', adjustable='box')
plt.grid(True)
ax.legend([line_approximation[0]], ["y = {:.3f}x".format(popt[0])],
loc='upper left',
numpoints=1,
fontsize=14)
plt.tick_params(labelsize=14)
plt.tight_layout()
plt.savefig(now + '_real_vs_pred.png')
plt.close()
plot_reg(predict_df, 'PRICE', 'Label')
本記事の雛形コードをJupyter Labで順番に実行して、機械学習の手順がわかる動画を下記リンク先に作成しました。
youtu.be
以上
<広告>
リンク