Multiple imputationについて①理論編 | 医療統計とStataプログラミングの部屋

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

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

さて今回の記事では、多重補完について、重要ポイントを説明した上でStataでの実装の仕方、ちょっとしたTipsについて解説したいと思います．長くなるので2回に分けてお送りします．

まず、本稿では理論編として多重補完の前提、入れるべき因子、欠損の割合などについて解説してみたいと思います．今回の記事を作成するにあたり、実際の研究をする際に統計家から直接おしえてもらったことを中心に解説します．

リンク

1．欠測のメカニズムと多重補完法（multiple imputation）
2．multiple imputationにおいて入れるべき因子
3．代入するアルゴリズムと変数の扱いについて
まとめ

1．欠測のメカニズムと多重補完法（multiple imputation）

1．欠測のメカニズム

欠損値がどのようなメカニズムで欠測しているか、そしてどんな場合に多重補完が適した方法なのか、という理論的な解説を少しだけしておきましょう．

missing completely at random（MCAR）：欠測値が完全にランダムに発生している場合
missing at random（MAR）：欠損値の発生が完全にはランダムではないものの、欠損しているか否かが欠損値そのものに依存していない場合
missing not at random （MNAR）：欠損値の発生が欠損値そのものに依存する場合．例えば婚姻状態についての申告者の割合がもしかすると未婚者・既婚者・離婚者で異なるかもしれない場合はそのような欠測のメカニズムはこれに相当します．

これを実際の論文でどのように記載していくか、というのは難しい問題です．研究分野ごとの専門知識に従って推定するしかないと思います．

2．欠測値を補う方法

欠損値を含む症例を全部除外して分析する方法をcomplete-case analysis（CCA）といいます．Stataでは欠測値を含むデータを多変量解析に投入すると自動的に出てくる結果はこのCCAの結果になります．MCARの場合や欠損値を含む症例の割合が少ない場合、完全ケース分析は許容されますが、それ以外の場合では、完全ケース分析はデータの推計結果に歪みをもたらすため、何らかの対処が必要とされます．

1つの対処方法として、「欠損値」という新たなカテゴリーを作ってモデルに組み込む方法がありますが、「欠損値」カテゴリーにおける推定値の解釈が困難です．

それ以外の方法としては、欠損値に何らかの数値を代入し、欠損値を補完する方法があります．

1）平均値代入法（mean imputation）
2）回帰代入法（regression imputation）
3）last observation carried forward（LOCF）　　　　など

しかし、いずれもデータの推計値に歪みが生じるため、推奨されません．そこで多重補完法（multiple imputation）の出番です．

ちなみに本来はMARにおいて多重補完法は適応すべきとされていますが、実際にはこのようなメカニズムは簡単に区別が付くわけではなく、しばしば妥当な方法としてこの多重補完法が採用され、解析を実施します．

3．多重補完法のポイント

多重補完法はシミュレーションを利用したモデル推定の手法であり、欠損値を復元するものではない
代入モデルが適切、かつ解析モデルが適切という強い仮定が必要
代入回数は5－20回が適切と言われるが、その目安となるような指標も後述する

以上のようなポイントがあります．

2．multiple imputationにおいて入れるべき因子

１．MIにはアウトカムをいれるべきなのか

では実際にどんな因子を多重補完においては投入すべきなのでしょうか．結論から言えば、「アウトカムも含めて全部」となります．以下参考文献です．

In principle, the imputation model should contain all variables that are included in the main analysis as well as those that may predict both the incomplete variable and its incompleteness. Most importantly, failure to include the outcome in the imputation model can weaken the apparent association between the outcome and the imputed covariate.
Furukawa et al. Statistical Methods in Medical Research 2017, Vol. 26(2) 707–723

根拠となった論文は、Moons KG, Donders RA, Stijnen T, et al. Using the outcome for imputation of missing predictor values was preferred. J Clin Epidemiol 2006; 59: 1092–1101. です．

２．生存時間分析のときは時間と打ち切りをモデルに入れるべき

（2021/2/27修正）

他の回帰モデルと同様、Cox回帰での多重補完においてもアウトカム変数を入れるべきとされています．方法は様々提案されているようです．

For Cox regression, it has been shown that an imputation model conditional on D and the cumulative baseline hazard H0(T) can lead to approximately valid results in the main analysis.
Furukawa et al. Statistical Methods in Medical Research 2017, Vol. 26(2) 707–723

より具体的に言及しているのがこちらです．

Nelson-Aalen estimate of the cumulative hazard to the survival time
White IR and Royston P. Imputing missing covariate values for the Cox model. Stat Med 2009; 28: 1982–1998.

現実的には、survival outcome としてT(時間)あるいはln(T) と D(打ち切り変数)を入れてもさほど結果はずれることはないのだそうです．

なお、このNelson-Aalen estimateについては算出方法が先日（2021/2/24）のLightstoneさんのWebinarで紹介されていましたので、そちらで実施されるのがより確実なのかもしれません．

具体的には、sts generate 新変数 = na とするとこの推定量が各症例において得られますのでこれを組み込めばよいということになります．

３．欠測の割合はどのくらいが妥当なのか？

