Statalist


[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

Re: st: stcox output: p-value and CI don't agree


From   jpitblado@stata.com (Jeff Pitblado, StataCorp LP)
To   statalist@hsphsun2.harvard.edu
Subject   Re: st: stcox output: p-value and CI don't agree
Date   Wed, 08 Aug 2007 16:09:22 -0500

Michael McCulloch <mm@pinest.org> 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 <philip.ryan@adelaide.edu.au>, and Richard Williams
<Richard.A.Williams.5@ND.edu> 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
jpitblado@stata.com
*
*   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/



© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index