LoginSignup
1
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

@maskot1977(Ikemen Mas Kot)

時系列の感度分析を滑らかに：Pythonで実装するEWMA × 平滑ラグ平均

Posted at

user:

EWMA/平滑ラグ平均について、Pythonを用いて詳しく教えてください。

assistant:

もちろんです。「EWMA（Exponentially Weighted Moving Average：指数加重移動平均）」や「平滑ラグ平均」は、時系列データの平滑化やラグの情報をまとめるのに使われるテクニックです。以下では、基本概念からPythonでの実装・応用まで丁寧に解説します。

✅ 1. EWMAとは？

● 定義

指数加重移動平均（EWMA）は、最近のデータほど重みを大きく、過去のデータほど重みを小さくする移動平均の一種です。

重みは指数関数的に減衰します。

$$
\text{EWMA}t = \alpha \cdot x_t + (1 - \alpha) \cdot \text{EWMA}{t-1}
$$

  • $\alpha$：平滑化係数（0 < α ≤ 1）
  • $x_t$：時点 $t$ の観測値
  • $\text{EWMA}_t$：時点 $t$ の指数加重移動平均

✅ 2. Pythonでの基本実装（pandas）

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

# サンプル時系列データ
np.random.seed(42)
n = 100
x = np.random.randn(n).cumsum() + 10
ts = pd.Series(x)

# EWMAの計算
ewma_ts = ts.ewm(span=10, adjust=False).mean()

# 可視化
plt.plot(ts, label='Original')
plt.plot(ewma_ts, label='EWMA (span=10)', color='orange')
plt.legend()
plt.title("指数加重移動平均")
plt.show()

● ewm(span=10)とは？

  • span: 移動平均における「期間」に相当。大きいほど滑らかに。
  • adjust=False: 重みを「オンライン更新方式」で計算（adjust=Trueは全履歴を重視）

✅ 3. 平滑ラグ平均（Smoothed Lag Average）

● 目的

  • ラグ変数 $x_{t-1}, x_{t-2}, \dots, x_{t-k}$ があるときに、それらを単に平均せず、重み付けして1つの特徴量にまとめたいという状況があります。
  • 例：因果推論・時系列回帰で「ラグを滑らかにまとめて使いたい」

✅ 4. 平滑ラグ平均をEWMAで作る例

● 例：過去3ラグの重み付き平均を作る（逐次的に）

def smoothed_lag_average(series, max_lag=3, alpha=0.5):
    """
    時系列のラグを指数加重平均で滑らかにまとめる
    """
    df = pd.DataFrame({'x': series})
    for lag in range(1, max_lag + 1):
        df[f'x_lag{lag}'] = df['x'].shift(lag)
    
    # NaN除去
    df.dropna(inplace=True)

    # EWMA的に加重平均を取る
    weights = np.array([(1 - alpha)**(i-1) for i in range(1, max_lag+1)])
    weights = weights / weights.sum()

    df['smoothed_lag'] = df[[f'x_lag{i}' for i in range(1, max_lag + 1)]].dot(weights)
    return df[['x', 'smoothed_lag']]

# 実行例
result = smoothed_lag_average(ts, max_lag=5, alpha=0.5)

# プロット
plt.plot(result['x'], label='Original')
plt.plot(result['smoothed_lag'], label='Smoothed Lag Avg', color='red')
plt.legend()
plt.title("平滑ラグ平均（EWMA風）")
plt.show()

✅ 5. 応用：特徴量としての活用

滑らかにしたラグ平均は、次のような用途に使えます。

  • 時系列予測（例：ARモデルの特徴量として）
  • 因果推論（ラグ変数を1つに集約しつつ、影響の強弱を残す）
  • 機械学習モデルの説明変数（XGBoostやLSTM等）
# 例：平滑ラグ平均を説明変数に
from sklearn.linear_model import LinearRegression

X = result['smoothed_lag'].values.reshape(-1, 1)
y = result['x'].values

model = LinearRegression().fit(X, y)
print("係数:", model.coef_)
print("切片:", model.intercept_)

