r - Defining a distribution for survival::survreg() -

i try fit survreg model using gamma distribution.

following ?survreg.distributions defined custom distribution this:

gamma <- list(name = 'gamma',           parms = c(2,2),           init = function(x, weights, ...){             c(median(x), mad(x))           },           density = function(x, parms){             shape <- parms[1]             scale <- parms[2]             cbind(pgamma(x, shape=shape, scale=scale),                   1-pgamma(x, shape=shape, scale=scale),                   dgamma(x, shape=shape, scale=scale),                   (shape-1)/x - 1/scale,                   (shape-1)*(shape-2)/x^2 - 2*(shape-1)/(x*scale) + 1/scale^2)           },           quantile = function(p, parms) {             qgamma(p, shape=parms[1], scale=parms[2])           },           deviance = function(...) stop('deviance residuals not defined') ) 

however can't run:

require(survival) survreg(surv(log(time), status) ~ ph.ecog + sex, lung, dist=gamma) #error in coxph.wtest(t(x) %*% (wt * x), c((wt * eta + weights * deriv$dg) %*%  :  #  na/nan/inf in foreign function call (arg 3) 

the error comes c-code think generated earlier...

any hints/suggestions/alternatives survreg?

i found flexsurv package implements generalized gamma distribution.

for weibull distribution estimates survreg , flexsurvreg similar (but note different parametrization:

require(survival) summary(survreg(surv(log(time), status) ~ ph.ecog + sex, data = lung, dist='weibull'))  call: survreg(formula = surv(log(time), status) ~ ph.ecog + sex, data = lung,      dist = "weibull")               value std. error      z         p (intercept)  1.7504     0.0364  48.13  0.00e+00 ph.ecog     -0.0660     0.0158  -4.17  3.10e-05 sex          0.0763     0.0237   3.22  1.27e-03 log(scale)  -1.9670     0.0639 -30.77 6.36e-208  scale= 0.14   weibull distribution loglik(model)= -270.5   loglik(intercept only)= -284.3     chisq= 27.62 on 2 degrees of freedom, p= 1e-06  number of newton-raphson iterations: 6  n=227 (1 observation deleted due missingness)   require(flexsurv) flexsurvreg(surv(log(time), status) ~ ph.ecog + sex, data = lung, dist='weibull')  call: flexsurvreg(formula = surv(log(time), status) ~ ph.ecog + sex,     data = lung, dist = "weibull")  maximum likelihood estimates:              est    l95%    u95% shape    7.1500  6.3100  8.1000 scale    5.7600  5.3600  6.1800 ph.ecog -0.0660 -0.0970 -0.0349 sex      0.0763  0.0299  0.1230  n = 227,  events: 164,  censored: 63 total time @ risk: 1232.1 log-likelihood = -270.5, df = 4 aic = 549 

with flexsurvreg can fit generalized gamma distribution data:

flexsurvreg(surv(log(time), status) ~ ph.ecog + sex, data = lung, dist='gengamma')  call: flexsurvreg(formula = surv(log(time), status) ~ ph.ecog + sex,     data = lung, dist = "gengamma")  maximum likelihood estimates:              est    l95%    u95% mu       1.7800  1.7100  1.8600 sigma    0.1180  0.0971  0.1440 q        1.4600  1.0200  1.9100 ph.ecog -0.0559 -0.0853 -0.0266 sex      0.0621  0.0178  0.1060  n = 227,  events: 164,  censored: 63 total time @ risk: 1232.1 log-likelihood = -267.57, df = 5 aic = 545.15 

the loglogistic distribution (in contrast survreg) not build in, can custumized (see examples of flexsurvreg).

i haven't tested much, flexsurv seems alternative survival.