Unboundedly

日々大学院で学んだこと、考えたことを更新

観察データを用いた因果推論で生じるバイアスの程度を考える:感度分析(Sensitivity analysis) & "E-value"入門

さて、今回からイデオロギー色の少ないブログ活動に戻ります。前回の受動喫煙に関する記事の中で、「感度分析("Sensitivity Analysis")」というテクニックを紹介しました。

どうも私の言葉足らずか、このテクニックに関して多くの方に誤解と混乱を招いたようなので今回改めて整理してみたいと思います。

なぜ感度分析が有効なのか?その意義は?

観察データを用いた分析から因果関係を言うためには、様々なバイアスに対処しなければいけません。数あるバイアスの中でも最もよくあるものが、交絡によるバイアスです。DAG上ではアウトカムと曝露因子の共通原因によって引き起こされるバイアスと整理されます。詳しくは以下の過去記事を参照ください。

さて、この交絡バイアスですが様々な対処方法があります。概要は上記の記事にまとめてあるので興味がある方は読んでみてください。一般的には、見たい因果関係を交絡する交絡因子に関するデータをすべて把握し、研究デザインや統計解析を通してその影響を制御する、というのが基本方針になります。

あらゆる交絡因子の影響を制御した場合、観察データからであってもドンピシャ正確な因果関係が得られるのです。これは統計的因果推論の数学的背景を確認すればすぐわかります。

問題は、「すべての交絡因子を把握、対処」という前提が成り立つことが実質不可能という点です。以前も書いたことなのですが、我々は神様ではありませんので注目している因果関係(例:喫煙と肺がん)に関係する交絡因子をすべて把握することなどできませんし、分かっていたとしてもその交絡因子に関するデータがとれるかどうかは別の議論です。金銭的・倫理的理由からデータがとれないことも多いでしょう。

このような未測定の交絡因子(Unmeasured Confounders)の影響が、観察研究から得られた結果を解釈するうえで非常に重要なポイントになります。ランダム化をしていない以上、多かれ少なかれ未測定交絡因子は必ず存在します。したがって、得られた推定結果は数学的に定義された平均因果効果(Average Causal Effect)とは異なる値になっています。つまり因果関係を”正確に”推定できていないことになりますので、よくこの辺りを批判されます。例えば、統計学の世界で著名なフィッシャーが喫煙の健康被害について懐疑的であったという話は有名ですが、フィッシャーの批判も上記の点を問題視したものでした。彼は「遺伝要因」による交絡の可能性を指摘していたそうです。以下のようなイメージ。タバコを吸いやすくなる遺伝子って何だよ、とつっこっみたくなるのですが、まあ仮にそういうものが存在するした場合、たしかに交絡によって計算された因果効果がゆがめられていそうです。

f:id:KRSK_phs:20170628090941p:plain

しかし、ここで考えるべきは「未測定交絡因子によってどのくらいの因果効果の推定がゆがめられているか(どの程度のバイアスがあるのか)」がありうるか、ということです。バイアスの程度に関する具体例はこちらを再読してみてください。未測定交絡因子によるバイアスの可能性を指摘するのであれば、どのような交絡因子が存在しうるのか、それによって観察されている関連は完全に説明できてしまうのか、逆に観察されている関連が完全に交絡によるものだというためにはどのような交絡因子が存在する必要があるのか、といった点をディスカッションしなければ建設的な批判とはいえません。

交絡によるバイアスを考える際には、①バイアスの方向(Direction)②その大きさ(Magnitude)の両方を検討する必要があります。この二つは①交絡因子とアウトカムの関連、②交絡因子と曝露因子の関連、③注目している集団における交絡因子の存在割合(Prevalence)という三つの要素によって決まります。これは超重要なポイントです。以下に詳細を説明します。

交絡バイアスの向き

バイアスにも向きがあります。例えば、喫煙が肺がんに与える”真の”因果効果がリスク比1.5だったとしましょう。未測定交絡因子によるバイアスがある場合、実際に観察される喫煙と肺がんの関連は1.5とは異なる値をとることになります。この時、観察されたリスク比が2.0だったとすると、真の因果効果よりも過大に因果効果を見積もっていることになります。これを上向きバイアス(upward bias)と呼ぶことにしましょう。逆に観察されたリスク比が1.2であった場合、真の因果効果よりも過少に因果効果を見積もっているので下向きバイアス(downward bias)が生じていることになります。

