Michael McCulloch <[email protected]> posted a question regarding a perceived
contradiction between a p-value and 95% CI from -bootstrap- of results from
the -stcox- command.
Michael, Philip Ryan <[email protected]>, and Richard Williams
<[email protected]> discussed the issue over several emails.
Before I offer my thoughts on this discussion, I note that there is nothing
wrong with -bootstrap- or -stcox-.
The perceived contradiction was due to the incorrect assumption that
-bootstrap- was reporting a p-value for the test of a hazard ratio against
one. The p-values in the tables reported by Stata's estimation commands,
including -bootstrap-, are always for tests against zero. -stcox- reports
exponentiated coefficients (hazard ratios) by default, but the reported
tests are of the null hypotheses that the true underlying coefficients are
zero. The default (hazard ratio) confidence intervals (CIs) for -stcox-
are obtained by exponentiating the limits of the CIs for the underlying
coefficient.
Now I'll summarize the thread and add a few of my own comments.
Michael was using the -bootstrap- prefix with the following user-written
program:
program boot_hr, rclass
stcox drug
indeplist, local
foreach var of varlist `X' {
return scalar `var' = exp(_b[`var'])
}
end
This program assumes the dataset is properly -stset- and the user-written
-indeplist- program has been installed (see -ssc describe indeplist-).
Michael shows an example of how he is using -bootstrap- and the above program:
. sysuse cancer
. stset studytime, failure(died)
. set seed 12358
. bootstrap drug=r(drug), reps(100) : boot_hr
The above can be done without using -boot_hr-, here's how
. bootstrap drug=exp(_b[drug]), reps(100) seed(12358) : stcox drug
Michael later pointed out his reason for using -boot_hr-. He needs
-bootstrap- to ignore the fact that -iweight-s were -stset-. Thus Michael
is using -boot_hr- to prevent -bootstrap- from figuring out that weights are
being used by the estimation command.
This assumes that resampling the observations and performing
estimation using the weights associated with the resampled
observations is appropriate. While this is reasonable for -aweight-s
and arguably reasonable for -pweight-s, I would strongly discourage
anyone from trying this with -fweight-s (instead consider expanding
the data and dropping the weights).
Michael indicated that his weights were the inverse probability of exposure,
something one might consider as -iweight-s as Michael did
The -force- option of -bootstrap- accomplishes the above, but I must
strongly warn that it is the user's responsibility to verify that
resampling observations and using the weights associated with the
resampled observations is appropriate.
When bootstrapping the coefficients -e(b)- from an estimation command,
-bootstrap- is able to properly label the reported values in its table.
-bootstrap- uses the estimation command to report its results when it can.
Thus something like
. bootstrap, reps(100) seed(12358) : stcox drug
will report hazard ratios, because -stcox- is being used to report the
results. However using -exp_list- or a user-written program that cannot
replay results, as in
. bootstrap drug=r(drug), reps(100) : boot_hr
implies that -bootstrap- must report its own results and labels the
bootstrapped point estimates using the standard 'Coef.' column header. I
believe this is what ultimately contributed to the confusion.
--Jeff
[email protected]
*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