✅ 6. パラメータ選択のコツ

パラメータ 意味 推奨設定の目安
alpha 重みの減衰速度 0.1〜0.5程度。小さいと過去を重視
max_lag 使用するラグの数 3〜20程度。時系列の性質に依存
span pandas.ewm()での平滑化幅 自由。トレンド除去なら10〜20程度

✅ 7. まとめ

項目 内容
EWMA 時系列の滑らか化に使う指数加重平均
平滑ラグ平均 複数ラグをEWMA風に加重して1変数に集約
利点 ノイズ除去・次元削減・ラグのまとめ
実装 pandas.ewm() or 自作関数（加重平均）

さらに応用的なケース（例えば、異なる系列で平滑化されたラグを組み合わせて因果効果を見る、など）についても興味があればご相談ください！

user:

ありがとうございます。EWMA/平滑ラグ平均に基づいて、次のPython関数を改変してください。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from sklearn.utils import resample

def sensitivity_analysis_with_ci_decay(
    df, date_col, target_var, explanatory_vars,
    intervention_scale=0.01, lags=12,
    horizons=range(1, 13), n_components=0.95,
    n_bootstrap=100, random_state=42,
    decay_rate=None, show_plot=True
):
    """
    Horizon別感度分析 + 95%信頼区間 + 減衰ウェイト付きラグ項
    
    decay_rate: float（例：0.95）を指定するとラグに減衰重みを掛ける。
    """

    np.random.seed(random_state)
    df = df.sort_values(date_col).copy()
    results_all = {}

    for horizon in horizons:
        df_lagged = df.copy()

        for lag in range(1, lags + 1):
            decay_weight = decay_rate**(lag - 1) if decay_rate is not None else 1.0
            df_lagged[f'{target_var}_lag{lag}'] = df_lagged[target_var].shift(lag) * decay_weight
            for var in explanatory_vars:
                df_lagged[f'{var}_lag{lag}'] = df_lagged[var].shift(lag) * decay_weight

        df_lagged['target_future'] = df_lagged[target_var].shift(-horizon)
        df_lagged.dropna(inplace=True)

        feature_cols = [col for col in df_lagged.columns if col not in [date_col, target_var, 'target_future']]
        X = df_lagged[feature_cols].copy()
        y = df_lagged['target_future'].copy()
        X_latest = X.iloc[[-1]].copy()

        baseline_preds = []
        intervened_preds_dict = {var: [] for var in explanatory_vars}

        for _ in range(n_bootstrap):
            X_resampled, y_resampled = resample(X, y)

            scaler = StandardScaler()
            X_scaled = scaler.fit_transform(X_resampled)
            pca = PCA(n_components=n_components)
            X_pca = pca.fit_transform(X_scaled)
            model = LinearRegression()
            model.fit(X_pca, y_resampled)

            X_latest_scaled = scaler.transform(X_latest)
            X_latest_pca = pca.transform(X_latest_scaled)
            y_pred_baseline = model.predict(X_latest_pca)[0]
            baseline_preds.append(y_pred_baseline)

            for var in explanatory_vars:
                X_modified = X_latest.copy()
                for lag in range(1, lags + 1):
                    col = f'{var}_lag{lag}'
                    if col in X_modified.columns:
                        X_modified[col] *= (1 + intervention_scale)
                X_modified_scaled = scaler.transform(X_modified)
                X_modified_pca = pca.transform(X_modified_scaled)
                y_pred_intervened = model.predict(X_modified_pca)[0]
                intervened_preds_dict[var].append(y_pred_intervened)

        baseline_preds = np.array(baseline_preds)
        result = {}
        for var in explanatory_vars:
            diffs = np.array(intervened_preds_dict[var]) - baseline_preds
            result[var] = {
                'mean': np.mean(diffs),
                'lower': np.percentile(diffs, 2.5),
                'upper': np.percentile(diffs, 97.5)
            }

        results_all[horizon] = result

    # 結果の整形
    df_result = pd.DataFrame({
        (var, metric): [results_all[h][var][metric] for h in horizons]
        for var in explanatory_vars
        for metric in ['mean', 'lower', 'upper']
    }, index=horizons)
    df_result.index.name = "horizon"

    # グラフ表示
    if show_plot:
        plt.figure(figsize=(12, 6))
        for var in explanatory_vars:
            mean = df_result[(var, 'mean')]
            lower = df_result[(var, 'lower')]
            upper = df_result[(var, 'upper')]
            plt.plot(mean.index, mean.values, label=f'{var} (mean)')
            plt.fill_between(mean.index, lower.values, upper.values, alpha=0.2)
        plt.xlabel("予測 Horizon（月）")
        plt.ylabel("1%の介入による目的変数の変化量")
        plt.title("感度分析：Horizonごとの影響と信頼区間（減衰付き）")
        plt.legend()
        plt.grid(True)
        plt.tight_layout()
        plt.show()

    return df_result