基本的には注目している因果関係が存在しない場合の値(Null value;リスク比=1またはリスク差=0)から離れていくようなバイアスのほうが問題だと考えます。上記の例では、交絡因子によるバイアスが上向きの場合は真の値が観察されている関連よりも小さいため、もしかすると真の因果効果がリスク比1(Null value)かもしれず注意が必要です。しかし、バイアスが下向きならば「真の因果効果は観察されているものよりももっと強い」といえるので「喫煙が肺がんのリスクをあげる」という結論自体は変わらなさそうです。

では、どうやってバイアスの向きを予想するのか。実は非常に簡単です。交絡因子とアウトカム、交絡因子と曝露因子それぞれの関連の向きに注目すれば、交絡バイアスの向きが分かります。以下の図において、は正の関連(片方が増えると、もう一方も増える)、は負の関連(片方が増えると、もう一方は減る)を表すとしてください。

f:id:KRSK_phs:20170628100554p:plain

 要するに、交絡因子と曝露因子の関連、交絡因子とアウトカムの関連の向きが同じならば上向きバイアス、向きが逆なら下向きバイアスということになります。フィッシャーの例をとるならば、ある遺伝子を持っている人で肺がんのリスクが高い(+)かつタバコを吸いやすい(+)、もしくは肺がんのリスクが低い(-)かつタバコを吸わない(-)ならば上向きバイアスが生じることになります。仮に、交絡バイアスを引き起こす遺伝要因があったとしても、遺伝子と喫煙・肺がんそれぞれとの関連の向きが一致していない場合は、因果効果が”ない”可能性の指摘としては的外れ(むしろ因果関係の存在をサポートする)といえます。

交絡バイアスの大きさ

交絡バイアスの向きが仮にNull valueに近づくものであったとして、でははたして考えうる未測定交絡因子は観察されている関連を交絡バイアスだけで引き起こしてしまうほど強力なものなのか、どのような未測定交絡因子が存在すればそのような状況になるのかを考える必要があります。

このようにバイアスの大きさを考えるためのテクニックとして(交絡バイアスのための)感度分析(Sensitivity Analysis)があります。一般に感度分析という言葉は、統計モデリングのやり方を変えたり、異なる定義をつかって変数を定義したりなどすることで結果の頑健性を確認するテクニックのことを指しますが、今回は交絡バイアスの大きさの程度に目星をつける手法のことだと思ってください。

”Without Assumptions”とは

これより先は、前回の記事の中で引用していた以下の論文で紹介されている感度分析を中心に書きます。

Ding, Peng, and Tyler J. VanderWeele. "Sensitivity analysis without assumptions." Epidemiology (Cambridge, Mass.) 27.3 (2016): 368.

最初に結論だけ言っておくと、交絡バイアスの大きさは①交絡因子とアウトカムの関連の強さ②交絡因子と曝露因子の関連の強さ③集団における交絡因子の存在割合によって決まります。

上記の論文のタイトルに”without assumptions"とあります。”前提を置いていない”ということです。以下紹介する感度分析のテクニックは①と②の条件のみを使い、③の条件を考慮していないため、このように書かれています。*1実は注目している集団において交絡因子が珍しいほど(割合0%に近づくほど)、もしくは逆に非常にありふれたものであるほど(100%に近づくほど)、交絡バイアスは弱くなることが知られています。

この交絡因子の存在割合に関する条件を指定していないため、以下の手法ではバイアスの大きさが過剰に見積もられています。言い換えると、①②の条件下で起こりうる最大のバイアスを推定しています。交絡因子のprevalenceによっては、推定されたよりも実際のバイアスは小さい可能性もありますので、これからご紹介するテクニックは非常に保守的(Conservative)な方法であると言われています。

また、以下では説明を簡単にするため未測定交絡因子、曝露因子、アウトカムすべて二値の変数であると仮定してください。今回は省略しますが、連続変数の交絡因子に対する感度分析も基本的には同じアイデアです。

