「回帰分析から学ぶ計量経済学」をPythonで写経 ～ 第５章「時系列分析」③コクラン・オーカット法
第５章「時系列分析」
書籍の著者 山澤成康 先生
この記事は、書籍「回帰分析から学ぶ計量経済学」第５章「時系列分析」の Python写経活動 を取り扱います。
今回はコクラン・オーカット法に挑戦します。
残差に自己相関がある場合の問題解決になるようです。
では書籍を開いて回帰分析の旅に出発です🚀
はじめに
書籍「回帰分析から学ぶ計量経済学」のご紹介
このシリーズは書籍「回帰分析から学ぶ計量経済学: Excelで読み解く経済のしくみ」（オーム社、「テキスト」と呼びます）の Python 写経です。
テキストは、2023年11月に発売され、副題「Excelで読み解く経済のしくみ」のとおり、主に Excel を用いて、計量経済学を平易に学べる素晴らしい書籍です。
テキストの「はじめに」に著者の先生が執筆の動機を書かれています。
社会人の統計リテラシーの向上をテーマの１つとした科研費プロジェクトの最終年度で、広く社会人に向けてわかりやすい経済分析の本を書きたかったのです。
私にとって計量経済学は高嶺の花ですが、このテキストでさまざまな回帰分析のアプローチを知ることができました。
また、書籍の Excel 処理を Python に置き換える「寄り道写経」の実践を通じて、回帰分析のお気持ちに少し近づけた感じがいたします。
回帰分析に慣れ親しむのに丁度良いレベル感と内容ですので、これはぜひともブログにしたい！と思って現在に至ります。
計量経済学の色を薄め、データ分析の色を濃いめに書いてまいります！
引用表記
この記事は、出典に記載の書籍に掲載された文章と配布データを引用し、適宜、掲載文章・配布データを改変して書いています。
【出典】
「回帰分析から学ぶ計量経済学: Excelで読み解く経済のしくみ」
第１版第１刷、著者 山澤成康、オーム社
記事中のイラストは、「かわいいフリー素材集いらすとや」さんのイラストをお借りしています。
ありがとうございます！
第５章 時系列分析
この記事は第５章の以下の節を取り扱います。
5.10 コクラン・オーカット法の実例（5.9節を含みます）
記事に用いるデータは、オーム社の書籍紹介サイトからダウンロードできる Excel ファイル内のデータをもとにしてCSVファイルを作成し、data フォルダに格納しています。
第５章で用いるライブラリをインポートします。
### インポート
# 数値計算
import numpy as np
import pandas as pd
# 統計処理
import scipy.stats as stats
import statsmodels.formula.api as smf # フォーミュラ構文
import statsmodels.tsa.api as tsa # 時系列
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf # コレログラム等
import pmdarima as pm # ARIMA次数自動探索
import lmdiag # 残差プロット
# 描画
import matplotlib.pyplot as plt
import seaborn as sns
import japanize_matplotlib
# ワーニング非表示
import warnings
warnings.simplefilter('ignore')
5.10 コクラン・オーカット法の実例
■ コクラン・オーカット法の概要 p.167～
残差に自己相関がある場合、最小二乗法の仮定「誤差は互いに無相関」を満たさないため、最小二乗法が望ましい推定値になるための条件を満たさないことになります。
テキストは、誤差項に相関がある場合は、誤差項が自己相関していることを明示した式を使うコクラン・オーカット法を用いることで問題が解決できる、と説明しています。
■ コクラン・オーカット法の手順 p.168～
私の理解では、次の４つのステップでコクラン・オーカット法を実践できます。
>>> ステップ１：
通常の回帰分析を行い、残差$${u_t}$$を得る
$$
Y_t = \alpha + \beta X_t + u_t
$$
>>> ステップ２：
次の式で残差の回帰分析を行い、残差の１次自己相関係数$${\rho}$$を得る
$$
u_t = \rho u_{t-1} + \epsilon_t
$$
>>> ステップ３：
自己相関係数$${\rho}$$を用いて、次の式で回帰分析を行い、回帰係数$${\alpha^{\prime}, \beta^{\prime}}$$を得る
$$
\begin{align*}
Y_t^{\prime} &= \alpha^{\prime} + \beta^{\prime} X_{t}^{\prime} + v_t \\
ただし \\
Y_t^{\prime} &= Y_t - \rho Y_{t-1} \\
X_t^{\prime} &= X_{t} - \rho X_{t-1} \\
v_t &= \varepsilon_t - \rho \varepsilon_{t-1} \\
\alpha^{\prime} &= \alpha - \rho \alpha
\end{align*}
$$
>>> ステップ４：
以下の式で定数項$${\alpha}$$を得る
$$
\alpha = \alpha^{\prime} / (1 - \rho)
$$
このステップでテキストの四半期消費関数の推定を行います。
■ コクラン・オーカット法による四半期消費関数の推定 p.169～
データを読み込みます。
### データの読み込み
# CSVファイルの読み込み
df5 = pd.read_csv('./data/05_04_consumption.csv', index_col=0)
# 四半期日付のインデックスを設定
df5.index = pd.date_range(start='1994-03-31', periods=df5.shape[0], freq='QE')
df5.index.name = '四半期'
# データフレームの表示
print('df5.shape:', df5.shape)
display(df5.head())【実行結果】
実質GDPを$${X}$$、実質民間最終消費支出を$${Y}$$とし、以後のモデルを検討します。
① ステップ１：通常の回帰分析
### コクラン・オーカット法 ステップ１－１ 通常の回帰を実行
# 回帰分析の実行
result1 = smf.ols(formula='実質民間最終消費支出 ~ 実質GDP', data=df5).fit()
display(result1.summary())【実行結果】
テキストの記載の通り、標本サイズは$${117}$$、自由度調整済み決定係数は$${0.877}$$、GDPに係る係数（限界消費性向）は$${0.4853}$$です。
残差のダービン・ワトソン比の近似値は$${0.167}$$であり、残差に高い自己相関があることがわかります。
続いて、残差$${u_t}$$と残差のラグ１を取得します。
### コクラン・オーカット法 ステップ１－２ 残差の取得
## 残差データの準備
# 残差の取得
resid_df = pd.DataFrame({'残差': result1.resid})
# 残差のラグ１を取得
resid_df['残差1期前'] = resid_df['残差'].shift(1)
# 欠損値を含む行を削除
resid_df = resid_df.dropna()
# データの表示
print('resid_df.shape: ', resid_df.shape)
display(resid_df)【実行結果】
残差、残差ラグ１のデータを次のステップの回帰分析で利用します。
② ステップ２：残差の回帰分析
残差の平均は０なので、切片なしでモデリングします。
フォーミュラ文に「 0 + … 」とすることで切片なしにできます。
### コクラン・オーカット法 ステップ２－１ 残差のAR(1)モデルの回帰分析を実行
# 回帰分析の実行 ※切片なしモデル
result2 = smf.ols(formula='残差 ~ 0 + 残差1期前', data=resid_df).fit()
display(result2.summary())【実行結果】
残差ラグ１（残差１期前）の自己相関係数$${\rho}$$は$${0.9174}$$になりました。
rho を取得します。
### コクラン・オーカット法 ステップ２－２ 残差のAR(1)モデルの自己相関係数を取得
# 自己相関係数を取得
rho = result2.params[0]
print('rho = ', rho)【実行結果】
③ ステップ３：$${\rho}$$を用いて回帰分析を実行
次の回帰分析に用いる$${Y^{\prime}, X^{\prime}}$$を作成します。
### コクラン・オーカット法 ステップ３－１ Y'、X'の作成
## データの準備
# 結果を格納するデータフレームの準備
df5_dash = pd.DataFrame()
# Y'_t = Y_t - ρY_t-1の算出
df5_dash['Y_dash'] = (df5['実質民間最終消費支出']
- rho * df5['実質民間最終消費支出'].shift(1))
# X'_t = X_t - ρX_t-1の算出
df5_dash['X_dash'] = df5['実質GDP'] - rho * df5['実質GDP'].shift(1)
# 欠損値を含む行を削除
df5_dash = df5_dash.dropna()
# 結果表示
print('df5_dash.shape: ', df5_dash.shape)
display(df5_dash)【実行結果】
では回帰分析を実行します。
### コクラン・オーカット法 ステップ３－２ 回帰分析の実行
# 回帰分析の実行
result3 = smf.ols(formula='Y_dash ~ X_dash', data=df5_dash).fit()
display(result3.summary())【実行結果】
ネタはすべて揃いました！
④ ステップ４：係数の取得
$${\alpha = \alpha^{\prime} / (1 - \rho)}$$で定数項の係数$${\alpha}$$を得ます。
実質GDPの係数$${\beta^{\prime}}$$は③の回帰で推定した$${0.4654}$$です。
### コクラン・オーカット法 ステップ４ 係数の取得
# 係数α、β'の算出
alpha = result3.params[0] / (1 - rho)
beta_dash = result3.params[1]
print(f'実質民間最終消費支出 = {alpha:.2f} + {beta_dash:.6f} × 実質GDP')【実行結果】
コクラン・オーカット法で推定した式が出来上がりました！
ステップ１の通常の回帰の残差と、ステップ３のコクラン・オーカット法の回帰の残差を可視化して比べます。
テキストの図表「消費関数の残差」に相当します。
### 消費関数の残差の可視化 p.170のグラフ
# 描画領域の設定
plt.figure(figsize=(10, 3))
# 残差（通常の回帰）の折れ線グラフの描画
plt.plot(result1.resid, lw=0.9, label='残差')
# 残差（コクラン・オーカット法）の折れ線グラフの描画
plt.plot(result3.resid, lw=0.9, color='tab:red',
label='残差(コクラン・オーカット法)')
# 修飾
plt.axhline(0, color='black', lw=0.5)
plt.ylim(-15000, 20000)
plt.xlabel('時間 [四半期]')
plt.ylabel('残差 [10億円]')
plt.grid(lw=0.5)
plt.legend(loc='upper left')
plt.show()【実行結果】
青い線が通常の回帰の残差、赤い線がコクラン・オーカット法による残差です。
テキストによると、残差の自己相関が少なくなっていることがわかります、とのことです。
コレログラム、偏自己相関係数を見てみましょう。
まずは通常の回帰の残差です。
### 通常の回帰の残差のコレログラム・偏自己相関係数
plot_acf_pacf(result1.resid.dropna(), 20)【実行結果】
ラグ１に高い偏自己相関係数が見られます。
続いてコクラン・オーカット法の残差です。
### コクラン・オーカット法の残差のコレログラム・偏自己相関係数
plot_acf_pacf(result3.resid.dropna(), 20)【実行結果】
若干残っている感は否めませんが、自己相関・偏自己相関はかなり小さくなっています。
今回の写経は以上です。
シリーズの記事
次の記事
前の記事
目次
ブログの紹介
note で７つのシリーズ記事を書いています。
ぜひ覗いていってくださいね！
１．のんびり統計
統計検定２級の問題集を手がかりにして、確率・統計をざっくり掘り下げるブログです。
雑談感覚で大丈夫です。ぜひ覗いていってくださいね。
統計検定２級公式問題集CBT対応版に対応しています。
Python、EXCELのサンプルコードの配布もあります。
２．実験！たのしいベイズモデリング１＆２をPyMC Ver.5で
書籍「たのしいベイズモデリング」・「たのしいベイズモデリング２」の心理学研究に用いられたベイズモデルを PyMC Ver.5で描いて分析します。
この書籍をはじめ、多くのベイズモデルは Ｒ言語＋Stanで書かれています。
PyMCの可能性を探り出し、手軽にベイズモデリングを実践できるように努めます。
身近なテーマ、イメージしやすいテーマですので、ぜひぜひPyMCで動かして、一緒に楽しみましょう！
３．実験！岩波データサイエンス1のベイズモデリングをPyMC Ver.5で
書籍「実験！岩波データサイエンスvol.1」の４人のベイジアンによるベイズモデルを PyMC Ver.5で描いて分析します。
この書籍はベイズプログラミングのイロハをざっくりと学ぶことができる良書です。
楽しくPyMCモデルを動かして、ベイズと仲良しになれた気がします。
みなさんもぜひぜひPyMCで動かして、一緒に遊んで学びましょう！
４．楽しい写経 ベイズ・Python等
ベイズ、Python、その他の「書籍の写経活動」の成果をブログにします。
主にPythonへの翻訳に取り組んでいます。
写経に取り組むお仲間さんのサンプルコードになれば幸いです🍀
５．RとStanではじめる心理学のための時系列分析入門 を PythonとPyMC Ver.5 で
書籍「RとStanではじめる心理学のための時系列分析入門」の時系列分析をPythonとPyMC Ver.5 で実践します。
この書籍には時系列分析のテーマが盛りだくさん！
時系列分析の懐の深さを実感いたしました。
大好きなPythonで楽しく時系列分析を学びます。
６．データサイエンスっぽいことを綴る
統計、データ分析、AI、機械学習、Pythonのコラムを不定期に綴っています。
統計・データサイエンス書籍にまつわる記事が多いです。
「統計」「Python」「数学とPython」「R」のシリーズが生まれています。
７．Python機械学習プログラミング実践記
書籍「Python機械学習プログラミング PyTorch & scikit-learn編」を学んだときのさまざまな思いを記事にしました。
この書籍は、scikit-learnとPyTorchの教科書です。
よかったらぜひ、お試しくださいませ。
最後までお読みいただきまして、ありがとうございました。