インパルス応答で何がわかる？VARモデルの基本と評価法

# **user**:

VARモデルからインパルス応答を使って効果を予測する方法をPythonで教えてください。

# **assistant**:

VAR（Vector Autoregression）モデルからインパルス応答関数（Impulse Response Function, IRF）を使って、ある変数が他の変数に与える影響を予測するには、Pythonの`statsmodels`ライブラリを使うのが一般的です。以下に、基本的な手順とコード例を示します。

---

## ✅ 前提：必要なライブラリのインストール
```bash
pip install statsmodels pandas matplotlib
```

---

## ✅ ステップ1：データの準備
```python
import pandas as pd
from statsmodels.tsa.api import VAR

# 適当な多変量時系列データ（例：経済データ）
data = pd.read_csv("your_data.csv", index_col=0, parse_dates=True)
data = data.dropna()  # 欠損値処理
```

---

## ✅ ステップ2：VARモデルのフィッティング
```python
model = VAR(data)

# ラグの次数はAICやBICで選ぶのが一般的
lag_order = model.select_order(maxlags=10)
print(lag_order.summary())

# 選ばれたラグ数でフィット
results = model.fit(lag_order.aic)
```

---

## ✅ ステップ3：インパルス応答関数（IRF）の計算と表示
```python
# 10期間先までのインパルス応答を計算
irf = results.irf(10)

# 可視化（すべての組み合わせ）
irf.plot(orth=False)  # orth=Trueで直交化（Cholesky分解によるショックの独立化）
```

---

## ✅ ステップ4：特定の変数に対するIRFを抽出
たとえば、変数 `"A"` のショックが `"B"` に与える影響：
```python
irf.plot(impulse='A', response='B')
```

---

## 補足
- `orth=True`にすると、ショックが直交化され、他の変数との共通ショックが除去されます。
- 信頼区間付きで表示したい場合は：
```python
irf.plot(impulse='A', response='B', stderr_type='mc')  # 'mc'はモンテカルロ法
```

---