assistant:

以下は、sensitivity_analysis_with_ci_decay 関数をEWMA（指数加重移動平均）/平滑ラグ平均に対応させた改良版です。具体的には、ラグ変数の代わりに「EWMAで平滑化された変数」を用いる形に変更します。

🔧 主な変更点

  • decay_rate 引数が指定されている場合、各変数のラグを計算する代わりに、EWMAで平滑化された特徴量を作成。
  • ewma_span 引数を新設（decay_rateの代わりに span 指定で平滑化できるように）。
  • 元のラグ作成部分を削除して、代わりに pandas.Series.ewm().mean() を使って最新値から過去方向に指数平滑。

🆕 改変後の関数（EWMA版）

def sensitivity_analysis_with_ewma(
    df, date_col, target_var, explanatory_vars,
    intervention_scale=0.01, lags=12,
    horizons=range(1, 13), n_components=0.95,
    n_bootstrap=100, random_state=42,
    ewma_span=6, show_plot=True
):
    """
    Horizon別感度分析 + 95%信頼区間 + EWMAによる平滑ラグ平均

    ewma_span: EWMAに使うspan（小さいほど最近の値を重視）
    """

    np.random.seed(random_state)
    df = df.sort_values(date_col).copy()
    results_all = {}

    for horizon in horizons:
        df_feat = df.copy()

        # EWMA処理（各列に対して、過去の値を使って滑らかに）
        for var in explanatory_vars + [target_var]:
            df_feat[f'{var}_ewma'] = df_feat[var].ewm(span=ewma_span, adjust=False).mean()

        df_feat['target_future'] = df_feat[target_var].shift(-horizon)
        df_feat.dropna(inplace=True)

        feature_cols = [f'{var}_ewma' for var in explanatory_vars + [target_var]]
        X = df_feat[feature_cols].copy()
        y = df_feat['target_future'].copy()
        X_latest = X.iloc[[-1]].copy()

        baseline_preds = []
        intervened_preds_dict = {var: [] for var in explanatory_vars}

        for _ in range(n_bootstrap):
            X_resampled, y_resampled = resample(X, y)

            scaler = StandardScaler()
            X_scaled = scaler.fit_transform(X_resampled)
            pca = PCA(n_components=n_components)
            X_pca = pca.fit_transform(X_scaled)
            model = LinearRegression()
            model.fit(X_pca, y_resampled)

            X_latest_scaled = scaler.transform(X_latest)
            X_latest_pca = pca.transform(X_latest_scaled)
            y_pred_baseline = model.predict(X_latest_pca)[0]
            baseline_preds.append(y_pred_baseline)

            for var in explanatory_vars:
                X_modified = X_latest.copy()
                col = f'{var}_ewma'
                if col in X_modified.columns:
                    X_modified[col] *= (1 + intervention_scale)
                X_modified_scaled = scaler.transform(X_modified)
                X_modified_pca = pca.transform(X_modified_scaled)
                y_pred_intervened = model.predict(X_modified_pca)[0]
                intervened_preds_dict[var].append(y_pred_intervened)

        baseline_preds = np.array(baseline_preds)
        result = {}
        for var in explanatory_vars:
            diffs = np.array(intervened_preds_dict[var]) - baseline_preds
            result[var] = {
                'mean': np.mean(diffs),
                'lower': np.percentile(diffs, 2.5),
                'upper': np.percentile(diffs, 97.5)
            }

        results_all[horizon] = result

    # 結果の整形
    df_result = pd.DataFrame({
        (var, metric): [results_all[h][var][metric] for h in horizons]
        for var in explanatory_vars
        for metric in ['mean', 'lower', 'upper']
    }, index=horizons)
    df_result.index.name = "horizon"

    # グラフ表示
    if show_plot:
        plt.figure(figsize=(12, 6))
        for var in explanatory_vars:
            mean = df_result[(var, 'mean')]
            lower = df_result[(var, 'lower')]
            upper = df_result[(var, 'upper')]
            plt.plot(mean.index, mean.values, label=f'{var} (mean)')
            plt.fill_between(mean.index, lower.values, upper.values, alpha=0.2)
        plt.xlabel("予測 Horizon（月）")
        plt.ylabel("1%の介入による目的変数の変化量")
        plt.title("感度分析：Horizonごとの影響と信頼区間（EWMA平滑付き）")
        plt.legend()
        plt.grid(True)
        plt.tight_layout()
        plt.show()

    return df_result

