[R] Error survreg: Density function returned an an invalid matrix
Israel Ortiz
isra4884 at gmail.com
Fri Nov 13 16:17:48 CET 2015
Thanks Terry but the error persists. See:
> library(foreign)> library(survival)> library(VGAM) > mypareto <- list(name='Pareto',+ init= function(x, weights,parms){+ alpha <- length(x)/(sum(log(x)))#this is a MLE for alpha+ c(media <-(alpha*1/(alpha-1)),varianza <- ((1/alpha)^2)*(alpha/(alpha-2)))},+ density= function(x,weights) {+ alpha <- length(x)/(sum(log(x)))+ cdf1 <- function(x, alpha) ifelse(x > 1 , 1 - (1/x)**alpha, 0 )+ cdf2 <- function(x, alpha) ifelse(x > 1, (1/x)**alpha ,0)+ distribution <- function(x, alpha) ifelse(x > 1 , alpha/(x**(alpha+1)), 0)+ firstdev <- function(x, alpha) ifelse(x > 1, -(alpha+x)/x, 0)+ seconddev <- function(x, alpha) ifelse(x > 1, (alpha+1)*(alpha+2)/x^2,0)+ cbind(cdf1(x,alpha),cdf2(x, alpha), distribution(x,alpha),firstdev(x,alpha),seconddev(x,alpha))},+ deviance=function(x) {stop('deviance residuals not defined')},+ quantile= function(p, alpha) ifelse(p < 0 | p > 1, NaN, 1*(1-p)**(-1/alpha)))> > survregDtest(mypareto, TRUE)[1] TRUE> set.seed(1)> a <- rpareto(100, 1, 1) > b <- rnorm(100,5,1)> c <- rep(1,100)> base <- cbind.data.frame(a,b,c)> mod1<-survreg(Surv(a, c) ~ b, base, dist = mypareto)Error in survreg.fit(X, Y, weights, offset, init = init, controlvals = control, :
Density function returned an invalid matrix
2015-11-04 7:52 GMT-06:00 Therneau, Terry M., Ph.D. <therneau at mayo.edu>:
> Hi, I want to perform a survival analysis using survreg procedure from
>> survival library in R for a pareto distribution for a time variable, so I
>> set the new distribution using the following sintax:
>>
>> library(foreign)
>> library(survival)
>> library(VGAM)
>>
>> mypareto <- list(name='Pareto',
>> init= function(x, weights,parms){
>>
> etc.
>
> The survreg routine fits location-scale distributions such that (t(y) -
> Xb)/s ~ F, where t is an optional transformation, F is some fixed
> distribution and X is a matrix of covariates. For any distribution the
> questions to ask before trying to add the distribution to survreg are
> - can it be written in a location-scale form?
> - if so, how do the parameters of the distribution map to the location
> (Xb) and scale (s).
>
> In fitting data we normally have per-subject location (X b) but an
> intercept-only model is of course possible.
>
> If y is Weibull then log(y) fits into the framework, which is how survreg
> fits it. The transformation of parameters location and scale parameters
> for log(y) back to the usual Weibull parameterization for y often trips
> people up (see comments in the Examples section of ?survreg).
>
> The log of a Pareto can be written in this form (I think?). The two
> parameters are the scale a and lower limit b, with survival function of
> S(x)= (b/x)^a, for x >= b. If y = log(x) the survival function for y is
> S(y) = (b/exp(y))^a = exp[-(y - log(b))/(1/a)], which has location log(b)
> and scale 1/a. But even if I am correct the discontinuity at b will cause
> the underlying Newton-Raphson method to fail.
>
> Terry Therneau
>
[[alternative HTML version deleted]]
More information about the R-help
mailing list