Mark Wardle <mark <at> wardle.org> writes: > I'm having difficulty getting access to data generated by survfit and > print.survfit when they are using with a Cox model (survfit.coxph). > > I would like to programmatically access the median survival time for > each strata together with the 95% confidence interval. I can get it on > screen, but can't get to it algorithmically. I found myself examining > the source of print.survfit to try to work out how it is done > internally, but is there a better way? > > An example (and I realise that estimating survival curses from an > average status and time is incorrect in this instance, but it keeps > this example simple): > > test1 <- list(time= c(4, 3,1,1,2,2,3), > status=c(1,NA,1,0,1,1,0), > x= c(0, 2,1,1,1,0,0), > sex= c(0, 0,0,0,1,1,1)) > c1 <- coxph( Surv(time, status) ~ x + strata(sex), test1) #stratified model > > f1 <- survfit(c1) > sf1 <- summary(f1) > str(f1) > print(f1) > print(sf1) > str(sf1)
(Disclaimer: there may be a better way got get it with library Design by Frank Harrell, but let's assume we have to do it the hard way) Looks like it is a bit hidden. f1 is of class(print.survfit), as str(f1) tells us. So let's try getAnyhwere(print.survfit). In the lower part you find line like the following: x1 <- pfun(nsubjects, stime, surv, x$n.risk, x$n.event, x$lower, x$upper) if (is.matrix(x1)) { if (is.null(x$lower)) dimnames(x1) <- list(NULL, plab) else dimnames(x1) <- list(NULL, c(plab, plab2)) } else { if (is.null(x$lower)) names(x1) <- plab else names(x1) <- c(plab, plab2) } if (show.rmean) print(x1) Make a copy of that function under a new name, and return x1. Dieter ______________________________________________ R-help@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.