欠測値の割合と推定の妥当性についてはよく知られている経験則として、全体の25％以上の欠測があるとよくないかもしれない、ということがあります．Applied Missing Data Analysisの本でもシミュレーションをしているところでは25％のデータを欠測させたりしています．

リンク

（Demirtas et al. An imputation strategy for incomplete longitudinal ordinal data. Stat Med. 2008 Sep 10;27(20):4086-93が根拠となっているようです．）

しかしながら、これに対しては、「欠測割合は気にしなくてよい」とする論文もあります(Journal of Clinical Epidemiology 110 (2019) 63-73)．代入モデルに十分に情報が与えられれば、大きな欠測割合の変数に対してもMIでバイアス、効率性の改善がある、ということが述べられています．

もう少し詳しく説明すると、欠測割合だけでは不十分であり、fraction of missing informationという考え方のほうが便利とのことです．そして一つ一つの変数に対する回帰式が十分な説明力（ここではR²）を有しているとき、MIによるばらつきは最小化できる、ということのようです．

The proportion of missing data is a common measure of how much information has been lost because of missing values in a data set. However, it does not reflect the information retained by auxiliary variables. Alternative measures such as the fraction of missing information (FMI) may be more useful as a tool for determining potential efficiency gains from MI. The FMI is a parameter-specific measure that is able to quantify the loss of information due to missingness, while accounting for the amount of information retained by other variables within a data set.
P. Madley-Dowd et al. Journal of Clinical Epidemiology 2019

もっとも、ここでいうモデルのR²は一番小さくて0.36、大きくて0.92とかなので、実は結構厳しいんじゃないかな～などと思ったりします．このFMIという指標もイマイチよくわかりませんので、統計の先生と相談しながらやるのが安全、というのが今のところの個人的な意見です．

3．代入するアルゴリズムと変数の扱いについて

1．代入するアルゴリズムの選択方法（2021/2/27　追記）

最も簡単なのが線形回帰モデルに基づいた代入法です．これは欠損のある変数が①連続変数であって②正規分布するという仮定に基づいています．正規性が成り立たない場合にはバイアスが生じてしまいます．

ここで、その他にもある代入モデルも含めて何を選択したらよいのかをきちんと整理しておきたいと思います．これもStataのWebinarから仕入れたネタになります．

モノトーンパターンかどうか（misstable nestedというコマンドで確認できる）
単一の変数の代入が必要か、複数の変数の代入が必要か
連続変数で正規性があるかどうか

最初にモノトーンパターンかどうかで分けるのは、モノトーンパターンだと計算が少なく楽だからです．しかし現実のデータではそこまでシンプルにはいかないようです．

単一の変数の場合には

①連続変数で正規分布ならmi impute regress

②正規分布でない場合には、その他（mi impute intreg, mi impute logit, mi impute mlogit, mi impute nbreg, mi impute ologit, mi impute pmm, mi impute poisson, mi impute truncreg）

複数の変数で代入が必要になった場合には

ベイズ統計を利用したデータ拡大アルゴリズム（mvn）もしくはchainedを用いることになります．

mi impute mvn を用いる場合には以下の条件が必要になります．

欠損データが連続変数で正規分布に従う
MCMCを使うので代入実行後に収束判定に関する診断をする必要がある．（ベイズ統計についてはまた項を改めて説明したいと思います．）

2．MICEかMVNか？

恐らく最もポピュラーな方法の1つが、Multivariate imputation by chained equations (MICE)というものです．多変量欠測データのimputationを変数ごとに回帰モデルを作って繰り返す、というものです．回帰モデルはregression, logistic, multinomial logisticなどの中から適切なモデルを変数ごとに選択します．これはもう一つのよく使われる方法であるmultivariate normal regression (MVN)と比べて次の様な利点があります．

MVNは正規分布する連続変数のみの代入に使うが、MICEは幅広くパラ・ノンパラ、バイナリー、その他カテゴリ変数につかえる
複数の変数で同時に欠測が生じているときにも多重補完が可能

MVNですと、０／１の変数もすべて連続値に置き換わりますし、マイナスの値すら取り得るためその後にlogistic 回帰を実行できないなどの問題が生じてしまいます．

MICEですと、問題としては、前提となる仮定として「各変数の欠測に対する補完モデルとして、欠測を補完したい変数以外の他のすべての変数によってお互いに説明可能」という強い仮定を必要とすること、「この仮定を満たせない場合に得られた最終的な推定値の分散が正しく推定できるか保証されていないこと」があります．非常に強い仮定で、MICEにおけるfully conditionally specifiedの仮定と言います．

とはいえ複数の変数の欠測に対処可能という柔軟性を持っているので、よく使用されている補完アルゴリズムの1つです．

それと非常に計算量が多く、時間がかかるということです．時間がかかった挙げ句に収束しない、というオチも経験しますので、時間に余裕を持ってスケジューリングしましょう(笑)

3．変数の扱い

基本的にはMIでは多変量正規分布を仮定しているため、正規分布でない連続変数については対数変換をするなど何らかの処理が必要になります．ただ正規性が崩れても大きな問題はないとも言われているようです．

また、交互作用を見るような状況では交互作用項もいれるとよいようです．

まとめ

多重補完を行うときには欠測値のメカニズム、入れるべき変数、アルゴリズム、変数の扱いに注意する必要があります．次回はStataでの実装方法についてまとめてみたいと思います．