Here are the model outputs.
Nathan
Survey package
ca.ATE.design <- svydesign(ids = ~ id, weights = ~ get.weights(ca.ATE.ps,
stop.method = 'ks.mean'), data = ca.dt)
Independent Sampling design (with replacement)
svydesign(ids = ~id, weights = ~get.weights(ca.ATE.ps, stop.method
"ks.mean"),
data = ca.dt)
> ca.ATE.dexmg.svy
Call:
svycoxph(formula = Surv(daysfromsurgerytodeath, as.logical(deceased)) ~
dexamethasonemg + paincontrol + histgrade + adjuvant + stage +
anesthetictransfusionunits, design = ca.ATE.design)
coef exp(coef) se(coef) z p
dexamethasonemg -0.0863 0.917 0.0339 -2.550 1.1e-02
paincontrolNot Epidural 0.6027 1.827 0.2370 2.543 1.1e-02
histgradeg2 0.9340 2.545 0.4307 2.168 3.0e-02
histgradeg3 1.2749 3.578 0.4453 2.863 4.2e-03
adjuvantyes -0.5810 0.559 0.2529 -2.298 2.2e-02
stageib -0.4394 0.644 0.6056 -0.726 4.7e-01
stageiia 1.6565 5.241 0.5193 3.190 1.4e-03
stageiib 1.6928 5.435 0.4902 3.453 5.5e-04
stageiii 1.8211 6.179 0.5130 3.550 3.9e-04
stageiv 2.3251 10.227 0.6940 3.350 8.1e-04
anesthetictransfusionunits 0.1963 1.217 0.0400 4.908 9.2e-07
Likelihood ratio test= on 11 df, p= n= 144, number of events= 102
rms package
> ca.ATE.dexmg.rms2
Cox Proportional Hazards Model
cph(formula = Surv(daysfromsurgerytodeath, as.logical(deceased)) ~
dexamethasonemg + paincontrol + histgrade + adjuvant + stage +
anesthetictransfusionunits + cluster(id), data = ca.dt,
weights = get.weights(ca.ATE.ps, stop.method = "ks.mean"),
robust = T, x = T, y = T, se.fit = T, surv = T, time.inc = 30)
Model Tests Discrimination
Indexes
Obs 144 LR chi2 117.80 R2 0.559
Events 102 d.f. 11 Dxy -0.459
Center 2.4016 Pr(> chi2) 0.0000 g 1.083
Score chi2 122.57 gr 2.953
Pr(> chi2) 0.0000
Coef S.E. Wald Z Pr(>|Z|)
dexamethasonemg -0.0863 0.0192 -4.49 <0.0001
paincontrol=Not Epidural 0.6027 0.1203 5.01 <0.0001
histgrade=g2 0.9340 0.2209 4.23 <0.0001
histgrade=g3 1.2749 0.2612 4.88 <0.0001
adjuvant=yes -0.5810 0.1741 -3.34 0.0008
stage=ib -0.4394 0.1899 -2.31 0.0207
stage=iia 1.6565 0.2097 7.90 <0.0001
stage=iib 1.6928 0.1979 8.55 <0.0001
stage=iii 1.8211 0.2411 7.55 <0.0001
stage=iv 2.3251 0.1886 12.33 <0.0001
anesthetictransfusionunits 0.1964 0.0214 9.17 <0.0001
From: Thomas Lumley <tlumley at uw.edu>
Date: Tuesday, February 25, 2014 at 3:09 PM
To: "Nathan Leon Pace, MD, MStat" <n.l.pace at utah.edu>
Cc: r help list <r-help at r-project.org>
Subject: Re: [R] SEs rms cph vs survey svycoxph
On Tue, Feb 25, 2014 at 2:51 PM, Nathan Pace
<n.l.pace at utah.edu> wrote:
I?ve used twang to get ATE propensity scores.
I?ve done multivariable, case weighted Cox PH models in survey using
svycoxph and in rms using cph with id(cluster) set to get robust estimates.
The model language is identical.
The point estimates are identical, but the CIs are considerably wider with
svycoxph estimates.
There is a note in the svycoxph help page stating the SEs should agree
closely unless the model fits poorly.
The actual note on the svycoxph help page says
"The standard errors agree closely with survfit.coxph for independent
sampling when the model fits well, but are larger when the model fits
poorly. "
That is, the note is for the survival curve rather than the coefficients.
It's still surprising that there's a big difference, but I think we need
more information.
-thomas
--
Thomas Lumley
Professor of Biostatistics
University of Auckland