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

From |
Jeph Herrin <junk@spandrel.net> |

To |
statalist@hsphsun2.harvard.edu |

Subject |
Re: st: repeat: predict after -streg- |

Date |
Sat, 14 Jul 2007 12:52:22 -0400 |

Steven, Thanks much for your careful reply. I was trying to get probabilities for each patient, which was going nowhere. cheers, Jeph Steven Samuels wrote:

LINE 0 Here is code that does what Jeph asks. On Jul 13, 2007, at 2:45 PM, Jeph Herrin wrote:Use the formula P{time<30} = P{ Z< log(30)-mu)/sigma) where Z has a normal distirbution and mu and sigma are the mean and variance for log(time)I can rank on the hospital effects, but how do I convert the model parameters to estimates of the probability of 30-day readmission?

/************************CODE BEGINS**********************************************************/(log survival time scale or log hazard scale). Check on the agreement of the rankings from the different specifications. You may be able to use BIC or some other criterion to choose a "best" model.

I did this initially, comparing LLs, AICs, and plots of Cox-Snell residuals. The log-normal model carried the day.

/* MAIL MAY ADD GREMLINS TO THIS TEXT. ZAP THEM IN YOUR EDITOR BEFORE RUNNING */

set more off

capture log close

log using pcox02, replace text

use http://www.stata-press.com/data/r9/catheter ,clear

stset time, fail(infect)

sum age, det

tab female, missing

//Generate fake hospital id

gen hosp_id=113

replace hosp_id=217 if _n>25

replace hosp_id=66 if _n>50

tab hosp_id

/* Generate a unique code numbered 1,2,3,... . Before this, delete any hospitals which will not appear in the regression */

egen hcode= group(hosp_id) //Give hospitals a temporary sequential code

sum hcode

local n_hosps= `r(max)' //Count Hospitals

di `n_hosps'

tab hcode, gen(xh) // Dummy variables for regression

sort hcode

tempfile t1

save `t1'

gen age45=age -45 //center age so that zero is for a real number

local covariates= "age45 female" //Covariate list

local ncovs: word count `covariates'

local n_predictors =`ncovs'+`n_hosps' //Needed to pick out the sigma

streg xh* `covariates', distribution(lnormal) nocons

matrix beta =e(b)

matrix list beta

di el(beta,1,`n_predictors'+1)

gen sigma=exp(el(beta,1,`n_predictors'+1)) //sigma for the normal distribution

di sigma

forvalues i=1/`n_hosps'{

di `i'

gen mu_`i'= el(beta,1,`i') // mean for hosp_idital i when other covariates =0

gen prob30`i'=normal((log(30)-mu_`i')/sigma) //Prob{T<=30} for hosp_idital i

di `i' " " round(prob30`i',.01)

}

// Now merge the probabilities back with the original hosp_id numbers

gen dum=1

keep in 1

keep dum prob30*

reshape long prob30, i(dum) j(hcode)

label var prob30 "Est Prob Readmit in <=30 days}"

drop dum

sort hcode

merge hcode using `t1'

bysort hosp_id: keep if _n==1

list hosp_id prob30

/**********************************CODE ENDS****************************************/

*

* 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/

* * 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/

**References**:**Re: st: repeat: predict after -streg-***From:*wgould@stata.com (William Gould, StataCorp LP)

**Re: st: repeat: predict after -streg-***From:*Jeph Herrin <junk@spandrel.net>

**Re: st: repeat: predict after -streg-***From:*Steven Samuels <ssamuels@albany.edu>

**Re: st: repeat: predict after -streg-***From:*Jeph Herrin <junk@spandrel.net>

**Re: st: repeat: predict after -streg-***From:*Steven Samuels <ssamuels@albany.edu>

- Prev by Date:
**st: repeat: predict after -streg-** - Next by Date:
**st: West Coast Stata Users Group meeting** - Previous by thread:
**Re: st: repeat: predict after -streg-** - Next by thread:
**Re: st: repeat: predict after -streg-** - Index(es):

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