LiNGAMで推定した因果グラフから2値データの変数への因果効果を計算する
はじめに
本記事では、LiNGAM(Linear Non-Gaussian Acyclic Model)を用いて推定した因果グラフから、2値データの変数への因果効果を計算する方法について解説します。
課題🤔
LiNGAMは主に連続値を対象とした因果探索手法で、データから因果関係を推定する強力なツールです。しかし、2値データの変数に対しては因果探索が難しいという課題があります。一方、実際の要因分析の現場では、正常・異常や良い・悪いといった2値データで表現できるカテゴリカル変数を因果の結果の変数として扱うケースも多く見られます。
この課題を克服するために、まず2値データの変数を除いたデータセットに対して、LiNGAMで因果グラフを推定します。その後、推定した因果グラフを基に、ロジスティック回帰を用いて各変数から2値データの変数への因果効果を計算します。
このアプローチにより、2値データの変数に対する要因分析が可能となります。
想定読者🧐
- Pythonでlingamパッケージを使って因果探索をしている人
- 2値データの変数を目的変数にして要因分析したい人
動作環境💻
- OS: Windows11
- Python: 3.10.9
- lingam: 1.9.0
- scikit-learn: 1.5.0
動作確認した日📅
- 2024年12月23日
2値データの変数への因果効果を計算する
使用したデータ📂
UCIのワインの品質データセットのうち白ワインのデータセットを使いました。🍷
quality(品質)変数を因果の結果の変数と仮定して分析します。
さらに、今回の手法の説明のためにquality変数を2値データに変換して分析します。
データの準備🐍
- 白ワインのデータセットをpandasデータフレームに読み込みます
- quality列の中央値を取得します
- quality列のデータを中央値以上で1(良い)、未満で0(悪い)に変換します
import pandas as pd
# データセットファイルを読み込み、qualityを中央値で2値に変換する
df = pd.read_csv("winequality-white.csv", encoding="UTF_8", sep=";")
median_quality = df["quality"].median()
df["quality"] = df["quality"].apply(lambda x:1 if x>=median_quality else 0)
2値データの変数を除いたデータセットでDirectLiNGAMを実行🐍
- quality列を削除したデータフレームを作成します
- DirectLiNGAMモデルを作成して、quality列を削除したデータフレームにフィットします
import lingam
# quality列を除いて、DirectLiNGAMを実行
_df = df.drop(columns=["quality"])
model = lingam.DirectLiNGAM()
model.fit(_df)
推定した因果グラフ
推定した因果グラフを基にロジスティック回帰を実行🐍
- quality列を削除したデータフレームの各列(変数)でループします
- 推定した因果モデルの隣接行列から、現在の変数の親の変数を取得します
- ロジスティック回帰モデルを作成し、quality変数にフィットさせます
- 現在の変数の値を操作し、介入後の予測確率を計算します
- 予測確率の差を因果効果として記録します
import numpy as np
from sklearn.linear_model import LogisticRegression
# 結果を格納するデータフレームを初期化
target_effects = pd.DataFrame(0, index=["quality"], columns=_df.columns.tolist())
# quality変数を除いたデータフレームの各列(変数)についてループ
for current_col in _df.columns:
# 元のデータフレームの列インデックスを取得
current_index = df.columns.get_loc(current_col)
# 隣接行列の該当変数の行で0以外のところから親の変数リストを取得
parents_indices = np.where(np.abs(model._adjacency_matrix[_df.columns.get_loc(current_col)]) > 0.0)[0]
# 隣接行列のインデックスをデータセットの列インデックスに変換
parents_indices = [df.columns.get_loc(df.columns[idx]) for idx in parents_indices]
# 説明変数に自分自身と親を追加
predictors = [current_index]
predictors.extend(parents_indices)
# ロジスティック回帰モデルを作成し、quality変数にフィットさせる
lr = LogisticRegression(solver='liblinear')
lr.fit(df.iloc[:, predictors], df.iloc[:, df.columns.get_loc("quality")])
# データフレームのコピーを作成し、現在の変数の値を平均値に設定
X_intervened = df.copy()
X_intervened.iloc[:, current_index] = df.iloc[:, current_index].mean() # do(x=E(x))
# 介入後の予測確率を計算
p1 = lr.predict_proba(X_intervened.iloc[:, predictors])
# 現在の変数の値を平均値+1に設定し、再度予測確率を計算
X_intervened.iloc[:, current_index] = df.iloc[:, current_index].mean() + 1 # do(x=E(x)+1)
p2 = lr.predict_proba(X_intervened.iloc[:, predictors])
# 二つの予測確率の平均の差を計算し、因果効果として記録
effect = p2[:, 1].mean() - p1[:, 1].mean()
target_effects.loc["quality", current_col] = effect
# 結果を表示用にソートする
target_effects = target_effects.T.sort_values(by="quality")
考察✏️
UCIの白ワインの品質データを用いて各変数がquality(品質)に与える因果効果を評価しました。
結果は以下の通りです:
quality変数への因果効果
変数 | 日本語訳 | 因果効果 |
---|---|---|
sulphates | 硫酸塩 | 0.168202 |
alcohol | アルコール度数 | 0.144802 |
pH | pH | 0.134204 |
citric acid | クエン酸 | 0.007925 |
free sulfur dioxide | 自由二酸化硫黄 | 0.003439 |
total sulfur dioxide | 総二酸化硫黄 | -0.000626 |
residual sugar | 残留糖 | -0.008489 |
fixed acidity | 固定酸度 | -0.042406 |
chlorides | 塩化物 | -0.263641 |
density | 密度 | -0.557571 |
volatile acidity | 揮発性酸度 | -0.648104 |
この結果から、以下の考察をしました:
- 正の因果効果:
硫酸塩(sulphates)、アルコール度数(alcohol)、およびpHは品質に対して正の影響を持つことが示されました。これらの変数は品質向上に寄与していると推測されます。 - 負の因果効果:
揮発性酸度(volatile acidity)と塩化物(chlorides)は品質に対して負の影響を持つことが示されました。これらの変数は品質低下の要因であると推測されます。 - その他の変数:
クエン酸(citric acid)、自由二酸化硫黄(free sulfur dioxide)、および総二酸化硫黄(total sulfur dioxide)は品質に対する影響が小さいことが示されました。これらの変数は品質に対して大きな影響を与えない可能性があります。
一般的に、高い揮発性酸度と過剰な塩化物は品質低下の要因とされており[1]、今回の結果もこれを裏付けるものとなっています。これらの知見を基に、品質向上のための具体的な対策を検討することが重要です。
おわりに
本記事では、LiNGAMを用いて推定した因果グラフから、2値データの変数への因果効果を計算する方法を紹介しました。この手法を用いることで、2値データの変数が因果の結果の変数である場合に、ロジスティック回帰を用いて因果効果を計算し、要因分析を行うことができます。
PR📺
弊社のCausalasでは、LiNGAMを用いた因果グラフの推定がノーコードで実現できます。ただし、2値データの扱いは今後のバージョンアップにご期待ください。
Discussion