user:
潜在因子を考えたくて、因子分析をしてみたんですが、よく考えると、取り扱うのが時系列データならば、ラグ変数を作成して因子分析をする必要があるんでしょうか?
assistant:
とても良い問いです。結論から言うと:
時系列データに因子分析を適用する場合、ラグ変数を導入するかどうかは、あなたが「どんな因子構造を仮定しているか」によって異なります。
以下で詳しく説明します。
🔹 通常の因子分析と時系列データ
因子分析(Factor Analysis)は、観測された多変量データを、少数の**潜在因子(共通因子)**と誤差に分解しようとするモデルです:
$$
\mathbf{x}_t = \mathbf{\Lambda} \mathbf{f}_t + \mathbf{\epsilon}_t
$$
- $\mathbf{x}_t$: 時点$t$での観測データ($p$次元)
- $\mathbf{f}_t$: 潜在因子($k$次元、通常 $k \ll p$)
- $\mathbf{\Lambda}$: 因子負荷行列
- $\mathbf{\epsilon}_t$: 特有因子(誤差)
🔹 時系列データにラグを取る必要があるケース
因子分析はもともと時点の独立性を前提としています。しかし時系列では自己相関があるため、単純な因子分析では時間的構造が無視されることになります。
このため、以下のような考慮が必要になります:
① 動的因子モデル(Dynamic Factor Model, DFM)
動的因子モデルでは、因子自体や観測値が過去の因子や観測値の影響を受けるというモデルです:
$$
\begin{align}
\mathbf{x}_t &= \mathbf{\Lambda} \mathbf{f}_t + \mathbf{\epsilon}_t \
\mathbf{f}t &= \sum{l=1}^L \mathbf{\Phi}l \mathbf{f}{t-l} + \mathbf{u}_t
\end{align}
$$
- 観測値は因子に依存するが、因子自体も自己回帰過程に従う。
- こうしたモデルでは、ラグを含めた因子分析が有効になります。
② ラグ付き因子分析(拡張アプローチ)
- 観測変数にラグを取ったもの(例:$x_{t}, x_{t-1}, \dots$)をすべて使って因子分析をする。
- データ行列の列数が増える(ラグ分だけ次元拡張)ことで、時間的情報を静的因子モデルに押し込むことができます。
# pandasでラグ変数を作る例(ラグ3まで)
lagged_df = pd.concat([df.shift(i) for i in range(4)], axis=1)
🔹 まとめ:ラグを使うかどうかの判断基準
| 条件 | ラグ変数を導入すべきか? | 説明 |
|---|---|---|
| 時系列の構造を無視して因子を取りたい場合 | ❌ 必要ない | 単純な因子分析で十分(例:月次データの共通要因) |
| 時系列の自己相関や動的構造も説明したい場合 | ✅ 必要あり | ラグ付き因子分析や動的因子モデルが有効 |
🔸 実務的アドバイス
- 観測変数が数個〜十数個程度なら:まずはラグなしで因子分析 → モデルがうまく説明できなければラグを試す。
- 高頻度 or 高次元データ(例:株式、センサーデータなど)なら:最初から動的因子モデルを検討してよいです。
必要であれば、Pythonでの「ラグ付き因子分析」や「動的因子モデル(DFM)」の実装もお手伝いできます!ご希望があればお知らせください。
user:
ありがとうございます。それでは、次の要件を満たす Python 関数を書いてください。
入力:
・ df : pandas dataframe 形式データ
・ explanatory_vars : 説明変数を表す列名のリスト
手順:
・ラグ付き因子分析を行う。
出力:
・結果を、pandas dataframe 形式で返す。
assistant:
以下は、ラグ付き因子分析を行うためのPython関数の実装です。
この関数では、指定された説明変数についてラグ変数を作成し、それらに対して**因子分析(Factor Analysis)を適用し、各サンプル時点における因子スコア(潜在因子)**を DataFrame として返します。
✅ 実装内容の要点
- ラグ数は
max_lagで指定可能(デフォルトは 3) - 欠損値処理のために
dropna()を実施 - 標準化(StandardScaler)を事前に実施
-
sklearn.decomposition.FactorAnalysisを使用 - 結果の出力は
pandas.DataFrame(各時点の因子スコア)
🔧 コード:lagged_factor_analysis
import pandas as pd
import numpy as np
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import FactorAnalysis
def lagged_factor_analysis(df, explanatory_vars, n_factors=1, max_lag=3):
"""
ラグ付き因子分析を行う関数
Parameters:
df (pd.DataFrame): 入力データフレーム
explanatory_vars (list): 説明変数の列名リスト
n_factors (int): 抽出する因子の数(デフォルト1)
max_lag (int): 使用する最大ラグ数(デフォルト3)
Returns:
pd.DataFrame: 各時点における因子スコア
"""
df = df.copy()
# ラグ変数の作成
lagged_vars = []
for lag in range(max_lag + 1):
lagged = df[explanatory_vars].shift(lag).add_suffix(f"_lag{lag}")
lagged_vars.append(lagged)
X_lagged = pd.concat(lagged_vars, axis=1).dropna()
# 標準化
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X_lagged)
# 因子分析の適用
fa = FactorAnalysis(n_components=n_factors)
factors = fa.fit_transform(X_scaled)
# 因子スコアのDataFrame化
factor_df = pd.DataFrame(factors,
index=X_lagged.index,
columns=[f"factor_{i+1}" for i in range(n_factors)])
return factor_df
🧪 使用例
# 説明変数
explanatory_vars = ["x1", "x2", "x3"]
# 因子スコアの抽出
factors = lagged_factor_analysis(df, explanatory_vars, n_factors=2, max_lag=3)
print(factors.head())
✅ 出力のイメージ
| factor_1 | factor_2 | |
|---|---|---|
| 2020-01-04 | 0.13 | -0.42 |
| 2020-01-05 | 0.18 | -0.39 |
| ... | ... | ... |
必要であれば、因子負荷行列や固有値、可視化なども追加可能です。さらに「動的因子モデル(状態空間モデルベース)」に拡張することもできますので、ご希望があれば教えてください。
user:
ありがとうございます。テスト用のデータを生成する関数も書いてください。
assistant:
こちらこそありがとうございます。
以下に、ラグ付き因子分析のテストに使える疑似データを生成する関数を提供します。
🔧 関数の設計方針
- 潜在因子(latent factor)を1つまたは複数生成
- 各観測変数は、この因子+ノイズから構成
- 必要に応じて、時系列構造(自己相関)も加味
-
pandas.DataFrameとして出力(時系列インデックス付き)
✅ コード:generate_synthetic_factor_data
import numpy as np
import pandas as pd
def generate_synthetic_factor_data(n_obs=200, n_vars=5, n_factors=1, seed=42):
"""
潜在因子に基づく時系列データを生成する関数
Parameters:
n_obs (int): 時系列の長さ(サンプル数)
n_vars (int): 観測変数の数
n_factors (int): 潜在因子の数
seed (int): 乱数シード
Returns:
pd.DataFrame: 説明変数からなる時系列データ
"""
np.random.seed(seed)
# 潜在因子(AR構造を持たせる)
factors = np.zeros((n_obs, n_factors))
for i in range(n_factors):
noise = np.random.normal(0, 1, size=n_obs)
for t in range(1, n_obs):
factors[t, i] = 0.8 * factors[t-1, i] + noise[t]
# 因子負荷行列(n_vars × n_factors)
loadings = np.random.uniform(0.5, 1.0, size=(n_vars, n_factors))
# 観測変数の生成: x = Λf + ε
epsilon = np.random.normal(0, 0.5, size=(n_obs, n_vars))
X = np.dot(factors, loadings.T) + epsilon
# DataFrame化
dates = pd.date_range(start="2020-01-01", periods=n_obs, freq="D")
columns = [f"x{i+1}" for i in range(n_vars)]
df = pd.DataFrame(X, index=dates, columns=columns)
return df
🧪 使用例
# データ生成
df_test = generate_synthetic_factor_data(n_obs=300, n_vars=6, n_factors=2)
# ラグ付き因子分析を適用
factor_scores = lagged_factor_analysis(df_test, explanatory_vars=df_test.columns.tolist(), n_factors=2, max_lag=3)
print(factor_scores.head())
📈 補足:データの可視化(例)
import matplotlib.pyplot as plt
df_test.plot(subplots=True, figsize=(10, 8), title="Synthetic Time Series")
plt.tight_layout()
plt.show()
必要であれば、潜在因子や因子負荷構造を可視化したり、系列の相関や構造も調整可能です。希望があれば教えてください。
user:
ありがとうございます。ここで出力される因子スコアの構造を見た限り、ラグ変数が考慮されているか分からない構造になっているのですが、されているんですか?
assistant:
鋭いご指摘、ありがとうございます。
結論から言うと、はい、ラグ変数は考慮されています。ただし、出力された因子スコア自体には「ラグの情報が明示的に見えない構造」になっています。
以下で詳しく解説します。
✅ なぜラグが「考慮されている」が「見えない」のか?
🔸 入力変数にはラグを含めています
# 例えば、explanatory_vars = ['x1', 'x2']
# max_lag = 3 のとき、次のような列が作成されます:
['x1_lag0', 'x2_lag0',
'x1_lag1', 'x2_lag1',
'x1_lag2', 'x2_lag2',
'x1_lag3', 'x2_lag3']
これを1つの大きな行列にして、因子分析(sklearnの FactorAnalysis)に渡しています。
つまり、元の観測変数 + ラグ変数が全部同じ因子分析モデルに投入されています。
🔸 しかし、出力される因子スコアは次元削減後の成分
sklearnの FactorAnalysis は、与えられた特徴行列(今回はラグ付き特徴)から、以下を求めます:
- 潜在因子(スコア):
f_tの推定 - 因子負荷行列($\Lambda$):あとで
.components_で取得可能
これにより、出力は「抽出された潜在因子」の時系列スコアになり、入力のどのラグがどの程度貢献したかの情報は見えなくなります。
🔍「ラグがどのように効いているか」を知る方法
出力された因子スコアだけではラグの影響はブラックボックス的になりますが、以下の方法で内部構造を確認できます:
方法①:因子負荷行列(loadings)を見てみる
loadings = fa.components_.T # shape: [n_features, n_factors]
columns = [f"factor_{i+1}" for i in range(n_factors)]
index = X_lagged.columns # ラグ付き変数の名前
loading_df = pd.DataFrame(loadings, index=index, columns=columns)
print(loading_df)
これを見ると、
- どの「ラグ付き変数」が
- どの「潜在因子」に
- どれくらい寄与しているか
が分かります。
方法②:可視化でラグの影響を見る
import seaborn as sns
import matplotlib.pyplot as plt
plt.figure(figsize=(10, 6))
sns.heatmap(loading_df, annot=True, cmap="coolwarm", center=0)
plt.title("因子負荷行列(ラグ付き変数 vs 因子)")
plt.tight_layout()
plt.show()
✅ まとめ
| 観点 | 説明 |
|---|---|
| ✅ ラグは考慮されているか? | はい。ラグ変数を特徴量に追加し、因子分析しています |
| ⚠ 出力の因子スコアにラグ情報はあるか? | 暗に含まれていますが、見えません |
| 🔍 どのラグが影響しているかを知るには? |
.components_ による因子負荷行列を確認する |
必要であれば、「因子スコア vs ラグ変数の寄与関係」をわかりやすく図示するアプローチや、因子負荷行列の自動解釈ロジックも作れます。興味があれば教えてください!
承知しました!
以下では、先ほどの lagged_factor_analysis 関数を拡張して、次のように機能を追加します:
✅ 拡張内容
- 因子スコア(各時点の潜在因子) に加えて
- 因子負荷行列(ラグ付き変数 vs 因子) も返す
- 因子負荷行列を ヒートマップで可視化
🔧 改良版:lagged_factor_analysis_with_loadings
import pandas as pd
import numpy as np
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import FactorAnalysis
import matplotlib.pyplot as plt
import seaborn as sns
def lagged_factor_analysis_with_loadings(df, explanatory_vars, n_factors=1, max_lag=3, plot_loadings=True):
"""
ラグ付き因子分析を行い、因子スコアと因子負荷行列を返す関数
Parameters:
df (pd.DataFrame): 入力データ
explanatory_vars (list): 説明変数の列名リスト
n_factors (int): 抽出する因子数(デフォルト1)
max_lag (int): ラグの最大数(デフォルト3)
plot_loadings (bool): 負荷行列をヒートマップで描画するか(デフォルトTrue)
Returns:
factor_scores (pd.DataFrame): 各時点の因子スコア
loading_df (pd.DataFrame): 因子負荷行列(ラグ付き変数 × 因子)
"""
df = df.copy()
# --- ラグ付き変数の作成 ---
lagged_vars = []
for lag in range(max_lag + 1):
lagged = df[explanatory_vars].shift(lag).add_suffix(f"_lag{lag}")
lagged_vars.append(lagged)
X_lagged = pd.concat(lagged_vars, axis=1).dropna()
# --- 標準化 ---
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X_lagged)
# --- 因子分析の適用 ---
fa = FactorAnalysis(n_components=n_factors)
factors = fa.fit_transform(X_scaled)
# --- 因子スコアを DataFrame に変換 ---
factor_scores = pd.DataFrame(factors,
index=X_lagged.index,
columns=[f"factor_{i+1}" for i in range(n_factors)])
# --- 因子負荷行列の作成 ---
loadings = fa.components_.T # shape: [n_features, n_factors]
loading_df = pd.DataFrame(loadings,
index=X_lagged.columns,
columns=[f"factor_{i+1}" for i in range(n_factors)])
# --- 可視化(ヒートマップ) ---
if plot_loadings:
plt.figure(figsize=(1.5 * n_factors + 4, 0.5 * len(loading_df) + 2))
sns.heatmap(loading_df, annot=True, fmt=".2f", cmap="coolwarm", center=0)
plt.title("因子負荷行列(ラグ付き変数 × 因子)")
plt.tight_layout()
plt.show()
return factor_scores, loading_df
🧪 使用例
# テストデータを作成
df_test = generate_synthetic_factor_data(n_obs=200, n_vars=5, n_factors=2)
# ラグ付き因子分析 + 因子負荷行列 + 可視化
factor_scores, loading_df = lagged_factor_analysis_with_loadings(
df_test, explanatory_vars=df_test.columns.tolist(),
n_factors=2, max_lag=3, plot_loadings=True
)
✅ 出力形式
factor_scores の例:
| 日付 | factor_1 | factor_2 |
|---|---|---|
| 2020-01-04 | -0.221 | 0.073 |
| 2020-01-05 | -0.315 | 0.034 |
| ... | ... | ... |
loading_df の例:
| 変数名 | factor_1 | factor_2 |
|---|---|---|
| x1_lag0 | 0.68 | -0.12 |
| x2_lag1 | 0.55 | 0.19 |
| ... | ... | ... |
必要に応じて、因子スコアの時系列プロット、因子ごとの寄与割合(固有値ベース)なども追加可能です。
ご希望があればお気軽にどうぞ!
user:
ありがとうございます。ところで、適切な因子の数を決めるための関数を作っていただけませんか?
assistant:
もちろんです!
因子の数(次元)を決めることは因子分析の重要なステップです。
以下では、情報量規準(AIC/BIC)や累積寄与率(固有値)に基づいて、最適な因子数を推定するPython関数をご紹介します。
✅ 提供する関数の特徴
- ラグ付き変数を含めて判断可能
- AIC, BIC(sklearnの
loglikeに基づく推定式) - 固有値に基づく「スクリープロット」「累積寄与率プロット」も作成可能
- 結果を
DataFrameとして返す
🔧 関数:select_optimal_n_factors
from sklearn.decomposition import FactorAnalysis
from sklearn.preprocessing import StandardScaler
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
def select_optimal_n_factors(df, explanatory_vars, max_factors=10, max_lag=0, plot=True):
"""
AIC/BIC/累積寄与率に基づいて因子数を選ぶ補助関数
Parameters:
df (pd.DataFrame): 入力データ
explanatory_vars (list): 説明変数名のリスト
max_factors (int): 最大因子数
max_lag (int): 最大ラグ数(0ならラグなし)
plot (bool): 結果をグラフ表示するか
Returns:
pd.DataFrame: 因子数ごとの AIC / BIC / 再構成誤差
"""
# --- ラグ付きデータの作成 ---
lagged_vars = []
for lag in range(max_lag + 1):
lagged = df[explanatory_vars].shift(lag).add_suffix(f"_lag{lag}")
lagged_vars.append(lagged)
X = pd.concat(lagged_vars, axis=1).dropna()
X_std = StandardScaler().fit_transform(X)
n_obs, n_features = X_std.shape
results = []
for n_factors in range(1, min(max_factors, n_features) + 1):
fa = FactorAnalysis(n_components=n_factors)
X_trans = fa.fit_transform(X_std)
X_recon = np.dot(X_trans, fa.components_) + fa.mean_
resid = X_std - X_recon
resid_var = np.var(resid)
# log-likelihood 近似(モデル自由度補正なし)
ll = -0.5 * n_obs * (n_features * np.log(2 * np.pi * resid_var) + n_features)
# パラメータ数 = loadings + unique variances
n_params = n_factors * n_features + n_features
aic = -2 * ll + 2 * n_params
bic = -2 * ll + n_params * np.log(n_obs)
results.append({
"n_factors": n_factors,
"AIC": aic,
"BIC": bic,
"resid_var": resid_var
})
result_df = pd.DataFrame(results).set_index("n_factors")
if plot:
result_df[["AIC", "BIC"]].plot(marker='o', title="AIC/BICによる因子数選択")
plt.xlabel("Number of Factors")
plt.grid(True)
plt.show()
return result_df
🧪 使用例
# テストデータ
df_test = generate_synthetic_factor_data(n_obs=200, n_vars=6, n_factors=2)
# 因子数選定(ラグ2まで含める)
result_df = select_optimal_n_factors(df_test, explanatory_vars=df_test.columns.tolist(), max_factors=6, max_lag=2)
print(result_df)
✅ 出力のイメージ
| n_factors | AIC | BIC | resid_var |
|---|---|---|---|
| 1 | 4221.9 | 4305.3 | 0.831 |
| 2 | 3921.5 | 4045.2 | 0.602 |
| 3 | 3918.4 | 4082.4 | 0.547 |
| ... | ... | ... | ... |
- AIC/BICが最小になる因子数が候補
-
resid_varが極端に変わらなくなる箇所も「ひざ点」として注目可能です
✅ 拡張案(ご希望あれば)
- 固有値のスクリープロット/累積寄与率を別関数で提供
- バリマックス回転付き因子分析
- クロスバリデーションでの復元誤差評価
必要に応じて、これらを統合して因子分析ワークフロー全体を整備することもできます。興味があればお知らせください!
user:
固有値を使って因子数を決める方法はできませんか?
assistant:
はい、固有値に基づく因子数の決定も非常に一般的で有力な方法です。
以下のようなアプローチが主に使われます:
✅ 固有値に基づく因子数決定法
| 方法 | 内容 |
|---|---|
| ① 固有値1超基準(Kaiserの基準) | 固有値 > 1 の主成分を因子として採用する(標準化前提) |
| ② スクリープロット(Scree Plot) | 固有値の減少カーブの「ひざ(エルボー)」で因子数を決定 |
| ③ 累積寄与率基準 | 例:累積寄与率が80%以上となる最小の因子数を採用 |
これらは通常、**主成分分析(PCA)**の固有値を使って「何因子でどれだけ分散が説明されるか」を見る形で行います。
🔧 固有値に基づく因子数判定関数
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
def select_n_factors_by_eigen(df, explanatory_vars, max_lag=0, threshold=1.0, cumvar_threshold=0.8, plot=True):
"""
固有値に基づいて因子数を決定する関数
Parameters:
df (pd.DataFrame): 入力データ
explanatory_vars (list): 説明変数の列名リスト
max_lag (int): 最大ラグ数(デフォルト0=ラグなし)
threshold (float): 固有値のしきい値(Kaiser基準:通常1.0)
cumvar_threshold (float): 累積寄与率の閾値(例:0.8)
plot (bool): スクリープロット表示
Returns:
dict: 以下の情報を含む辞書
- 'eigenvalues': 固有値一覧(np.ndarray)
- 'n_factors_kaiser': 固有値1超の数
- 'n_factors_cumvar': 累積寄与率で満たされる最小因子数
"""
# ラグ付きデータの作成
lagged_vars = []
for lag in range(max_lag + 1):
lagged = df[explanatory_vars].shift(lag).add_suffix(f"_lag{lag}")
lagged_vars.append(lagged)
X = pd.concat(lagged_vars, axis=1).dropna()
X_std = StandardScaler().fit_transform(X)
# PCAによる固有値の取得
pca = PCA()
pca.fit(X_std)
eigenvalues = pca.explained_variance_
cumvar = np.cumsum(pca.explained_variance_ratio_)
# 基準に基づく因子数
n_factors_kaiser = np.sum(eigenvalues > threshold)
n_factors_cumvar = np.argmax(cumvar >= cumvar_threshold) + 1
if plot:
plt.figure(figsize=(8, 5))
plt.plot(np.arange(1, len(eigenvalues) + 1), eigenvalues, marker='o', label='Eigenvalues')
plt.axhline(y=threshold, color='r', linestyle='--', label=f'λ > {threshold}')
plt.axvline(x=n_factors_cumvar, color='g', linestyle='--', label=f'Cumulative > {cumvar_threshold*100:.0f}%')
plt.xlabel("Component Number")
plt.ylabel("Eigenvalue")
plt.title("スクリープロット(固有値による因子数判定)")
plt.legend()
plt.grid(True)
plt.tight_layout()
plt.show()
return {
"eigenvalues": eigenvalues,
"n_factors_kaiser": n_factors_kaiser,
"n_factors_cumvar": n_factors_cumvar
}
🧪 使用例
# テストデータ
df_test = generate_synthetic_factor_data(n_obs=200, n_vars=6, n_factors=2)
# 固有値による因子数推定
result = select_n_factors_by_eigen(df_test, df_test.columns.tolist(), max_lag=2)
print(f"固有値 > 1 の因子数: {result['n_factors_kaiser']}")
print(f"累積寄与率80%以上に必要な因子数: {result['n_factors_cumvar']}")
✅ 出力例(dict)
{
'eigenvalues': array([...]),
'n_factors_kaiser': 3,
'n_factors_cumvar': 2
}
🔍 おすすめの実務フロー
- この関数で因子の数を候補として決める
- AIC/BICでもう一度確認
- 最終的にラグ付き因子分析を実行
ご希望があれば、AIC/BIC・固有値・累積寄与率を統合的に表示するダッシュボード的な関数も作成可能です。
いつでもどうぞ!
user:
ありがとうございました。ここまでの内容をブログにまとめたいので、タイトル案をいくつか提案してください。
assistant:
もちろんです!以下は、時系列データにラグを付加して因子分析を行うプロセスと、**因子数の決定方法(固有値やスクリープロット)**を含んだ内容にふさわしいブログ記事タイトル案です:
🔍 解説系タイトル案(学術・技術寄り)
- 「時系列データに因子分析を適用する:ラグ変数と因子数の選び方」
- 「ラグ付き因子分析のすすめ:時系列構造を活かす多変量解析」
- 「因子分析 × 時系列:ラグ付加と固有値による因子数の選定法」
- 「時系列因子分析入門:遅行変数と固有値から見る潜在構造」
- 「ラグ付き因子分析で時系列の共通因子を探る:Python実装付き」
🧪 実践・コード重視スタイル
- 「Pythonで実践するラグ付き因子分析:時系列データの潜在構造を解明」
- 「因子数は固有値で決める?Pythonで学ぶ時系列因子分析の手順」
- 「ラグを考慮した因子分析の実装と因子数決定法(Python付)」
- 「スクリープロットで見る時系列因子の構造:実装と分析ガイド」
🧠 読み物・問いかけスタイル
- 「時系列データに因子分析ってアリ?潜在因子を追うラグ付きアプローチ」
- 「共通因子は未来を語るか?ラグ付き因子分析で探る構造」
気に入ったスタイルがあれば、それに沿ってさらにアレンジも可能です!
user:
「時系列データに因子分析を適用する:ラグ変数と因子数の選び方」でいきたいと思います。このタイトルに沿って、今回の内容を詳しくまとめてください。
assistant:
以下に、ブログ記事「時系列データに因子分析を適用する:ラグ変数と因子数の選び方」にふさわしい構成と内容を、詳しくまとめました:
時系列データに因子分析を適用する:ラグ変数と因子数の選び方
はじめに
因子分析(Factor Analysis)は、多変量データの背後にある共通因子(潜在因子)を抽出する手法として知られています。しかし、もし取り扱っているデータが時系列だったらどうでしょうか?
本記事では、時系列データに因子分析を適用する際のポイントとして、
- ラグ(遅行)変数の導入
- 適切な因子数の選び方(固有値やスクリープロットの利用)
に焦点を当て、Pythonによる実装付きで詳しく解説します。
因子分析とは?
因子分析は、観測変数を少数の潜在因子で表現するモデルです。通常の因子分析は時系列構造を考慮しませんが、時系列に内在する動きや遅延効果を捉えたい場合には、ラグを取って構造を拡張する必要があります。
ラグ変数を加える理由
時系列データでは、現在の値が過去の値に依存することが多く見られます。例えば、経済指標、気象データ、株価などは、**時間遅れの影響(ラグ効果)**を反映する必要があります。
そのため、因子分析にラグ変数を加えることで、時系列的な共通因子の構造を捉えることができるようになります。
Pythonによるラグ付き因子分析の実装
ステップ1:ラグ付き因子分析関数の定義
import pandas as pd
import numpy as np
from sklearn.decomposition import FactorAnalysis
def lagged_factor_analysis(df, explanatory_vars, n_lags=1, n_factors=2):
"""
ラグ変数を含めた因子分析を実行する。
Parameters:
- df (pd.DataFrame): 入力の時系列データ
- explanatory_vars (list): 説明変数の列名リスト
- n_lags (int): ラグの数(何期前まで使うか)
- n_factors (int): 抽出する因子の数
Returns:
- pd.DataFrame: 因子スコア(各時点ごとの因子値)
"""
# ラグ付きデータの作成
lagged_data = []
for lag in range(n_lags + 1):
lagged = df[explanatory_vars].shift(lag)
lagged.columns = [f"{col}_lag{lag}" for col in explanatory_vars]
lagged_data.append(lagged)
combined = pd.concat(lagged_data, axis=1).dropna()
# 因子分析の実行
fa = FactorAnalysis(n_components=n_factors)
scores = fa.fit_transform(combined)
# 因子スコアの DataFrame 化
score_df = pd.DataFrame(scores, index=combined.index,
columns=[f"factor_{i+1}" for i in range(n_factors)])
return score_df
ステップ2:テストデータの作成
def generate_test_data(n_samples=200, seed=42):
"""
テスト用の時系列データを作成する関数。
潜在因子に基づく構造を持った系列を生成。
Returns:
- pd.DataFrame: テスト用の時系列データ
"""
np.random.seed(seed)
t = np.arange(n_samples)
common_factor = np.sin(2 * np.pi * t / 50) + np.random.normal(0, 0.2, size=n_samples)
df = pd.DataFrame({
'X1': 0.8 * common_factor + np.random.normal(0, 1, n_samples),
'X2': 0.6 * common_factor + np.random.normal(0, 1, n_samples),
'X3': 0.4 * common_factor + np.random.normal(0, 1, n_samples),
})
return df
因子数の決め方:固有値とスクリープロット
因子分析では、因子の数(次元数)をどう決めるかが非常に重要です。
よく使われる基準は次の2つです:
1. 固有値基準(Kaiser基準)
- 固有値が1以上の成分のみを因子とする
- 直感的で使いやすいが、過剰な因子数を選びやすい
2. スクリープロット
- 固有値のグラフを描き、「肘のように折れ曲がる点」までを因子とする
- 視覚的に適切な因子数を判断できる
実装例(PCAベース)
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA
def plot_scree_and_get_eigenvalues(df, explanatory_vars, n_lags=1):
"""
ラグ付きデータに対してスクリープロットを作成し、固有値を返す。
Returns:
- np.ndarray: 固有値の配列
"""
# ラグ付きデータの構築
lagged_data = []
for lag in range(n_lags + 1):
lagged = df[explanatory_vars].shift(lag)
lagged.columns = [f"{col}_lag{lag}" for col in explanatory_vars]
lagged_data.append(lagged)
combined = pd.concat(lagged_data, axis=1).dropna()
# PCA による固有値の計算
pca = PCA()
pca.fit(combined)
eigenvalues = pca.explained_variance_
# プロット
plt.figure(figsize=(8, 4))
plt.plot(range(1, len(eigenvalues) + 1), eigenvalues, marker='o')
plt.axhline(1.0, color='r', linestyle='--', label='Eigenvalue = 1')
plt.title("Scree Plot (固有値プロット)")
plt.xlabel("成分番号")
plt.ylabel("固有値")
plt.legend()
plt.grid(True)
plt.tight_layout()
plt.show()
return eigenvalues
おわりに
本記事では、時系列データに因子分析を適用する際に重要な、
- ラグ変数の導入(時間的依存性を考慮)
- 適切な因子数の選定(固有値・スクリープロット)
について、実装コード付きで解説しました。
時系列構造を意識した因子分析は、共通トレンドや潜在的なダイナミクスを理解するための強力なアプローチです。ぜひ自身のデータでも試してみてください。
ご希望があれば、Streamlit などで可視化できるアプリ化や、より動的な因子モデル(動学的因子モデルなど)への展開もご案内できます!
user:
ありがとうございました!
assistant:
どういたしまして!いつでもお手伝いしますので、ブログ執筆や分析でまた質問があれば遠慮なくどうぞ。
素敵な記事になりますように。執筆、応援しています!📚✨