[R] survreg() with frailty

Prof Brian Ripley ripley at stats.ox.ac.uk
Thu Apr 17 15:31:58 CEST 2008


On Thu, 17 Apr 2008, Dimitris Rizopoulos wrote:

> Dear R-users,
>
> I have noticed small discrepencies in the reported estimate of the
> variance of the frailty by the print method for survreg() and the
> 'theta' component included in the object fit:

This is the print method for class "survreg.penal".  It is using

> fit1$printfun[[1]]
function (coef, var, var2, df, history)
{
     if (!is.null(history$history))
         theta <- history$history[nrow(history$history), 1]
     else theta <- history$theta
     clog <- history$c.loglik
     if (is.matrix(var))
         test <- coxph.wtest(var, coef)$test
     else test <- sum(coef^2/var)
     df2 <- max(df, 0.5)
     list(coef = c(NA, NA, NA, test, df, 1 - pchisq(test, df2)),
         history = paste("Variance of random effect=", format(theta),
             "  I-likelihood =", format(round(clog, 1), digits = 10)))
}

that is fit1$history[[1]]$history[10, 1]

Short answer: your fit1 did not converge, and the others only rouighly 
converged.

>
> # Examples in R-2.6.2 for Windows
>
> library(survival) # version 2.34-1 (2008-03-31)
>
> # discrepancy
> fit1 <- survreg(Surv(time, status) ~ rx + frailty(litter), rats)
> fit1
> fit1$history[[1]]$theta
>
> # OK
> fit2 <- survreg(Surv(time, status) ~ rx + frailty(litter, df = 13),
> rats)
> fit2
> fit2$history[[1]]$theta
>
> # discrepancy
> fit3 <- survreg(Surv(time, status)~ age + sex + frailty(id), kidney)
> fit3
> fit3$history[[1]]$theta
>
> # OK
> fit4 <- survreg(Surv(time, status)~ age + frailty(id), kidney)
> fit4
> fit4$history[[1]]$theta
>
>
> Am I missing something? Thanks in advance for any pointers!
>
> Best,
> Dimitris
>
> ----
> Dimitris Rizopoulos
> Biostatistical Centre
> School of Public Health
> Catholic University of Leuven
>
> Address: Kapucijnenvoer 35, Leuven, Belgium
> Tel: +32/(0)16/336899
> Fax: +32/(0)16/337015
> Web: http://med.kuleuven.be/biostat/
>     http://www.student.kuleuven.be/~m0390867/dimitris.htm
>
>
> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595



More information about the R-help mailing list