Matplotlibでstatsmodelsの線形回帰（OLS）を美しくプロットする方法

statsmodelsの線形回帰（OLS）の結果は、非線形の曲線としてもプロットできます。データ自体は線形モデルに基づいていますが、回帰曲線や信頼区間を組み合わせることで、より見やすく洗練されたグラフを作成できます。この記事では、Matplotlibを使ってOLS回帰結果を美しく可視化する具体的な手順を解説します。

実装の手順

• 図のサイズを設定し、figure.autolayoutでサブプロット間および周囲の余白（パディング）を自動調整します。

• np.random.seed()メソッドで乱数シードを固定し、実行ごとに同じ結果が得られるようにします。

• サンプル数（nsample）とノイズの大きさ（sig）の変数を初期化します。

• numpyを使って、線形データポイント x、デザイン行列 X、係数 beta、真の値 y_true、観測値 y を作成します。

• res は最小二乗法（Ordinary Least Squares）クラスのインスタンスで、fit()により回帰モデルをフィッティングします。

• wls_prediction_std()を使って予測の標準偏差を計算します。なお、予測の信頼区間はWLSとOLSに適用可能ですが、独立ではあるものの同一分布に従わない観測値を扱う一般的なGLSには適用できません。

• plt.subplots()メソッドで、図（Figure）とサブプロット（Axes）のセットを作成します。

• plot()メソッドで、(x, y) の散布図、(x, y_true) の真の曲線、(x, res.fittedvalues) のOLS近似曲線、さらに (x, iv_u) と (x, iv_l) の信頼区間上限・下限をそれぞれ描画します。

• legend(loc='best')で、重なりを避けた最適な位置に凡例を配置します。

• 最後にshow()メソッドを呼び出して、図を表示します。

コード例

import numpy as np
from matplotlib import pyplot as plt
from statsmodels import api as sm
from statsmodels.sandbox.regression.predstd import wls_prediction_std

plt.rcParams["figure.figsize"] = [7.50, 3.50]
plt.rcParams["figure.autolayout"] = True

np.random.seed(9876789)
nsample = 50
sig = 0.5
x = np.linspace(0, 20, nsample)
X = np.column_stack((x, np.sin(x), (x - 5) ** 2, np.ones(nsample)))
beta = [0.5, 0.5, -0.02, 5.]
y_true = np.dot(X, beta)
y = y_true + sig * np.random.normal(size=nsample)

res = sm.OLS(y, X).fit()
prstd, iv_l, iv_u = wls_prediction_std(res)

fig, ax = plt.subplots()
ax.plot(x, y, 'o', label="data")
ax.plot(x, y_true, 'b-', label="True")
ax.plot(x, res.fittedvalues, 'r--.', label="OLS")
ax.plot(x, iv_u, 'r--')
ax.plot(x, iv_l, 'r--')
ax.legend(loc='best')

plt.show()

出力結果

上記のコードを実行すると、観測データ（青い点）、真の回帰曲線（青い実線）、OLSによる近似曲線（赤い破線）、そして赤い点線で示される95%信頼区間の上下限が1つのグラフに表示されます。これにより、モデルの当てはまり具合と予測の不確実性をひと目で確認できます。