✅ 補足

引数名 内容
ewma_span 小さいほど最近の値を重視（例：6で「6期間」程度の重み）
lags 形式上は残してありますが、EWMAでは使用していません。将来的に削除してもOKです
target_var もEWMA平滑する理由 target_var_lag* があった元コードに倣い、目的変数の履歴も特徴量に含める想定です

🔁 旧関数との使い分け

方式 特徴
sensitivity_analysis_with_ci_decay ラグ付き、かつ減衰係数で加重平均
sensitivity_analysis_with_ewma（新） ラグを使わずEWMAで平滑化（より滑らか・次元圧縮向き）

ご希望があれば、両方の方式を1つの関数で切り替え可能にする統合版もご提供できます！

はい、sensitivity_analysis_with_ci_decay と sensitivity_analysis_with_ewma はどちらも時系列の**感度分析（変数の微小な変化が将来の目的変数に与える影響）**を行うための関数ですが、ラグ変数の取り扱いと重みづけの方法が異なり、実用上の意味や用途も違ってきます。以下にその違いを詳しく整理します。

🔍 1. アプローチの違い：ラグの使い方

関数名 特徴 使用する情報
sensitivity_analysis_with_ci_decay 明示的にラグ変数を多数作成し、それぞれに減衰係数を掛ける 各ラグ（例：1ヶ月前, 2ヶ月前, ...）の変数
sensitivity_analysis_with_ewma ラグを使わず、指数加重移動平均（EWMA）で過去情報を1変数に圧縮 過去の値を重み付きで滑らかにした1変数

🧠 2. モデルの直感的な違い

◉ with_ci_decay: 「各ラグの効果を個別に考慮」

  • 各時点のラグがどれだけ目的変数に影響しているかを明示的に加味。
  • 例：1ヶ月前の気温は強く、6ヶ月前の気温は弱く効く → そういう構造が学習できる。

◉ with_ewma: 「過去の情報を滑らかに統合して1つの特徴量で表現」

  • 各変数の「過去からの影響の合計」を1変数で要約。
  • 例：気温が最近高めかどうかだけを捉える（時系列の過去形状は問わない）

📊 3. 特徴量の次元・モデルの安定性

観点 with_ci_decay with_ewma
特徴量の次元 高くなる（変数数 × ラグ数） 低い（変数数）
学習データ量が必要 多め（高次元なので） 少なめでも安定しやすい
過学習リスク 高め（特に小規模データでは） 低め（圧縮されてるため）
PCAの必要性 高い（多変量なので） 低い（次元が小さいため）

📈 4. 時系列のダイナミクス表現力

