survreg penalized likelihood?

	  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.

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).

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

Thanks very much.

After I distributed that example, I noticed that the documtation for the
survival package says, "Survival analysis, including penalised
likelihood", so I thought that some penalty function might be invoked
under certain circumstances.