[R] Survreg with gamma distribution

Thomas Lumley tlumley at u.washington.edu
Mon Oct 18 17:48:40 CEST 2004

On Sun, 17 Oct 2004, Kuan-Ta Chen wrote:

> Hi, all:
> I find survreg {survival} has provided many distributions such as weibull,
> lognormal, etc. But I wonder why it doesn't have the support for gamma
> distribution since it should be a good distr. in lifetime analysis. Can
> anybody figure out the reason?

Presumably the actual reason is because Terry Therneau didn't need to use 
the Gamma model.  However, all the distributions in survreg are 
location-scale families, which the Gamma is not, so the basic algorithm 
would have to be different.

> I've tried to implement the likelihood function of progressively censored
> data for gamma distr. and use optim() to solve the paramemters. The
> log-likelihood function L contains some integrations.

It shouldn't have to: we do have pgamma() built in (and digamma, trigamma, 
etc for derivatives).

>							I use tryCatch() to
> capture the error when integration lead to divergence and return Inf.
> But if consequent two calls to the objective function return Inf, optim()
> will raise errors:
> Error in optim(c(ga, 1/la), fr, method = "BFGS") :
>        non-finite finite-difference value [1]
> What can I do except for choosing better initial values?

Choose better initial values.  You should be able to get quite good 
initial values for regression coefficients by using survreg on a lognormal 
distribution, since Gamma and lognormal distributions agree pretty well 
except in the extreme tails.  You could then try getting the shape 
parameter by matching the variance of the Gamma to the variance of the 
fitted lognormal.

> The last question, by its name "survreg", survreg does its job by
> regression,
> but why p.75 in Tableman, Kim (2004) said that "We use the S function
> survReg to fit parametric models (with the MLE approach)...". Does survreg
> use regression or MLE approach?

Its job *is* regression. It uses maximum likelihood to fit a regression 


More information about the R-help mailing list