第 92 章 Regression Methods with binary outcomes 結果變量爲二分類變量
92.1 二分類結果變量的因果被估計量 (causal estimand):
Average causal effect (因果邊際危險度差，marginal causal risk difference), ACE :
\[ \text{Pr}\{Y(1) = 1\} - \text{Pr}\{Y(0) = 1\} \]
因果邊際危險度比 (marginal causal risk ratio):
\[ \frac{\text{Pr}\{Y(1) = 1\}}{\text{Pr}\{Y(0) = 1\}} \]
或，因果邊際比值比 (marginal causal odds ratio):
\[ \frac{[\frac{\text{Pr}\{Y(1) = 1\}}{1-\text{Pr}\{Y(1) = 1\}}]}{[\frac{\text{Pr}\{Y(0) = 1\}}{1-\text{Pr}\{Y(0) = 1\}}]} \]
或，因果邊際對數危險度比/比值比:
\[ \log\{\text{Pr}\{Y(1) = 1\}\} - \log\{\text{Pr}\{Y(0) = 1\} \} \]
\[ \log[\frac{\text{Pr}\{Y(1) = 1\}}{1-\text{Pr}\{Y(1) = 1\}}] - \log[\frac{\text{Pr}\{Y(0) = 1\}}{1-\text{Pr}\{Y(0) = 1\}}] \]
和上一章節一樣，我們可能實際上還會關心這些被估計量的調整後的條件平均因果效應 (conditinal ACE):
\[ \text{Pr}\{Y(1) = 1 | \mathbf{V=v} \} - \text{Pr}\{Y(0) = 1 | \mathbf{V=v}\}\\ \log\{\text{Pr}\{Y(1) = 1| \mathbf{V=v}\}\} - \log\{\text{Pr}\{Y(0) = 1| \mathbf{V=v}\} \}\\ \log[\frac{\text{Pr}\{Y(1) = 1 | \mathbf{V=v}\}}{1-\text{Pr}\{Y(1) = 1 | \mathbf{V=v}\}}] - \log[\frac{\text{Pr}\{Y(0) = 1 | \mathbf{V=v}\}}{1-\text{Pr}\{Y(0) = 1 | \mathbf{V=v}\}}] \]
92.1.1 比值比的不可壓縮性 non-collapsibility of the odds ratio
即便是沒有效應修飾，在 GLM 章節我們也學到過，由於邏輯回歸模型的不可壓縮性質，一般地，選擇的條件變量不同的話，(對數)比值比的大小都會發生變化。所以沒辦法用回歸系數的變化來推斷是否有明顯的混雜效應。
92.2 鑑定 identification - conditional effects
如果我們對暴露 \(X\) 和結果 \(Y\) 之間的條件因果危險度差 (conditional causal risk difference):
\[ \begin{aligned} \text{Pr}\{ Y(1) = 1 & | \mathbf{C=c}\} - \text{Pr}\{ Y(0) = 1 | \mathbf{C=c}\} \\ & = \text{Pr}\{ Y(1) = 1 | X=1, \mathbf{C=c}\} - \text{Pr}\{ Y(0) = 1 | X=1, \mathbf{C=c}\} \\ & \text{By conditional exchangeability given }\mathbf{C} \uparrow \\ & = \text{Pr}\{ Y = 1 | X=1, \mathbf{C=c}\} - \text{Pr}\{ Y = 1 | X=1, \mathbf{C=c}\} \\ & \text{By consistency } \uparrow \\ \end{aligned} \]
相似地，條件因果對數危險度比 (conditional causal log risk ratio):
\[ \begin{aligned} \log[\text{Pr}\{ Y(1) = 1 &| \mathbf{C=c} \}] - \log[\text{Pr}\{ Y(0) = 1 | \mathbf{C=c} \}] \\ & = \log\{\text{Pr}( Y = 1 | X=1, \mathbf{C=c})\} - \log\{\text{Pr}( Y = 1 | X=1, \mathbf{C=c})\} \\ \end{aligned} \]
條件因果對數比值比 (conditioanl causal log odds ratio):
\[ \begin{aligned} & \log[\frac{\text{Pr}\{Y(1) = 1 | \mathbf{C=c}\}}{1-\text{Pr}\{Y(1) = 1 | \mathbf{C=c}\}}] - \log[\frac{\text{Pr}\{Y(0) = 1 | \mathbf{C=c}\}}{1-\text{Pr}\{Y(0) = 1 | \mathbf{C=c}\}}] \\ & =\log[\frac{\text{Pr}\{Y = 1 | X = 1, \mathbf{C=c}\}}{1-\text{Pr}\{Y = 1 |X = 1, \mathbf{C=c}\}}] - \log[\frac{\text{Pr}\{Y = 1 |X = 0,\mathbf{C=c}\}}{1-\text{Pr}\{Y = 1 |X = 0, \mathbf{C=c}\}}] \\ \end{aligned} \]
92.3 鑑定 identification - marginal effects
92.3.1 Marginal causal risk difference (ACE)
\[ \begin{aligned} \text{Pr}\{ Y(1) =1 \} & - \text{Pr}\{ Y(0) =1 \} \\ = & \sum_c\text{Pr}\{ Y(1)=1 |C = c \}\text{Pr}(C=c) \\ & - \sum_c\text{Pr}\{ Y(0)=1 | C=c\}\text{Pr}(C=c) \\ & (\text{by the law of total probability } \uparrow) \\ = &\sum_c\text{Pr}\{ Y(1)=1|X=1, C=c \} \text{Pr}(C=c) \\ & - \sum_c\text{Pr}\{ Y(0)=1|X=1, C=c \} \text{Pr}(C=c) \\ & (\text{by conditional exchangeability } \uparrow) \\ = &\sum_c\text{Pr}( Y=1|X=1, C=c ) \text{Pr}(C=c) \\ & - \sum_c\text{Pr} (Y=1|X=1, C=c) \text{Pr}(C=c) \\ & (\text{by consistency } \uparrow) \\ = & \sum_c\{ \text{Pr}( Y=1|X=1, C=c) - \\ & \;\;\;\;\; \text{Pr}( Y=1|X=1, C=c) \}\text{Pr}(C=c) \end{aligned} \]
92.3.2 Marginal causal log risk ratio
\[ \begin{aligned} & \log[\text{Pr}\{ Y(1) = 1 \}] - \log[\text{Pr}\{ Y(0) =1 \}] \\ & = \log[\sum_c\text{Pr}(Y = 1|X=1, C=c)\text{Pr}(C=c)] \\ & \;\;\;\; - \log[\sum_c\text{Pr}(Y = 1|X=0, C=c)\text{Pr}(C=c)] \\ \end{aligned} \] ### Marginal causal log odds ratio (cannot be calculated)
\[ \begin{aligned} & \log[\frac{\text{Pr}\{Y(1) = 1\}}{1-\text{Pr}\{Y(1) = 1\}}] - \log[\frac{\text{Pr}\{Y(0) = 1\}}{1-\text{Pr}\{Y(0) = 1\}}] \\ & = \log\{ \frac{\sum_c\text{Pr}(Y = 1|X=1, C=c)\text{Pr}(C=c)}{1-\sum_c\text{Pr}(Y = 1|X=1, C=c)\text{Pr}(C=c)} \} \\ & \;\;\;\; - \log\{ \frac{\sum_c\text{Pr}(Y = 1|X=0, C=c)\text{Pr}(C=c)}{1-\sum_c\text{Pr}(Y = 1|X=0, C=c)\text{Pr}(C=c)} \} \end{aligned} \]
92.4 通過邏輯回歸估計這些被估計量
Log_lbw <- glm(lbweight ~ as.factor(mbsmoke) + mage + as.factor(fbaby) + as.factor(prenatal), family = binomial(link = "logit"), data = cattaneo2)
summary(Log_lbw)##
## Call:
## glm(formula = lbweight ~ as.factor(mbsmoke) + mage + as.factor(fbaby) +
## as.factor(prenatal), family = binomial(link = "logit"), data = cattaneo2)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.15417 -0.33142 -0.30949 -0.29641 2.60453
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.470282 0.402678 -1.1679 0.24285
## as.factor(mbsmoke)1 0.775782 0.137371 5.6474 1.629e-08 ***
## mage -0.021202 0.012607 -1.6817 0.09263 .
## as.factor(fbaby)1 -0.088391 0.136386 -0.6481 0.51692
## as.factor(prenatal)1 -1.950799 0.274350 -7.1106 1.155e-12 ***
## as.factor(prenatal)2 -1.912776 0.299804 -6.3801 1.770e-10 ***
## as.factor(prenatal)3 -2.101474 0.430343 -4.8833 1.043e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 2115.30 on 4641 degrees of freedom
## Residual deviance: 2026.64 on 4635 degrees of freedom
## AIC: 2040.64
##
## Number of Fisher Scoring iterations: 5但是，在邏輯回歸的模型下，即使是滿足無互相幹擾，一致性和條件可置換性的前提，並且你就算是 100% 自信地認爲你的模型絕對正確，計算獲得的條件比值比 (conditional odds ratio) 總也無法被賦予因果關系的含義。這是由於邏輯回歸的不可壓縮性 (non-collapsiblity)，這也是越來越多的人傾向與不使用比值比作爲評價治療效果 (treatment effect) 的指標的原因之一。
也因此，STATA 裏的 teffects ra 中即使你用的結果模型中加入 logit 的選項，它計算的是因果平均危險度差 (Marginal causal risk difference (ACE))。
##
## . use "backupfiles/cattaneo2.dta"
## (Excerpt from Cattaneo (2010) Journal of Econometrics 155: 138-154)
##
## . teffects ra (lbweight mage i.fbaby i.prenatal, logit) (mbsmoke)
##
## Iteration 0: EE criterion = 2.950e-26
## Iteration 1: EE criterion = 2.710e-34
##
## Treatment-effects estimation Number of obs = 4,642
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## lbweight | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATE |
## mbsmoke |
## (smoker |
## vs |
## nonsmoker) | .0583827 .0117242 4.98 0.000 .0354036 .0813617
## -------------+----------------------------------------------------------------
## POmean |
## mbsmoke |
## nonsmoker | .0503207 .0036151 13.92 0.000 .0432353 .0574061
## ------------------------------------------------------------------------------92.5 Average causal/treatment effect in the exposed/treated (ATET)
這是爲了回答公共衛生學上一個抽象的政策性問題: 對那些真的接受了治療/暴露/衛生政策幹預的人來說，他們身上發生的治療效果是怎樣的？因爲有些情況下你無法讓“所有人”都接受治療或幹預。
It is often of public health interest to ask “what is the effect of this exposure on those who choose to take it?” rather than “what would be its effect on everyone?”
\[ E\{ Y(1) - Y(0) | X=1 \} \]
此時，條件可置換性的前提發生了微妙變化:
\[ Y(0) \perp\perp X|\mathbf{C} \]
對於一個簡單的分類型條件變量 \(C\) 來說，它的 ATET 的鑑定過程如下:
\[ \begin{aligned} E\{ Y(1) -Y(0) |X =1 \} = & \sum_cE\{ Y(1) |X=1, C=c \}\text{Pr}(C=c|X=1) \\ & - \sum_cE\{ Y(0) |X=1, C=c \}\text{Pr}(C=c|X=1) \\ & \text{(by the law of total probability } \uparrow) \\ = & \sum_cE\{ Y(1) |X=1, C=c \}\text{Pr}(C=c|X=1) \\ & - \sum_cE\{ Y(0) |X=0, C=c \}\text{Pr}(C=c|X=1) \\ & \text{(by conditional exchangeability } \uparrow) \\ = &\sum_cE (Y |X=1, C=c)\text{Pr}(C=c|X=1) \\ & - \sum_cE(Y |X=0, C=c)\text{Pr}(C=c|X=1) \\ & \text{(by consistency } \uparrow) \\ = & \sum_c\{ E (Y |X=1, C=c) \\ & \;\;\;- E(Y |X=0, C=c) \}\text{Pr}(C=c|X=1) \end{aligned} \]
這時，我們只關心那些真正暴露的人 (predicted potential outcomes are predicted only for the exposed)。
在 STATA 的 teffects ra 後面加上 atet 的選項即可:
teffects ra (lbweight mage i.fbaby i.prenatal, logit) (mbsmoke), atet##
## . use "backupfiles/cattaneo2.dta"
## (Excerpt from Cattaneo (2010) Journal of Econometrics 155: 138-154)
##
## . teffects ra (lbweight mage i.fbaby i.prenatal, logit) (mbsmoke), atet
##
## Iteration 0: EE criterion = 2.950e-26
## Iteration 1: EE criterion = 1.149e-34
##
## Treatment-effects estimation Number of obs = 4,642
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## lbweight | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATET |
## mbsmoke |
## (smoker |
## vs |
## nonsmoker) | .0537169 .0114307 4.70 0.000 .0313131 .0761206
## -------------+----------------------------------------------------------------
## POmean |
## mbsmoke |
## nonsmoker | .0562368 .0045959 12.24 0.000 .047229 .0652447
## ------------------------------------------------------------------------------這裏你看到的是 ATET 和 ACE 很接近的數字，還有別的情況下，你會發現，某種治療方案對於接受治療的人來說是有好處的，但是對其他人是有害/沒有用的。
92.6 Practical04 - causal inference
注意: 這裏的練習使用的是 STATA 因爲，我在 R 裏找不到像 STATA 的 teffects 這樣靈活且方便的命令，如果你知道，歡迎告訴我: abelardccwang@gmail.com。
數據來自一個觀察性研究，樣本量是 3551 的一個肺癌患者數據庫，四家醫院的肺癌患者正準備選用 1) 常規手術 或者 2) 射頻消蝕 (radiofrequency ablation, RFA) 兩種方法取出肺部的轉移腫塊。該數據的變量如下:
##
## . use "backupfiles/RFA.dta"
##
## . describe
##
## Contains data from backupfiles/RFA.dta
## obs: 3,551
## vars: 14 5 Nov 2013 15:05
## size: 198,856
## -------------------------------------------------------------------------------
## storage display value
## variable name type format label variable label
## -------------------------------------------------------------------------------
## id float %9.0g Patient ID
## age float %9.0g
## gender float %9.0g gender
## smoke float %9.0g smoke Smoking status
## hospital float %9.0g Hospital ID
## nodules float %9.0g Number of nodules
## mets float %9.0g Number of other metastatic sites
## duration float %9.0g Duration of disease (in months)
## maxdia float %9.0g Diameter of largest nodule (in
## cm)
## primary float %22.0g primary Location of primary cancer
## position float %9.0g position Ease with which nodules can be
## reached
## coag float %9.0g coag Coagulopathy
## rfa float %23.0g rfa Treatment variable: RFA or
## standard surgery
## dodp float %9.0g dodp Outcome variable: death or
## disease progression within 36
## months
## -------------------------------------------------------------------------------
## Sorted by: id其中，主要的混雜因子是:
hospital: 有些醫院可能本身更傾向於/不傾向於使用 RFA，或者有些醫院的患者整體症狀較輕/較重;maxdia: 如果腫塊太大，那就不適合使用 RFA，而且腫塊較大的患者，生存的概率一般來說比較低;position: 腫塊位置，相比較傳統常規手術摘除的方法，RFA 能夠治療那些手術難以摘除的腫塊的部位。
這三個主要的變量被認爲非常重要，需要在分析中被調整。
其他的變量被認爲不太會是混雜因子，但是醫生認爲對患者的預後有很強的預測效果: age, gender, smoke, nodules, mets, duration, primary。
最後一點，對於有凝血障礙的患者 coag 來說，RFA 是不安全的。
92.6.1 在STATA裡打開數據，初步分析和熟悉數據
##
## . use "backupfiles/RFA.dta"
##
## .
## . summarize age
##
## Variable | Obs Mean Std. Dev. Min Max
## -------------+---------------------------------------------------------
## age | 3,551 51.99944 4.379916 33 84
##
## .
## . tab gender
##
## gender | Freq. Percent Cum.
## ------------+-----------------------------------
## male | 2,124 59.81 59.81
## female | 1,427 40.19 100.00
## ------------+-----------------------------------
## Total | 3,551 100.00
##
## .
## . tab hospital
##
## Hospital ID | Freq. Percent Cum.
## ------------+-----------------------------------
## 1 | 570 16.05 16.05
## 2 | 631 17.77 33.82
## 3 | 1,395 39.28 73.11
## 4 | 955 26.89 100.00
## ------------+-----------------------------------
## Total | 3,551 100.00
##
## .
## . summarize maxdia
##
## Variable | Obs Mean Std. Dev. Min Max
## -------------+---------------------------------------------------------
## maxdia | 3,551 1.816981 .5713256 .7 4
##
## .
## . tab position
##
## Ease with |
## which |
## nodules can |
## be reached | Freq. Percent Cum.
## ------------+-----------------------------------
## easy | 857 24.13 24.13
## moderate | 1,820 51.25 75.39
## difficult | 874 24.61 100.00
## ------------+-----------------------------------
## Total | 3,551 100.00
##
## .
## . tab dodp
##
## Outcome |
## variable: |
## death or |
## disease |
## progression |
## within 36 |
## months | Freq. Percent Cum.
## ------------+-----------------------------------
## no | 2,604 73.33 73.33
## yes | 947 26.67 100.00
## ------------+-----------------------------------
## Total | 3,551 100.0092.6.2 用標準邏輯回歸模型分析 rfa (暴露) 和 dodp (結果) 之間的關係
##
## . use "backupfiles/RFA.dta"
##
## .
## . *(a)
## . logit dodp rfa
##
## Iteration 0: log likelihood = -2059.3462
## Iteration 1: log likelihood = -2030.2625
## Iteration 2: log likelihood = -2030.123
## Iteration 3: log likelihood = -2030.123
##
## Logistic regression Number of obs = 3,551
## LR chi2(1) = 58.45
## Prob > chi2 = 0.0000
## Log likelihood = -2030.123 Pseudo R2 = 0.0142
##
## ------------------------------------------------------------------------------
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## rfa | -.5881255 .0777773 -7.56 0.000 -.7405661 -.4356849
## _cons | -.7496965 .0498312 -15.04 0.000 -.8473639 -.652029
## ------------------------------------------------------------------------------
##
## .
## . *(b)
## . logit dodp rfa i.hospital maxdia i.position
##
## Iteration 0: log likelihood = -2059.3462
## Iteration 1: log likelihood = -1782.1086
## Iteration 2: log likelihood = -1770.8713
## Iteration 3: log likelihood = -1770.8481
## Iteration 4: log likelihood = -1770.8481
##
## Logistic regression Number of obs = 3,551
## LR chi2(7) = 577.00
## Prob > chi2 = 0.0000
## Log likelihood = -1770.8481 Pseudo R2 = 0.1401
##
## ------------------------------------------------------------------------------
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## rfa | -.0215334 .0976595 -0.22 0.825 -.2129425 .1698757
## |
## hospital |
## 2 | -.1087261 .147873 -0.74 0.462 -.3985519 .1810998
## 3 | .4063856 .1222275 3.32 0.001 .1668242 .645947
## 4 | -.1398633 .1340884 -1.04 0.297 -.4026717 .1229451
## |
## maxdia | 1.426368 .0891745 16.00 0.000 1.251589 1.601147
## |
## position |
## moderate | .072014 .1096704 0.66 0.511 -.142936 .286964
## difficult | 1.226249 .1191946 10.29 0.000 .9926321 1.459866
## |
## _cons | -4.234223 .2365292 -17.90 0.000 -4.697812 -3.770634
## ------------------------------------------------------------------------------92.6.3 比較上面(a)和(b)兩個邏輯回歸模型的結果，你認為混雜因素對暴露和結果的關係的影響是怎樣的？
總體來說，混雜的模式 (pattern) 應該是正混雜 (positive confounding)。模型 (a)，不經過任何調整，RFA 似乎比標準手術療法好很多 (幾乎減少一半三年內死亡的對數比值 log-odds)。但是調整了其他混雜因素的模型 (b)，結果暗示兩者之間對於患者預後沒有太大的影響 (treatment effect reduced to suggest very little evidence to departure from no effect)。所以調整的這些變量對暴露和結果之間的關係的混雜是正向的(傾向於把關係改變成為接近另假設 tend to change the association to null)。
92.6.4 在怎樣的前提假設條件下，上面模型 (b) 可以被賦予因果關係的解釋？
這些假設包括：
- 無相互幹擾 no interference：某個病人接受的療法，不影響另一個病人療法的結果。
- 一致性 consistency：對於真的接受了 RFA 療法的病人來說，他/她的療效，和潛在暴露 (potential exposure) 為接受 RFA 療法時的潛在療效 (potential outcome) 是一致的。接受標準手術療法的患者中也是需要一樣的假設。
- 條件可置換性 conditional exchangeability：對於同一所醫院，腫塊大小相同，腫塊位置相同的患者來說，他/她的兩種潛在治療結果 (potential outcome)，和該病人最終到底是接受了常規手術治療，還是接受 RFA 之間是相互獨立的。用更通俗的話說是，暴露變量
rfa和結果變量dodp之間的所有可能的混雜，都被模型中加入的hospital, maxdia, position囊括進去了。 - 正確的模型結構 correct specification of the model：因為模型 (b) 中不包括任何交互作用的相乘項，要給這個模型擬合的回歸係數以因果關係的解釋，我們需要認為模型中的變量之間沒有任何較互作用，也就是說
rfa的療效，不因為醫院，腫塊位置，和腫塊大小不同而不同。
92.6.5 在前面提出的所有前提假設都滿足的情況下，請給模型 (b) 的回歸係數賦予一個因果效應的解釋。
當四個前提假設都可以滿足，模型 (b) 的rfa的回歸係數 -0.022 是一個條件對數比值比 (conditional log-odds ratio)。該條件對數比值比調整了醫院，腫瘤大小，和腫瘤位置。患者在三年內死亡或者疾病加重的對數比值 (log odds) 被估計為 0.022 (95% CI: -0.17, 0.21)。這個對數比值比較的是，情形 A.所有患者都被實施 RFA 手術，和情形 B. 所有患者都被實施常規手術摘除腫塊，兩種潛在情形的潛在結果。用潛在結果的數學語言來解釋就是：
\[ \log\{ \frac{\text{Pr}[Y(1) = 1|\mathbf{C = c}]}{1-\text{Pr}[Y(1) = 1|\mathbf{C = c}]} \} - \log\{ \frac{\text{Pr}[Y(0) = 1|\mathbf{C = c}]}{1-\text{Pr}[Y(0) = 1|\mathbf{C = c}]} \} \]
其中，\(\mathbf{C}=\) {hospital, maxdia, position}
92.6.6 用 STATA 的 teffects ra 擬合上面兩個模型
##
## . use "backupfiles/RFA.dta"
##
## .
## . *(a)
## . teffects ra (dodp, logit) (rfa)
##
## Iteration 0: EE criterion = 5.296e-17
## Iteration 1: EE criterion = 6.380e-32
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATE |
## rfa |
## (radiofre.. |
## vs |
## standard..) | -.113019 .0146495 -7.71 0.000 -.1417316 -.0843064
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .3208874 .0108592 29.55 0.000 .2996039 .342171
## ------------------------------------------------------------------------------
##
## .
## . *(b)
## . teffects ra (dodp i.hospital maxdia i.position, logit) (rfa)
##
## Iteration 0: EE criterion = 6.343e-18
## Iteration 1: EE criterion = 9.535e-33
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATE |
## rfa |
## (radiofre.. |
## vs |
## standard..) | .0261149 .0164826 1.58 0.113 -.0061903 .0584201
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .290968 .0110295 26.38 0.000 .2693506 .3125854
## ------------------------------------------------------------------------------92.6.7 在怎樣的假設前提條件下，前一步擬合的模型 (b) 結果中的 ATE 可以被賦予因果關係的解釋？
這些假設包括：
- 無相互幹擾 no interference，解釋同前。
- 一致性 consistency，解釋同前。
- 條件可置換性 conditional exchangeability，解釋同前。
- 正確的模型結構 correct specification of the mode：調整了醫院，腫塊大小，和腫塊位置以後，患者的死亡或者疾病加重的對數比值 (log odds of death or diseae progression) 和腫塊大小，醫院，腫塊位置不再有任何依賴性(independent)，但是接受 RFA 療法和常規手術療法之間的療效差，被允許在不同的醫院，腫塊大小，以及腫塊位置的不同而有所不同。
92.6.8 前一問和你擬合完簡單的邏輯回歸之後做的模型假設的回答，有什麼不同？
用 STATA 的 teffects ra 命令的時候，我們允許了療效差在不同的醫院，不同腫塊的大小，和不同腫塊的位置之間有所不同。The effect of treatment is allowed to differ nby hospital, position of the nodule and its diameter. 這種不同在擬合簡單邏輯回歸模型時是被忽略掉的。
92.6.9 用因果關係語言解釋 teffects ra 擬合的模型 (b) 的結果
模型 (b) 擬合的結果是邊際因果危險度差 (marginal causal risk difference)。情形 A. 所有患者都被實施 RFA 療法，和情形 B. 所有患者都被實施常規手術療法相比，患者中三年內死亡或者病情加重的比例 (proportion) 被估計要高出 0.026 (95%CI: -0.01, 0.06):
\[ E\{ Y(1) \} - E\{ Y(0) \} \]
92.6.10 如果模型中加入 age, gender, smoke, nodules, mets, duration, primary 等和預後相關但是和決定療法並不太有關係的變量，結果會有什麼不同呢？
##
## . use "backupfiles/RFA.dta"
##
## .
## . teffects ra (dodp age gender i.smoke i.hospital nodules mets duration ///
## > maxdia i.primary i.position, logit) (rfa)
##
## Iteration 0: EE criterion = 1.417e-06
## Iteration 1: EE criterion = 1.527e-07
## Iteration 2: EE criterion = 1.726e-08
## Iteration 3: EE criterion = 3.724e-09
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATE |
## rfa |
## (radiofre.. |
## vs |
## standard..) | .0378445 .0137717 2.75 0.006 .0108525 .0648366
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .2856634 .0100116 28.53 0.000 .266041 .3052858
## ------------------------------------------------------------------------------
##
## .ACE 的估計結果沒有發生非常劇烈的變化，但是，它的標準誤被大大降低了，有效地提高了療效估計的精確度。而且此時的結果已經提示平均因果危險度差是有統計學意義的 (p = 0.006)。這時候，對於整體患者來說，如果全部實施了 RFA，那麼和全部實施標準手術療法相比較會有略差的結果，這個相差是有統計學意義的。
92.6.11 如果再向模型中加入和暴露變量相關，和預後沒什麼關係的變量 coag，結果該怎麼解讀？
##
## . use "backupfiles/RFA.dta"
##
## . teffects ra (dodp age gender i.smoke i.hospital nodules mets duration ///
## > maxdia i.primary i.position coag, logit) (rfa)
##
## Iteration 0: EE criterion = 1.414e-06
## Iteration 1: EE criterion = 1.516e-07
## Iteration 2: EE criterion = 1.711e-08
## Iteration 3: EE criterion = 3.634e-09
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATE |
## rfa |
## (radiofre.. |
## vs |
## standard..) | .0401436 .0151604 2.65 0.008 .0104297 .0698574
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .2815725 .0101806 27.66 0.000 .2616189 .3015261
## ------------------------------------------------------------------------------ACE 的估計量的標準誤因為調整了只和暴露變量相關的變量 coag 變得比之前大了一些。但是此時的結果依然提示全部實施 RFA 療法的話結果會比全部實施常規手術治療要差。這裡應該考慮的是，因為 coag 本身不是暴露和結果變量之間的混雜因子，我們本不該調整這個變量，一旦調整了只和暴露變量相關的變量，我們反而會降低療效估計的精確度，所以不是說模型中想加多少變量就加多少變量的。(Violation of positivity)
另外，已經有證據證明模型中調整了和暴露相關，但是並不是混雜因子的變量會嚴重增加估計量的偏倚 (Pearl 2011)。
92.6.12 使用 atet 的選項重新擬合上面的因果效應模型，解釋結果發生的變化，並作出相應的結論。
##
## . use "backupfiles/RFA.dta"
##
## .
## .
## . teffects ra (dodp, logit) (rfa), atet
##
## Iteration 0: EE criterion = 5.296e-17
## Iteration 1: EE criterion = 2.524e-32
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATET |
## rfa |
## (radiofre.. |
## vs |
## standard..) | -.113019 .0146495 -7.71 0.000 -.1417316 -.0843064
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .3208874 .0108592 29.55 0.000 .2996039 .342171
## ------------------------------------------------------------------------------
##
## .
## . teffects ra (dodp i.hospital maxdia i.position, logit) (rfa), atet
##
## Iteration 0: EE criterion = 6.343e-18
## Iteration 1: EE criterion = 6.099e-33
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATET |
## rfa |
## (radiofre.. |
## vs |
## standard..) | -.0506326 .0164035 -3.09 0.002 -.0827828 -.0184824
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .2585011 .0141292 18.30 0.000 .2308084 .2861938
## ------------------------------------------------------------------------------
##
## .
## . teffects ra (dodp age gender i.smoke i.hospital nodules mets duration ///
## > maxdia i.primary i.position, logit) (rfa), atet
##
## Iteration 0: EE criterion = 1.118e-06
## Iteration 1: EE criterion = 9.719e-08
## Iteration 2: EE criterion = 2.682e-09
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATET |
## rfa |
## (radiofre.. |
## vs |
## standard..) | -.0396473 .0137245 -2.89 0.004 -.0665468 -.0127479
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .2475158 .0130394 18.98 0.000 .221959 .2730726
## ------------------------------------------------------------------------------
##
## .
## .
## .
## . teffects ra (dodp age gender i.smoke i.hospital nodules mets duration ///
## > maxdia i.primary i.position coag, logit) (rfa), atet
##
## Iteration 0: EE criterion = 1.108e-06
## Iteration 1: EE criterion = 9.419e-08
## Iteration 2: EE criterion = 3.077e-09
##
## Treatment-effects estimation Number of obs = 3,551
## Estimator : regression adjustment
## Outcome model : logit
## Treatment model: none
## ------------------------------------------------------------------------------
## | Robust
## dodp | Coef. Std. Err. z P>|z| [95% Conf. Interval]
## -------------+----------------------------------------------------------------
## ATET |
## rfa |
## (radiofre.. |
## vs |
## standard..) | -.031106 .0143989 -2.16 0.031 -.0593273 -.0028847
## -------------+----------------------------------------------------------------
## POmean |
## rfa |
## standard .. | .2389745 .0137425 17.39 0.000 .2120396 .2659094
## ------------------------------------------------------------------------------模型的最後加了 atet 之後最大的不同是導致結果分析的結論完全改變了。這時候 RFA 的療效反而變得有利起來。每個模型都給出了小於零的治療組因果危險度差 (causal risk difference in the treated)。這主要是因為接受 RFA療法和常規手術治療的兩類患者本身有太大的不同 (strong effect heterogeneity)。因為患者如果腫塊尺寸小， RFA 可以給出很好的預後，如果腫塊的位置太難找無法使用標準療法的時候又只能選用RFA。這些都是決定醫生最終給患者使用哪種療法的重要參考因素。所以如果 ATT 是小於零的，意味著，在那些(被)選擇了 RFA療法的患者中，療效是良好的，可以有效的減少三年內死亡和病情加重的概率。
另外值得一提的是，這裡調整了 coag 對ATT 估計的影響較小。這主要是因為，如果一個患者有凝血功能障礙，他/她就不可能接受 RFA 療法:
##
## . use "backupfiles/RFA.dta"
##
## . tab coag rfa, col chi
##
## +-------------------+
## | Key |
## |-------------------|
## | frequency |
## | column percentage |
## +-------------------+
##
## | Treatment variable:
## | RFA or standard
## Coagulopat | surgery
## hy | standard radiofreq | Total
## -----------+----------------------+----------
## no | 1,466 1,701 | 3,167
## | 79.33 99.88 | 89.19
## -----------+----------------------+----------
## yes | 382 2 | 384
## | 20.67 0.12 | 10.81
## -----------+----------------------+----------
## Total | 1,848 1,703 | 3,551
## | 100.00 100.00 | 100.00
##
## Pearson chi2(1) = 388.2057 Pr = 0.000但是一個患者如果沒有凝血功能障礙，那麼他/她接受哪種療法就變得不受 coag 的影響。而此時我們使用 etat 估算的是所有接受 RFA 療法的患者中的療效。這裡受到模型變量共線性影響就較小。
References
Pearl, Judea. 2011. “Invited Commentary: Understanding Bias Amplification.” American Journal of Epidemiology 174 (11). Oxford University Press: 1223–7.