古典的Cornfieldの条件

以下、簡単な数学表現も使うので一度用語を定義しておきます。

A:曝露因子(例:喫煙)
Y:アウトカム(例:肺がん)
U:未測定交絡因子(例:遺伝子)
RR: A-Yの関連
UA:U-Aの関連
UY:U-Yの関連

いわゆるCornfield conditionsと呼ばれる古典的な条件は次のようなものです。

”UA>RRかつUY>RRでないならば、Uだけで観察されている関連(RR)を説明することができない”

言い換えると、RRがUAまたはUYの少なくともどちらか一方よりも大きい場合は、交絡バイアスだけで観察されている関連を説明できないということです。

フィッシャー大先生の具体例を使いましょう。フィッシャーは喫煙Aと肺がんYのリスク比を10.73とする研究(Hammond and Horn, 1958)に対して、遺伝子Uによる交絡の存在を問題視しました。Cornfield条件を使うと、

UA > 10.73

UY> 10.73

を満たす遺伝子Uでなければ喫煙と肺がんの因果関係は覆されないと言えます。

Cornfield条件にはもう一つ条件があって、UA=10.73で観測されている関連を説明するためには、UとYの関連がPerfect(UY=∞)でなければいけないというものです。UY=10.73のときにはUA=∞です。

このような遺伝子が存在するかどうかは、その他の研究結果を踏まえて考えるのですが、基本的には分野の常識から言ってありえない条件だろうと言えます。したがってフィッシャー大先生の指摘する遺伝要因による交絡バイアスだけでは、喫煙と肺がんの因果関係を否定することは難しそうです。

仮想の交絡因子によるバイアス:Bounding factorをもとめる

Cornfield条件ではUA、UYの大きさから交絡バイアスで観察されている関連が説明されてしまうかが分かりました。実はUA、UYが分かっていればでバイアスの程度を直接求めることもできます。

B = UY * UA / (UY + UA -1)*2

上記の式から得られるBのことを筆者らはBounding factorと呼んでいます。このBはバイアスの程度を表す数値であり、観測された関連をBで割ることで、交絡バイアスがない場合の真の因果効果を計算することができます。

再びRR=10.73の場合を考えます。UY=10.73、UA=10.73だったとしましょう。この時、

B=10.73*10.73 / (10.73 + 10.73 -1) = 5.63

観測されたリスク比10.73 ÷ 5.63 = 1.91 >> 1 (null value)なので、これほど強力な交絡因子であってもそれだけで観察された関連を説明することができないと分かります。

ちなみにB=UY * UA / (UY + UA -1)の式で

UY→∞でB=UA

UA→∞でB=UY

となりますのでCornfield条件と一致します。

Bounding factorを解釈するうえで、注意すべきポイントが2点あります。

1.得られるバイアスは最大(possible maximum)

さきほども書きましたが、交絡因子のPrevalenceに関する条件を指定していないので、B=5.63はUY=UA=10.73を満たすUが起こしうる最大のバイアスです。逆に言うと、UY=UA=10.73でも、実際のバイアスはもっと小さいことになります。これらの感度分析がConservativeと言われるゆえんです。

2.UY、UAはconditional on measured confounders

これも非常に重要なポイントです。元論文を見れば詳しく数学的背景も書いてあるのですが、実は上記の式に含まれているパラメータは”すでに調整されている交絡因子を考慮したうえでの”値を採択します。どういうことでしょうか。多くの研究では単純に曝露因子とアウトカムの関連を見たりせず、回帰分析などを使って性、年齢、学歴など主な交絡因子の影響を調整したうえで注目している二つの因子の関連を見ていると思います。Bの計算に使うUA、UYはこれら交絡因子の影響を除いたうえでの値になります。つまり、

交絡因子が多く観測され、調整されているほどUA、UYが大きい値をとるのが難しい

→Bが大きい値を取りにくくなる

→未測定交絡因子によるバイアスだけで説明がつきにくくなる

ということです。

どんな交絡因子ならば観測された関連を説明できるか

実際にはUA、UYの値はわかりません(Uが未測定なので当然です)。なので、ほかの研究から得られた知見からUAとUYのおおよその値を予想したり、さまざまな値を入れてみたときのBの一覧を表にするなどの方法がとられます。