表現力 with_ci_decay with_ewma
短期 vs 長期効果の区別 明確に可能（ラグの番号で） 困難（1つの変数に圧縮されている）
反転・周期性の検出 可能（複数ラグを保持するため） 困難（滑らか化されて平均的に）
シグナルが鋭い変化を含むとき 有利（各ラグを保持） 不利（スムージングにより情報喪失）
ノイズが多いとき 不利（ノイズも含めて高次元化） 有利（ノイズ除去される）

🛠 5. 適用例・向いている用途

タスク with_ci_decay が向いている with_ewma が向いている
因果時系列モデルの分析 ✅ YES △ 部分的
各ラグの重要性評価 ✅ YES ❌ 不可（圧縮済）
小規模データで安定分析 ❌ 不安定 ✅ YES
ノイズ除去・トレンド抽出 △ ✅ YES
モデルの説明性 高い（ラグ単位で可視化可能） 低い（要約されている）

⚙ 補足：decay_rate vs ewma_span の違い

パラメータ 使われる関数 解釈
decay_rate with_ci_decay 指数的なラグ重み（例：0.9 → 1ヶ月前は1.0、2ヶ月前は0.9、3ヶ月前は0.81）
ewma_span with_ewma EWMAの平滑化幅（小さいほど最近の値を重視）

✅ 結論・どちらを使うべき？

状況 推奨関数
各ラグの影響を定量的に分析したい sensitivity_analysis_with_ci_decay
過去の傾向全体が目的変数に効いている（と思われる） sensitivity_analysis_with_ewma
データが少ない／次元削減したい with_ewma（低次元で安定）
特徴量の効果を明確に可視化したい with_ci_decay

必要であれば、両者を統一して method='decay' or 'ewma' のように切り替えられる関数にしても構いません。ご希望があれば統合版も提供します。

user:

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

assistant:

こちらこそありがとうございました！
本日の内容（EWMA・平滑ラグ平均・感度分析関数の改良）をふまえて、ブログタイトル案をいくつかご提案します。

📝 ブログタイトル案（技術寄り）

  1. 「時系列の感度分析を滑らかに：Pythonで実装するEWMA × 平滑ラグ平均」
  2. 「指数加重平均でラグを圧縮！感度分析関数の実践的アップグレード」
  3. 「ラグか？EWMAか？Pythonで比べる感度分析アプローチ」
  4. 「ブートストラップで信頼区間も：平滑ラグ平均による時系列介入分析」
  5. 「感度分析 × 時系列 × EWMA：次元圧縮と滑らかさのトレードオフに迫る」

🧑‍💻 読み物・学習寄りタイトル案

  1. 「Pythonで学ぶ：感度分析におけるEWMAとラグ平均の違いと実装」
  2. 「時系列の未来を予測するには？ラグより滑らかさ？EWMAで考える」
  3. 「減衰付きラグとEWMA、どちらが有効？Pythonで検証してみた」
  4. 「小さな介入が未来にどう効く？時系列の感度分析入門」
  5. 「感度分析の次の一手：ラグ情報をなめらかに扱うには？」

ご希望に合わせて絞り込んだり、学術寄り／初心者向け／実践Tips風に調整も可能です。お気に入りの方向性があればお知らせください！

user:

「時系列の感度分析を滑らかに：Pythonで実装するEWMA × 平滑ラグ平均」でいきたいと思います。タイトルに沿って、ここまでの内容を詳しくまとめてください。

assistant:

承知しました。以下に、ブログ記事「時系列の感度分析を滑らかに：Pythonで実装するEWMA × 平滑ラグ平均」の構成と内容を詳しくまとめました。

📝 時系列の感度分析を滑らかに：Pythonで実装するEWMA × 平滑ラグ平均

はじめに

時系列データの解析において、「過去の変数が将来にどの程度影響を与えるのか？」を調べる感度分析（sensitivity analysis）は、因果推論やシナリオシミュレーションにとって不可欠です。

本記事では、次の2つのアプローチをPythonで実装し、その違いや使い分け方を詳しく解説します：

  • 減衰ウェイト付きラグ変数（Decayラグ）
  • EWMA（指数加重移動平均）による平滑ラグ平均