## 参考文献
- [statsmodels.tsa.api.VAR](https://www.statsmodels.org/stable/generated/statsmodels.tsa.api.VAR.html)
- VAR分析の教科書や時系列解析の基礎理論

---

必要であれば、ダミーデータを使った実行可能なサンプルコードも提供できます。ご希望ですか？

# **user**:

ダミーデータを使った実行可能なサンプルコードをお願いします

# **assistant**:

以下に、**ダミーデータを使ったVARモデルとインパルス応答分析の実行可能なPythonコード**を示します。すべてのコードをコピペしてそのまま実行できます。

---

## ✅ 実行可能なサンプルコード（VARモデル + インパルス応答）

```python
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.tsa.api import VAR

# ダミーデータの作成（2変数、100期間分）
np.random.seed(42)
n_obs = 100

# ランダムウォーク + ノイズ（適当な構造）
x = np.cumsum(np.random.normal(size=n_obs))  # 変数X
y = np.cumsum(np.random.normal(size=n_obs)) + 0.5 * x  # 変数YはXに依存

# データフレームにまとめる
data = pd.DataFrame({'X': x, 'Y': y})

# 時系列データとしてのインデックス
data.index = pd.date_range(start='2000-01-01', periods=n_obs, freq='M')

# 欠損除去（今回は不要だが一般的な処理）
data = data.dropna()

# VARモデルのフィッティング
model = VAR(data)

# ラグ次数の自動選択
lag_order = model.select_order(maxlags=10)
print("選ばれたラグ次数（AIC）:", lag_order.aic)

# モデルフィッティング（選ばれたラグを使う）
results = model.fit(lag_order.aic)

# インパルス応答の計算（10期間）
irf = results.irf(10)

# すべての変数ペアのインパルス応答をプロット
irf.plot(orth=False)  # orth=Trueだと直交化ショック
plt.tight_layout()
plt.show()

# 特定のインパルス応答：XのショックがYに与える影響
irf.plot(impulse='X', response='Y', orth=False)
plt.title("IRF: XのショックがYに与える影響")
plt.tight_layout()
plt.show()
```

---

### 🔍 説明

- `X`: 外生的なショック（ランダムウォーク）
- `Y`: `X` に依存して動く変数（因果的に影響されるイメージ）
- `VAR` モデルをフィットして、`X` が `Y` に与える影響をインパルス応答で可視化
- `irf.plot(impulse='X', response='Y')`で、**「Xが1単位ショックを受けたらYがどう変化するか」**を見られます

---

必要であれば、**信頼区間付きプロット**や**直交化ショック（Cholesky分解）**も追加可能です。ご希望ですか？

# **user**:

ある変数が他の変数に与える影響の大きさを実数で得られないですか？

# **assistant**:

はい、**インパルス応答の影響の大きさを実数（数値）で取得することは可能**です。Pythonの`statsmodels`では、`.irf()`の返すオブジェクトから、インパルス応答値を**数値配列（NumPy配列）やPandasデータフレーム**として取得できます。

---

## ✅ 実数の影響値を取得する方法

```python
# インパルス応答オブジェクト（すでに fit 済み）
irf = results.irf(10)

# 変数 'X' のショックが 'Y' に与えるインパルス応答値（配列として取得）
response_values = irf.irfs[:, data.columns.get_loc('Y'), data.columns.get_loc('X')]

# 表形式で確認（10期間分）
df_response = pd.DataFrame({
    'Period': range(1, len(response_values) + 1),
    'X→Y Impact': response_values
})

print(df_response)
```

---

### 🔍 出力例（意味）

| Period | X→Y Impact |
|--------|------------|
|   1    |  0.123     |
|   2    |  0.245     |
|   3    |  0.301     |
|  ...   |   ...      |

- これは、**Xに1単位のショックを与えたとき、Yが各期間にどれだけ変化するか（単位で）**を示します。
- 単に合計の影響を知りたい場合は、次のように合計できます：

```python
total_effect = response_values.sum()
print("XのショックによるYへの累積影響:", total_effect)
```

---

## ✅ その他の方法（累積インパルス応答）

```python
# 累積インパルス応答（積み上げ効果）
cumulative_irf = irf.cum_effects[:, data.columns.get_loc('Y'), data.columns.get_loc('X')]
```

---

### ご参考までに
`irf.irfs`の形は `(期間, 応答変数, ショック変数)` です。

---

必要であれば、**全変数ペアの影響をDataFrameで一覧表示するコード**も提供できます。見たいですか？

# **user**:

irf = results.irf(10) としたときの結果と、irf = results.irf(20)したときの先頭１０レコードは同一ですか？

# **assistant**:

はい、**`results.irf(10)` の結果と `results.irf(20)` の結果の**「先頭10レコード（期間1〜10のインパルス応答）」**は**同一になります**。

---

### 🔍 理由

`results.irf(n)` は **VARモデルに基づいて、n期間先までのインパルス応答を一括計算する関数**です。これは「n期間分を一度に計算する」だけであり、1〜10期間目の値は `results.irf(10)` でも `results.irf(20)` でもまったく同じVARパラメータを用いて算出されます。

---

### ✅ 実証コード（確認用）

```python
irf10 = results.irf(10).irfs
irf20 = results.irf(20).irfs

# 比較：先頭10レコードが等しいか（全変数ペアについて）
np.allclose(irf10, irf20[:10])
```

結果が `True` なら、**先頭10期間のインパルス応答は完全に一致**しています。

---

### ✅ 注意点

- `results.irf(10)` は **10期間までしか持たない**ため、11期間以降の応答は得られません。
- 長い期間の応答を分析する場合は `irf = results.irf(n)` の `n` を大きめにしておくと良いです。

---

必要であれば、**`irfs` 配列の構造を可視化するコード**も提供できます。ご希望ですか？

# **user**:

irfsは、説明変数・目的変数の区別なく、任意の２変数間の影響を見ることができるという理解であってますか？

# **assistant**:

はい、その理解は**概ね正しい**です。

---

## ✅ 結論
`irfs`（インパルス応答関数）は、**VARモデルに含まれるすべての変数間の影響を対等に扱い**ます。したがって：

> **「説明変数・目的変数の区別なく、任意の2変数間で“ある変数のショックが他の変数にどう影響するか”を見る」**  
ことが可能です。

---

## 🔍 なぜ「説明変数・目的変数の区別がない」のか？

- VAR（Vector Autoregression）モデルは、**すべての変数を同時に内生変数（＝互いに影響し合う変数）**として扱います。
- 回帰のように「従属変数（Y）」「独立変数（X）」という明確な役割分担はしません。
- よって、`irfs[i, j]` は「**変数 j に1単位のショックがあったとき、変数 i がどう反応するか**」を意味します。

---

## ✅ 具体例（2変数 X, Y の場合）

```python
irfs[:, i, j]
```

| 項目 | 意味 |
|------|------|
| `j` | **ショックを受けた変数（インパルス）** |
| `i` | **その影響を受ける変数（レスポンス）** |
| `:` | 時間（ステップ）方向 |

---

## ✅ 実務上の補足

- **因果解釈が目的**なら、LiNGAMなどの因果モデルの方が適しています。
- VARは、**時系列内での動的な影響を数量的に評価する手法**です（因果性ではなく、**動学的応答性**に注目）。

---

### ご参考
もし因果性も考慮したい場合には、以下の選択肢があります：

- Granger因果性検定（VARに基づく）
- LiNGAM（非ガウス性を仮定した構造モデル）
- SVAR（構造VAR、経済学的制約を導入）

---

必要であれば、**VARモデルとLiNGAMとの違い**についても整理できます。興味ありますか？

# **user**:

はい、VARモデルとLiNGAMとの違いについて整理してください。

# **assistant**:

もちろんです。以下に、**VARモデルとLiNGAMの違い**を目的・前提・推定方法・得られる結果などの観点で整理しました。

---

## ✅ VARモデル vs LiNGAMモデルの違いまとめ

| 観点 | VAR（Vector Autoregression） | LiNGAM（Linear Non-Gaussian Acyclic Model） |
|------|------------------------------|--------------------------------------------|
| **目的** | **時系列データにおける動的な影響（時間の流れに沿った依存関係）**の分析 | **因果構造（変数間の因果関係）**の推定 |
| **対象** | 時系列データ（同時点での影響は扱わない） | 時系列 or 横断データ（特に同時点での因果性に強い） |
| **変数の扱い** | すべての変数を内生変数として対等に扱う | 変数に**因果的な順序（非巡回性 DAG）**を仮定 |
| **前提** | - 線形性<br>- 定常性<br>- ラグ構造での自己回帰関係 | - 線形性<br>- **非ガウス分布**（重要）<br>- **非巡回的（ループなし）な因果構造** |
| **モデル式** | $$ X_t = A_1 X_{t-1} + A_2 X_{t-2} + \dots + \varepsilon_t $$ | $$ X = B X + e $$（$B$は因果構造、$e$は非ガウスな独立ノイズ） |
| **因果性の推定** | Granger因果性で**時間遅れの影響**を検定可能だが、因果と呼ぶには限定的 | **構造方程式から明示的に因果方向を推定**可能 |
| **インパルス応答** | あり（ショックの時間的影響を評価） | なし（時間的変化よりも構造的関係に注目） |
| **代表的な利用分野** | 経済学（マクロ経済、金融の時系列分析） | 社会科学、心理学、医療、マーケティングなど幅広い分野での**因果探索** |
| **主なPythonライブラリ** | `statsmodels` | `lingam`（`causal-learn`, `py-causal` なども） |

---

## 🔍 補足：選ぶ基準のポイント

| 目的 | 推奨モデル |
|------|------------|
| 時系列の動学的な変化を見る（例：政策ショックの影響） | VARモデル |
| 「何が原因で何が結果か？」を知りたい（因果構造を解明したい） | LiNGAM |

---

## ✅ 両者を組み合わせる考え方もある

例えば、**VAR-LiNGAM**という応用モデルがあります：

- 時系列データを使ってVAR構造をモデル化しつつ、
- **その残差構造にLiNGAMを適用して因果構造を推定**する方法。

これにより、時間遅れの影響（VAR）と**即時的な因果構造（LiNGAM）**を両立できます。

---

必要であれば、**LiNGAMをPythonで使う実行サンプルコード**も提供できます。ご希望ですか？

# **user**:

VARモデルの良し悪しはどのように評価しますか？また、LiNGAMモデルの良し悪しはどのように評価しますか？

# **assistant**:

良い質問です。VARモデルとLiNGAMモデルは目的が異なるため、**評価基準も異なります**。それぞれについて、以下に整理して説明します。

---

## ✅ VARモデルの良し悪しの評価方法

### 1. **ラグ次数の選定**
- AIC（赤池情報量基準）、BIC（ベイズ情報量基準）、FPE、HQIC などで評価。
- `model.select_order()` により自動選択可能。

### 2. **残差診断（モデルの妥当性）**
- **残差の自己相関がないこと（白色ノイズ）**を確認：
  - Ljung-Box検定（`statsmodels.stats.diagnostic.acorr_ljungbox`）
  - 残差の自己相関関数（ACF）プロット
- **残差の正規性**（Jarque-Bera検定など）
- **定常性チェック**（ユニットルート検定やルートの絶対値 < 1）

### 3. **予測性能**
- クロスバリデーションや直近データによる予測精度（RMSE, MAE）などで評価可能。
- `results.forecast()` を用いて未来予測 → 実測値と比較。

### 4. **インパルス応答の解釈可能性**
- 経済・業務の知識と一致しているか（ドメイン知識と照らし合わせて解釈可能か）

---

## ✅ LiNGAMモデルの良し悪しの評価方法

### 1. **推定された因果構造の信頼性**
- 結果として得られる**因果順序と因果係数の解釈性**
- ドメイン知識と一致するか？

### 2. **統計的仮定の妥当性**
- **非ガウス性が前提**なので、変数や誤差項が非正規分布であるかを確認：
  - Q-QプロットやShapiro-Wilk検定などで正規性を否定

### 3. **再現性と頑健性**
- **ブートストラップ法による信頼区間評価**（`BootstrapLiNGAM`）
- 少しのデータ変動で推定結果が安定しているか

### 4. **評価指標（専用）**
- `lingam` ライブラリでは以下のスコアが得られる：
  - **推定係数の有意性（p値）**
  - **有向非巡回グラフ（DAG）におけるエッジのスコアや数**
  - グラフ構造の正解が分かっている場合は、**構造的Hamming距離（SHD）**などで比較可能

---

## ✅ 評価観点の違い（まとめ）

| 観点 | VARモデル | LiNGAMモデル |
|------|-----------|--------------|
| 主な目的 | 時系列の動的予測 | 因果構造の発見 |
| 評価軸 | ラグ次数、残差診断、予測精度、IRFの妥当性 | 非ガウス性の仮定、構造の安定性、ブートストラップ結果、ドメイン妥当性 |
| 指標例 | AIC/BIC、RMSE、ACF、IRF解釈 | ブートストラップ信頼区間、因果順序の解釈性、非正規性検定 |

---

必要であれば、LiNGAMの「ブートストラップ付き因果係数の評価コード」も提供できます。ご希望ですか？

# **user**:

ありがとうございます。今日はこの辺で大丈夫です。今回の会話の内容をブログに書き残したいのですが、タイトル案をいくつか提案していただけますか？

# **assistant**:

どういたしまして！以下に、今回の会話内容（VARモデルとLiNGAMの違いや評価方法、インパルス応答分析など）を反映した**ブログタイトル案**をいくつか提案します。

---

## ✨ タイトル案一覧

### VARモデルと因果探索にフォーカスしたもの
1. **「時系列か因果か？VARモデルとLiNGAMを徹底比較してみた」**
2. **「インパルス応答で何がわかる？VARモデルの基本と評価法」**
3. **「VARとLiNGAMの違いをChatGPTに聞いてみた」**
4. **「動的影響と因果関係をどう見極める？VARモデルとLiNGAMの使い分け」**

### 分析＋実践寄りのトーン
5. **「Pythonで始めるVARモデルとLiNGAM：分析・評価・使いどころ」**
6. **「インパルス応答×因果探索：VARとLiNGAMをPythonで実践」**
7. **「モデル選びに迷ったら：VARかLiNGAMか、その違いと見極め方」**

### ちょっとキャッチーなタイトル
8. **「全部『因果』じゃない！？VARとLiNGAMの違いを直感的に理解する」**
9. **「時間遅れか構造依存か？VARとLiNGAMで世界の見え方が変わる」**

---

ご希望に応じて、記事の**構成案**や**導入文のドラフト**もお手伝いできます。書き進める際はぜひお声がけください！

記事の主な読者層（例：データ分析初心者／経済系研究者など）が決まっていれば、それに合わせたタイトルにも調整できますが、いかがですか？

# **user**:

了解しました！以下に、「**インパルス応答で何がわかる？VARモデルの基本と評価法**」というタイトルにふさわしい、**今回の会話内容の詳細なまとめ**をご用意しました。ブログ記事の本文としてそのまま使えるように整理しています。

---

# インパルス応答で何がわかる？VARモデルの基本と評価法

時系列データに含まれる複数の変数が、互いにどのように影響し合っているのかを分析する方法として、「**VAR（Vector Autoregression）モデル**」は広く使われています。今回は、Pythonを使ってVARモデルを構築し、インパルス応答の結果から得られる示唆や、モデルの評価方法について整理しました。また、因果関係の発見に特化したLiNGAMモデルとの違いについても触れています。

---

## VARモデルとは何か？

VARモデルは、複数の時系列変数が**互いにラグ（過去の値）を通じて影響を与え合う**ことを前提に構築されます。特徴は以下の通りです。

- **すべての変数が内生変数（説明・目的の区別なし）**
- **ラグ構造**（過去の値が現在の変動を説明）
- 多変量時系列に特化

モデル形式は以下のようになります：

$$
X_t = A_1 X_{t-1} + A_2 X_{t-2} + \dots + A_p X_{t-p} + \varepsilon_t
$$

---

## インパルス応答とは？

インパルス応答（Impulse Response Function）は、**ある変数に1単位の外的ショックを加えたときに、他の変数が時間とともにどう反応するか**をシミュレーションするものです。

たとえば、2変数モデルで `irfs[:, i, j]` を見ると：

- `j`: ショックを受けた変数（インパルス）
- `i`: その影響を受ける変数（レスポンス）
- `:`: 時間の経過（0期先、1期先、...）

となり、**任意の変数ペア間で「ショック→反応」の関係を確認できます**。

---

## Pythonでの実装と注意点

`statsmodels` ライブラリを使うことで、Pythonで簡単にVARモデルとインパルス応答分析が行えます。以下は重要なポイントです：

- `results.irf(n)` は n期先までの応答を返すが、**1〜10期目の値は irf(10) でも irf(20) でも同じ**
- `irfs` はすべての変数ペア間の応答を出力するので、**説明変数・目的変数の区別は不要**
- 解釈にはドメイン知識が不可欠（モデルは因果を保証しない）

---

## VARモデルの評価方法

1. **ラグ次数の選定**  
   → AIC、BIC、FPE、HQIC などの情報量基準を使用

2. **残差診断**
   - Ljung-Box検定：残差に自己相関がないか
   - 正規性検定（Jarque-Beraなど）
   - 定常性（ユニットルートが存在しないか）

3. **予測精度**
   - RMSE、MAEなどで将来値の予測性能を評価

4. **インパルス応答の妥当性**
   - 経済的・実務的な知識と照らし合わせて違和感がないかを確認

---

## LiNGAMモデルとの違い

VARモデルは**「時間の流れによる影響」**を分析するのに適しており、因果性を直接推定するわけではありません。一方、**LiNGAM（Linear Non-Gaussian Acyclic Model）**は、以下のような特徴を持ちます：

| 特徴 | VARモデル | LiNGAM |
|------|-----------|--------|
| 主目的 | 時系列内の影響分析 | **因果構造の推定** |
| 変数関係 | 全て内生的 | 因果順序あり（非巡回グラフ） |
| 分布仮定 | 正規でもOK | **非ガウス性必須** |
| 時系列前提 | 必須 | あってもなくても可 |
| インパルス応答 | あり | なし |

LiNGAMは、非ガウス分布のノイズという前提をもとに、**変数間の因果的な順序と係数を明示的に推定**します。ブートストラップを使った信頼区間の評価なども可能です。

---

## まとめ

- VARモデルは、**時系列における変数間の動的な関係性**を把握するために有効
- インパルス応答を用いることで、**ある変数のショックが他変数に与える影響を時間軸で視覚化**できる
- 良いVARモデルには、**残差診断や予測精度の高さ、解釈の妥当性**が必要
- 因果性を直接推定したい場合は、LiNGAMなど**構造モデル系のアプローチ**が適している

---

次回は、実データを使ったLiNGAMモデルの実装や、VARとLiNGAMを組み合わせる「VAR-LiNGAM」的アプローチにも触れてみたいと思います。

---

ブログの導入文やコード付きの具体例が必要な場合も、お気軽にご相談ください！

# **user**:

ありがとうございました。今日はこの辺にしましょう。次回またよろしくお願いします！

# **assistant**:

こちらこそ、ありがとうございました！とても深い内容でしたね。次回もいつでもお待ちしています。どうぞ良い一日をお過ごしください😊