発想を逆転させて、「観察されている関連を完全に交絡バイアスだけで説明するためには、未測定の交絡因子はどれくらい強力である必要があるか」を求めることができます。それには次のような条件を使います。

UYとUAがどちらも

RR + sqrt[RR * (RR-1)]*3

よりも大きい場合は交絡バイアスだけで観測されている関連を説明できるが、これより”弱い”交絡では説明できない。

RR=10.73の場合、10.73+sqrt[10.73 * (10.73-1)] =  20.95なのでUY=UA=20.95かそれ以上の交絡因子UならばRR=10.73をすべて交絡バイアスで説明できるが、これより”弱い”交絡では説明できないということです。

これはUYとUAが必ず20.95以上でないといけないという意味ではありません。UYが20.95を下回る場合、UAが20.95よりも大きくなければUY=UA=20.95のUが起こす交絡バイアスよりも”弱い”交絡となってしまいます。UYが小さくなればUAがさらに大きい必要があるということです。

 ①最大のバイアスであるためConservativeであること

②調整済交絡因子の影響を考慮したうえでのUA,UYであること

以上の二点に注意すべきという点はBounding factorと同じです。

 

ちなみに、この論文の共著者の一人であるVanderweele教授はこの値のことをp-valueにちなんでE-value (evidence for causality)と呼んでいます。

RによるBounding factor、E-valueの可視化

Bounding factorとE-valueの関係をグラフにしました。使用したRコードは以下の通り。Rコードのうまい貼り付け方法ってあるんですかね?なんか味気ない、、、詳しい人いたら教えてください。

library(boot)
#観測されたリスク比=10.73
RR <- 10.73
# X軸とY軸の最小値、最大値の設定
a <- 0
b <- 50
# E-valueを計算
evalue.estimate <- RR + sqrt(RR*(RR-1))
min <- signif(RR+sqrt(RR*(RR-1)), digits=3)
# UA, UYの値を変化させたときのBounding factorの変化の動きをプロット
# UA=UY=20.9と”同じ”強さの交絡を引き起こすUAとUYの組み合わせ
# UA=xとする
curve((RR*(1-x))/(RR-x), from = RR+0.001, to = b, xlim = c(a, b), ylim = c(a, b),
      xlab=expression("UA"), ylab=expression('UY'))
points(evalue.estimate, evalue.estimate, pch=19)
text(min, min, pos=4,
     bquote(paste("(", .(min), ', ', .(min), ')' )), cex=0.75)
legend('bottomleft', as.expression(bquote('UA'*'UY'/('UA'+'UY'-1)
                                            == .(RR))),
       lty=1:2, cex=0.75)

 アウトプットが以下のグラフです。

f:id:KRSK_phs:20170628130148p:plain

上記の線は、UがBounding factor = 10.73(すなわち、観測されたリスク比をすべて説明できてしまう)の未測定交絡因子となるようなUY=UAの組み合わせをプロットしたものです。UY=UA=E-value=20.9ならば、交絡バイアスだけで説明がつくことがわかります。

さらに注目すべき点はUYが20.9より小さくなるほど(グラフの右のほうにくほど)、UAが大きくなっていることです。逆もまたしかりです。これが「UY=UA=20.9の交絡よりも”弱い”交絡では観測された関連を説明できない」の意味するところです。

UY=∞でUA=10.73に収束していきますし、UA=∞でUY=10.73となりますのでCornfield条件が成立していることも確認できます。

 

以上が未測定の交絡因子によるバイアスの程度を考える際に便利な感度分析の概要です。疫学の世界でも比較的新しい概念で、実際の論文で使われることはまだまだ稀です。しかし、これらのテクニックを使って個々の観察研究のバイアスの程度を考えることでより建設的な議論・批判的吟味が行われ、分野全体で研究の質が上がっていくのではと期待されています

 

 

*1:実はそれ以外にもいくつか前提条件を置いていないのですが、詳しくは論文を読んでください

*2:これはUが一つだけの時の式ですが、Uが複数ある場合もUYとUAを定義し、Bを求めることができます

*3:RRが1より小さい場合は逆数をとります