######## Week 4, RQ1 ## propensity scores fit and ATE weights formed in Week 3 CC > m2full = matchit(abcix ~ stent + height + female + diabetic + acutemi + ejecfrac + ves1proc, data = lindner, method = "full") > # weights for estimating ATE (according to lore) > # for abcix = 1, 1/(m2full.dat$distance) > # for abcix = 0, 1/(1 - m2full.dat$distance) > m2full.dat$weight.ATE = ifelse(m2full.dat$abcix == 1, 1/m2full.dat$distance, 1/(1 - m2full.dat$distance)) ################################################################################################## > ## IPTW ATE for lifepres outcome > table(lifepres) lifepres 0 11.6 26 970 > str(m2full.dat) 'data.frame': 996 obs. of 16 variables: $ lifepres : num 0 11.6 11.6 11.6 11.6 11.6 11.6 11.6 11.6 11.6 ... $ cardbill : int 14301 3563 4694 7366 8247 8319 8410 8517 8763 8823 ... $ abcix : int 1 1 1 1 1 1 1 1 1 1 ... $ stent : int 0 0 0 0 0 0 0 0 0 0 ... $ height : int 163 168 188 175 168 178 185 173 152 180 ... $ female : int 1 0 0 0 1 0 0 1 1 0 ... $ diabetic : int 1 0 0 1 0 0 0 0 0 0 ... $ acutemi : int 0 0 0 0 0 0 0 0 0 0 ... $ ejecfrac : int 56 56 50 50 55 50 58 30 60 60 ... $ ves1proc : int 1 1 1 1 1 1 1 1 1 1 ... $ distance : num 0.408 0.578 0.524 0.473 0.493 ... $ weights : num 1 1 1 1 1 1 1 1 1 1 ... $ subclass : num 1 6 85 140 176 209 229 238 244 66 ... $ weight.ATE : num 2.45 1.73 1.91 2.12 2.03 ... $ rpropen : num 0.622 0.622 0.622 0.622 0.622 ... $ boost_propen: num 0.695 0.695 0.695 0.695 0.695 ... # make a 0,1 outcome, a factor should be ok also > table(lifepres/11.6) 0 1 26 970 > m2full.dat$life = m2full.dat$lifepres/11.6 > table(m2full.dat$life) 0 1 26 970 # still get warning cuz of weights > glm.ATE = glm(life ~ abcix, family = binomial, data = m2full.dat, weights = (weight.ATE)) Warning message: In eval(expr, envir, enclos) : non-integer #successes in a binomial glm! > summary(glm.ATE) Call: glm(formula = life ~ abcix, family = binomial, data = m2full.dat, weights = (weight.ATE)) Deviance Residuals: Min 1Q Median 3Q Max -10.9754 0.1951 0.2151 0.5773 1.7473 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) 2.4250 0.1154 21.021 < 2e-16 *** abcix 1.7471 0.2837 6.158 7.36e-10 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 775.82 on 995 degrees of freedom Residual deviance: 723.82 on 994 degrees of freedom AIC: 707.97 Number of Fisher Scoring iterations: 6 > exp(coef(glm.ATE)) (Intercept) abcix 11.302532 5.737676 ## From glmer using fullmatch subclasses for lifepres, much larger effect Fixed effects: Estimate Std. Error z value Pr(>|z|) (Intercept) 2.9389 0.2659 11.053 < 2e-16 *** abcix 6.0279 0.8381 7.192 6.39e-13 *** ---------------------------------------------------------------- # compare glm without weighting > glm.ATE_noW = glm(life ~ abcix, family = binomial, data = m2full.dat) > exp(coef(glm.ATE_noW)) (Intercept) abcix 18.866667 3.310312 > summary(glm.ATE_noW)) Error: unexpected ')' in "summary(glm.ATE_noW))" > summary(glm.ATE_noW) Call: glm(formula = life ~ abcix, family = binomial, data = m2full.dat) Deviance Residuals: Min 1Q Median 3Q Max -2.8811 0.1782 0.1782 0.3214 0.3214 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) 2.9374 0.2650 11.086 < 2e-16 *** abcix 1.1970 0.4032 2.969 0.00299 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 240.89 on 995 degrees of freedom Residual deviance: 232.04 on 994 degrees of freedom AIC: 236.04 Number of Fisher Scoring iterations: 7 > ?glm