Mediation analysisをStataで実行する | 医療統計とStataプログラミングの部屋

このブログでは、統計解析ソフトStataのプログラミングのTipsや便利コマンドを紹介しています．

Facebook groupでは、ちょっとした疑問や気づいたことなどを共有して貰うフォーラムになっています．ブログと合わせて個人の学習に役立てて貰えれば幸いです．

さて、今回は「媒介効果分析」について、それをStataで実行するためのコマンドを紹介したいと思います．媒介効果、というのは、「ある要因がある結果を引き起こすときに、何らかの因子を媒介するかどうか」ということを指し、その分析を行うことが媒介効果分析になります.

概要を理解するためには最も単純な「要因」「結果」「媒介因子」の３つを矢印でつないだ図を理解するところから始めると思いますが、もう少し込み入った因果効果の中で理解するためにはCausal mediation analysisというものを理解する必要があります.

本稿では実践の部分をStataのコードを示しながら解説することとし、理論的な詳しい説明は別の方がとてもわかりやすく記述してくださっているので、そちらを参考にしていただければよいと思います．

1．媒介効果分析の概要
2．コマンド”paramed”でオーソドックスな媒介分析
3．コマンド”med4way”を実装してみる
4．まとめ

1．媒介効果分析の概要

詳細は省きますが、最も重要なエッセンス部分だけは簡単に解説しておきます.

下の図3つを順に解説しますので、因果効果の表現方法に関する概観を何となくつかんでみてください.

図Aでは媒介因子を考慮しないモデルです．話をシンプルにするために交絡因子を考慮しないことにしています．この矢印の関係がtotal effectです．

図Bでは媒介因子Mを考慮したモデルです．図Aのtotal effectは図BにおいてはX→Yのcontrolled (あるいはconditioned) direct effectとX→M→Yのindirect effectに分けられます．

最も基本的な媒介効果分析はBK法という方法で分析しますが、以下のステップを経ます．

X→Mの回帰モデルで有意な関連を持つことを確認
X、Mを含めてYをアウトカムにした回帰モデルにおいて、Mが有意であることを確認．このとき、Xが有意でなくなったら完全媒介、有意なままであれば部分媒介
間接効果の有意性検定（ここでブートストラップを行うのが最も一般的）

この関係を示すのに最もシンプルなコマンドが”paramed“というコマンドです．

図3は、今日のメインである”med4way“を走らせるために必要な最低限の知識のベースです．通常のDAG（Directed acyclic graph）と異なり、augmented DAGと呼ばれるものなのですが、mediationとinteractionをよりわかりやすく表現したもの、ととりあえず理解しておきます．

交互作用の部分をXMで表現して別の箱にいれることでわかりやすくなります．X→Yに向けて行く道は下図のように4通りあります．なのでこれをFour-way decomposition と呼びます．この道筋に沿って媒介分析を行います．

この４つの経路を通ってXからYに行くのでfour-way decompositionといいます

最も優れてるのは、このコマンドを使えば生存時間分析もそのまま使えるという点です．告白してしまえば、その一点のみで飛びついた、というのが正直なところです．

これを理解する前に最もオーソドックスな媒介分析をStataで行ってみましょう．

2．コマンド”paramed”でオーソドックスな媒介分析

このコマンドはlinear, logistic, loglinear, Poisson, Negative binomial の回帰モデルで適用できますので、横断的な関係をみるのには使えます．デフォルトでは入っていませんので、いつものようにhelpあるいはfinditで取りに行きましょう．

syntaxは以下の様になっています．

paramed depvarname, avar(indepvarname) mvar(medvarname) a0(#) a1(#) m(#) ///
 yreg(regressname) mreg(regressname) [cvars(varlist) nointeraction casecontrol ///
 fulloutput c(numlist) bootstrap reps(integer1000) level(cilevel) seed(passthru)]

paramedの後すぐにアウトカム変数を置き、その後のオプションには、

avar: 説明変数
mvar: 媒介変数
a0: natural level of treatmentなので、通常は0/1なら0
a1: 曝露状態なので1
m: mediatorが存在している状態のときの変数の値
yreg: outcomeに対する回帰モデル
mreg: mediatorに対する回帰モデル
cvar: 調整変数群
c: 条件付き効果が計算される共変量に指定する固定値

具体的な使用例を交えた細かい説明は次のmed4wayのところで解説します．この2つのコマンドはよく似ていますので．

3．コマンド”med4way”を実装してみる

さて、先ほども述べたように、このコマンドは生存時間分析にも使用できる点が特徴です．

入手先がちょっと調べにくいのですが、以下をコピペしてもらえれば取りに行けます．（ネットワーク環境下でお願いしますね）

net install med4way, from("https://raw.githubusercontent.com/anddis/med4way/master/") replace

Discacciati et al. Med4way: a Stata command to investigate mediating and interactive mechanisms using the four-way effect decomposition. International Journal of Epidemiology 2019; 15–20によれば、このコマンドを著者達はGithubに公開してくれています．素敵すぎます．惚れます．

実際にこれを実装した論文も見つけました．日本の研究者の書いた論文で、炎症が非心臓手術後のAKI予測因子であり、術後AKI後の死亡の媒介因子である、という内容の論文です（Murashima et al. Infammation as a predictor of acute kidney injury and mediator of higher mortality after acute kidney injury in non-cardiac surgery. Scientific Reports 2019; 9: 20260）．オープンアクセスなので内容を確認できます．

Abstract: The proportions explained by mediating effects of lnCRP and albumin were 29.3% and 39.2% and mediation effects were statistically significant.
「媒介因子が主要因とアウトカムの間の関係性を～％説明するとわかりました．」
Result: By mediation analyses, regression coefficients (95% CI) for pure indirect effect (mediation effect) were 0.09 (0.05–0.14) (p < 0.001) and 0.16 (0.10–0.23) (p < 0.001) for lnCRP and albumin levels, respectively. The proportion explained by mediating effect of lnCRP and albumin was 29.3% and 39.2%, respectively.
Murashima et al. Sci Rep. 2019; 9: 20260

という感じで結果には記載されていますので、そういう使い方ができます．わかりやすいですね．

さてこのコマンド、”med4way”ですが、syntaxは先ほどのparamedとよく似たような見た目です．

med4way [yvar] avar mvar [cvars] [if] [in], a0(#) a1(#) m(#) yreg(string) mreg(string) [ options ]

yvar: 目的変数（生存時間分析のときには不要）
avar: 説明変数
mvar: 媒介変数
a0: natural level of treatmentなので、通常は0/1なら0
a1: 曝露状態なので1
m: four-way decompositionが計算される媒介変数の値
yreg: outcomeに対する回帰モデル
mreg: mediatorに対する回帰モデル
cvar: 調整変数群
c: 条件付き効果が計算される共変量に指定する固定値

exampleのデータから実例を紹介します．

use med4way_example_2.dta
stset y_cens, failure(fail) noshow 　　/* 生存分析の宣言 */
med4way treat m_bin cvar1, a0(0) a1(1) m(0) yreg(aft, e) mreg(logistic) c(1)
    /* パラメトリックな生存分布(exponential survival distribution), mediatorに対してはlogistic regression*/
test _b[ereri_cde] = _b[ereri_pie]
   /* controlled direct effectによるexcess relative riskがpure indirect effectによるexcess relative riskと統計学的に異なるかを検定（P大きい） */
med4way treat m_bin cvar1, a0(0) a1(1) m(0) yreg(aft, e) mreg(logistic) c(1) fulloutput   /* 最終結果 */