感度分析の目的と前提

感度分析では、以下のような問いに答えることを目的とします：

「もし説明変数を1%だけ変化させたら、将来の目的変数はどれくらい変わるか？」

本記事では、以下のような前提に基づき感度分析を行います：

  • 入力は時系列のDataFrame
  • target_var：目的変数（例：売上、温度など）
  • explanatory_vars：説明変数（例：広告費、CO2濃度など）
  • Horizon：将来の何ヶ月後を予測するか（1〜12ヶ月など）

アプローチ①：減衰ウェイト付きラグ（Decay）

✔ コンセプト

各説明変数のラグ（過去の値）を複数生成し、過去になるほど減衰重みを掛けて特徴量とする方法です。

✔ 実装概要

for lag in range(1, lags + 1):
    decay_weight = decay_rate ** (lag - 1)
    df[f'{var}_lag{lag}'] = df[var].shift(lag) * decay_weight

このラグ変数を使って、ブートストラップ・PCA・線形回帰を用いて感度分析を行い、信頼区間付きで出力します。

✔ メリット・注意点

利点 課題
各ラグの影響を個別に分析できる 特徴量が増えて高次元になりやすい
時系列的な効果を精密に捉えられる 小データには過学習のリスクあり

アプローチ②：EWMA（指数加重移動平均）

✔ コンセプト

ラグを複数持たず、過去の情報をEWMAで滑らかに1変数にまとめるアプローチです。

$$
\text{EWMA}t = \alpha \cdot x_t + (1 - \alpha) \cdot \text{EWMA}{t-1}
$$

または pandas.Series.ewm(span=6, adjust=False).mean() を使います。

✔ 実装ポイント

df[f'{var}_ewma'] = df[var].ewm(span=ewma_span, adjust=False).mean()

これにより各変数の履歴情報を低次元かつ滑らかに特徴量化できます。

✔ メリット・注意点

利点 課題
データが少なくても安定 各ラグの効果を区別できない
ノイズ除去に優れる 鋭い変化への反応は鈍くなる

実装：感度分析関数（EWMA版）

以下のような関数を実装して、介入効果の平均と信頼区間を推定できます。

def sensitivity_analysis_with_ewma(...):
    # 入力：時系列DataFrame, 目的変数, 説明変数, span, horizonなど
    # 処理：EWMA平滑、PCA、回帰、ブートストラップ
    # 出力：Horizon別の感度と95%信頼区間
    ...

グラフ出力例：

  • 横軸：将来の予測 Horizon（月）
  • 縦軸：1%の介入による目的変数の平均変化量
  • エラーバー：95%信頼区間

比較まとめ：Decayラグ vs EWMA

観点 減衰ラグ方式 EWMA方式
各ラグの影響を見たい ✅ 可能 ❌ 不可
特徴量数 多くなりがち 少なく済む
学習安定性 要データ量 少ないデータでも可
説明性 高い（ラグを追跡可能） 低い（圧縮されている）
適用例 因果探索、構造把握 ノイズ除去、短期トレンドの要約

おわりに

時系列の感度分析では、過去情報の扱い方が分析の精度や安定性に大きく影響します。本記事で紹介した 減衰付きラグ平均 と EWMA平滑 のアプローチは、目的やデータ構造に応じて選択するのが理想的です。

今後はこのアプローチを拡張して、

  • 異なる変数間の非線形効果の検出
  • 多変量系列のクラスタリング前処理
  • 因果推論との統合分析

などにも応用していきたいと考えています。

🔗 ソースコード・リポジトリ（任意）

GitHubリンク or Colabノートブックなどがあれば紹介

ご希望であれば、Markdown形式や公開用のJupyter Notebook版に整形することも可能です。ご自由にカスタマイズしてお使いください！

user:

ありがとうございました！

assistant:

どういたしまして！こちらこそ、ありがとうございました😊
またいつでもお気軽にご相談ください。素敵な記事になりますように、応援しています！

よい一日を！

1
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
1
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?