Re: [R] Gaussian frailty leads to segmentation fault

From: Thomas Lumley <tlumley_at_u.washington.edu>
Date: Fri 30 Jul 2004 - 08:53:59 EST

The crash happens when the penalty on the frailty term gets so small that the model is singular and there are no degrees of freedom left for the treatment effect. The attached revision for survival/R/coxph.wtest.s prevents the crash, but the model is still singular.

If you use method="aic" instead of the default method="reml" the model converges. It is interesting to note that the estimated frailty variance is then very small, so the two methods of estimating the frailty variance disagree completely on these data.

This may well indicate that the penalized likelihood approach to the frailty should be approached with strong distrust in these data.

        -thomas

On Thu, 29 Jul 2004, Christian Lederer wrote:

> Dear Thomas,
>
> attached you find a data frame which produces the error.
> I am using survival 2.11-5 under R 1.9.1-1 and 1.9.0-1.
>
> By the way, if i randomly omit 50% of the data, i usually
> get no crash, but a warning message like this:
> Inner loop failed to coverge for iterations 1 2 3 in: coxpenal.fit(X,
> Y, strats, offset, init = init, control, weights = weights,
>
> Maybe, the model is not appropriate for this kind of data.
> But on the other hand, as soon the treatment group (study == 1,
> treatment == 1) is smaller than the randomized placebo group
> (study == 1, treatment == 0), the warnings disappear.
> and the model gives reasonable results in my first simulations
> with normally distributed study effects.
>
> Christian
>
>
>
>
> Thomas Lumley wrote:
> > We really need a reproducible example to find segmentation faults. Can
> > you make one?
> >
> > -thomas
> >
> >
> > On Wed, 28 Jul 2004, Christian Lederer wrote:
> >
> >
> >>Dear R gurus,
> >>
> >>for a simulation concerning study effects and historical controls
> >>in survival analysis, i would like to experiment with a gaussian
> >>frailty model.
> >>
> >>The simulated scenario consists of a randomized trial
> >>(treatment and placebo) and historical controls (only placebo).
> >>
> >>So the simulated data frames consist of four columns
> >>$time, $cens, $study, $treat.
> >>$time, $cens are the usual survival data.
> >>For the binary thretment indicator we have
> >>$treat == 0 or 1, if $study == 1,
> >>$treat == 1 if $study > 1
> >>
> >>Typical parameters for my simulations are:
> >>sample sizes (per arm): between 100 and 200
> >>number of historical studies: between 7 and 15
> >>hazard ratio treatment/placebo: between 0.7 and 1
> >>variance of the study effekt: between 0 and 0.3
> >>
> >>Depending on the sample sizes, the following call sometimes leads to
> >>a segmentation fault:
> >>
> >>coxph(Surv(time,cens) ~
> >> as.factor(treatment) + frailty(study, distribution="gaussian"),
> >> data=data)
> >>
> >>I noticed, that this segmentation fault occures most frequently, if the
> >>number of randomized treatment patients is higher than the number of
> >>randomized placebo patients, and the number of historical studies is
> >>large.
> >>There seems to be no problem, if there are at least as many randomized
> >>placebo patients as treated patients. Unfortunately, this is not the
> >>situation i want to investigate (historical controls should be used
> >>to decrease the number of treated patients).
> >>
> >>Is there a way to circumwent this problem?
> >>
> >>Christian
> >>
> >>P.S.
> >>Is it allowed, to attach gzipped sample data sets in this mailing list?
> >>
> >>______________________________________________
> >>R-help@stat.math.ethz.ch mailing list
> >>https://www.stat.math.ethz.ch/mailman/listinfo/r-help
> >>PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
> >>
> >
> >
> > Thomas Lumley Assoc. Professor, Biostatistics
> > tlumley@u.washington.edu University of Washington, Seattle
> >
> > ______________________________________________
> > R-help@stat.math.ethz.ch mailing list
> > https://www.stat.math.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
> >
> >
>
>

Thomas Lumley			Assoc. Professor, Biostatistics
tlumley@u.washington.edu	University of Washington, Seattle

______________________________________________

R-help@stat.math.ethz.ch mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Received on Fri Jul 30 09:01:13 2004

This archive was generated by hypermail 2.1.8 : Wed 03 Nov 2004 - 22:55:22 EST