survreg penalized likelihood?
On Sat, 19 Apr 2003, Spencer Graves wrote:
What objective function is maximized by survreg with the default Weibull model? I'm getting finite parameters in a case that has the likelihood maximzed at Infinite, so it can't be a simple maximum likelihood.
The objective function *is* the loglikelihood.
Consider the following: #############################
> set.seed(3) > Stress <- rep(1:3, each=3) > ch.life <- exp(9-3*Stress) > simLife <- rexp(9, rate=1/ch.life) > Data <- data.frame(Stress=factor(Stress),
+ Time=pmin(simLife, 50), dead = (simLife<=50))
> Data[1:3,]
Stress Time dead 1 1 50 FALSE 2 1 50 FALSE 3 1 50 FALSE # All three observations at Stress=1 are censored.
> Fit <- survreg(Surv(Time, dead)~Stress-1, data=Data) > Fit
Call:
survreg(formula = Surv(Time, dead) ~ Stress - 1, data = Data)
Coefficients:
Stress1 Stress2 Stress3
41.3724178 2.2819204 -0.5281005
Scale= 0.2577886
Loglik(model)= -6.6 Loglik(intercept only)= -22.1
Chisq= 30.88 on 2 degrees of freedom, p= 2e-07
n= 9
# Even though 100% of observations at Stress=1 are censored,
# I still get a finite estimate for log(characteristics life).
I suppose technically 41.3 is a finite estimate of log(life), but since your data are censored at 50, a life expectancy of nearly 10^18 isn't terribly finite. The likelihood at this value is very very very close to the likelihood at the true mle, and that's all a numerical optimisation technique can really be expected to give you (and all that statistical theory, frequentist or Bayesian, demands). It's quite difficult to come up with a reliable test for ridges in the loglikelihood -- coxph() tries, but gives too many false positives. Incidentally, if you relax the assumption that the scale parameter is the same for each Stress you end up with much larger finite predictions for Stress=1, as the scale parameter ends up around 10^15 in that stratum. -thomas