[R] Survreg with gamma distribution
Thomas Lumley
tlumley at u.washington.edu
Mon Oct 18 19:37:45 CEST 2004
On Tue, 19 Oct 2004, Kuan-Ta Chen wrote:
>
> One million thanks to Prof. Ripley and Prof. Lumley. I think I now have more
> understanding regarding survreg with gamma distribution. But one of my
> problems is still there: in the text of Lee, Wang (2003), there are two
> "kinds" of parametric fitting: 1) fitting of survival distributions (like
> regular probabillity distribution fitting) 2) regression model fitting
> (mostly assume an accelerated failure time model). Survreg {survival}
> provides model fitting of (2). But I still have one problem regarding (1):
> try to estimate the parameters of gamma distributions for some data.
There aren't really two separate kinds: 1 is a special case of 2, so
survreg() can do 1.
> For regular gamma distr. fitting, we could use fitdistr (mass) or use
> optim()/mle() with log-likelihood composed by dgamma()/pgamma(). But because
> the data contains (randomly) censored observations, the log-likelihood
> function must be modified to include the effect of duration of censored
> observations.
Yes. The loglikelihood is
pgamma(x,shape,scale=scale,lower.tail=FALSE,log.p=TRUE)
for a censored observation and
dgamma(x,shape,scale=scale,log=TRUE)
for an uncensored observation. No integration necessary.
You might want to work with log(shape) and log(scale) instead, to avoid
the boundaries at 0.
eg if your data were in variables "times" and "status"
ll <-function(logshape,logscale){
-sum(ifelse(status,
pgamma(times,exp(logshape),scale=exp(logscale),
lower.tail=FALSE,log.p=TRUE),
dgamma(times,exp(logshape),scale=exp(logscale),log=TRUE)
))
}
This works in mle() without too much sensitivity to starting values.
-thomas
More information about the R-help
mailing list