[R] survreg.distributions() error
    Nick Reich 
    nick.reich at gmail.com
       
    Mon Jun  1 05:47:37 CEST 2009
    
    
  
Hi there.
I am receiving an unexpected error message when creating a new
distribution for the survreg() function in the survival package.  I
understand the survival.distributions() function and have been
following the Cauchy example provided in the help file.
My goal is to use survreg to fit a gamma distribution to interval
censored data.
Here is a simple example of what I'm trying to do.  This code
recreates the error on my machine.
require(survival)
## create the gamma distribution for survreg()
## paramter vector is parms = c(shape, scale)
mygamma <- list(name='gamma',
		init = function(x, weights, ...) {
			c(median(x), mad(x)) },
		density = function(x, parms) {
			shape <- parms[1]  ## "k" in wikipedia definition
			scale <- parms[2]  ## "theta" in wikipedia
			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=parms2)},
		deviance= function(...) stop('deviance residuals not defined'),
		parms=c(2,2)
		)
## create some interval censored data
xx <- Surv(time=rep(c(.5, 1, 1.5), each=5),
	   time2=rep(2, 15),
	   event=rep(3, 15),
	   type="interval")
## lognormal model runs fine
survreg(xx~1, dist="lognormal")
## mygamma model gives an error
survreg(xx~1, dist=mygamma)
This is the error message that I get:
Error in density(z2) : argument "parms" is missing, with no default
I've tracked this down to a single line in the survreg.fit() function.
 This function defines another function derfun() whose first few lines
look like this:
derfun <- function(y, eta, sigma, density, parms) {
        ny <- ncol(y)
        status <- y[, ny]
        z <- (y[, 1] - eta)/sigma
        dmat <- density(z, parms)  ## line A: THIS COMMAND DOES NOT
GIVE AN ERROR
        dtemp <- dmat[, 3] * dmat[, 4]
        if (any(status == 3)) {
            z2 <- (y[, 2] - eta)/sigma
            dmat2 <- density(z2)    ## line B: THIS IS THE COMMAND
THAT GIVES THE ERROR
        }
        else {
            dmat2 <- matrix(0, 1, 5)
            z2 <- 0
        }
>From playing with debug() a little bit, it looks to me like if line B
(see labels above) called "density(z2, parms)", similarly to line A,
then the error message might be avoided.  Instead, it calls
"density(z2)".  Maybe there's a fundamental error in my setting up of
the gamma distribution in the first place, but this seems like
possibly just as simple as a typo mistake.
So...  two questions:
1.  Can anyone verify that this is indeed a bug that's as simple as a
typo in the package code?
2.  I tried cutting and pasting the code for survreg.fit() into a new
file, making the edit myself and then re-sourcing the function code to
overwrite the "bad" function.  But this isn't working.  Is there a
simple way to test out a new version of this function?
Thanks,
Nick
    
    
More information about the R-help
mailing list