summary.coxph. print.summary.coxph is needed so that printing of the objects
that are thus created has the same effect as calling the existing summary.coxph
It needs a modified help file to go with it. I have been meaning to try to get this,
or something like it, incorporated into the survival package.
"sum.coxph" <-
function (object, table = TRUE, coef = TRUE, conf.int = 0.95,
scale = 1, ...)
{
cox <- object
coxsum <- list(call=cox$call, fail=cox$fail, na.action=cox$na.action,
n=cox$n, icc=cox$icc, naive.var=cox$naive.var)
if (!is.null(cox$fail)) {
return()
}
omit <- cox$na.action
if (is.null(cox$coef)) {
return()
}
beta <- cox$coef
nabeta <- !(is.na(beta))
beta2 <- beta[nabeta]
if (is.null(beta) | is.null(cox$var)){
coxsum$invalid <- TRUE
return()
}
else coxsum$invalid<-FALSE
se <- sqrt(diag(cox$var))
if (!is.null(cox$naive.var))
nse <- sqrt(diag(cox$naive.var))
if (coef) {
if (is.null(cox$naive.var)) {
tmp <- cbind(beta, exp(beta), se, beta/se, 1 -
pchisq((beta/se)^2, 1))
dimnames(tmp) <- list(names(beta), c("coef", "exp(coef)",
"se(coef)", "z", "p"))
}
else {
tmp <- cbind(beta, exp(beta), nse, se, beta/se, 1 -
pchisq((beta/se)^2, 1))
dimnames(tmp) <- list(names(beta), c("coef", "exp(coef)",
"se(coef)", "robust se", "z", "p"))
}
} else tmp <- NULL
if (conf.int) {
z <- qnorm((1 + conf.int)/2, 0, 1)
beta <- beta * scale
se <- se * scale
tmpc <- cbind(exp(beta), exp(-beta), exp(beta - z * se),
exp(beta + z * se))
dimnames(tmpc) <- list(names(beta), c("exp(coef)", "exp(-coef)",
paste("lower .", round(100 * conf.int, 2), sep = ""),
paste("upper .", round(100 * conf.int, 2), sep = "")))
} else tmpc <- NULL
logtest <- -2 * (cox$loglik[1] - cox$loglik[2])
sctest <- cox$score
df <- length(beta2)
Rsquare <- 1 - exp(-logtest/cox$n)
maxRsquare <- 1 - exp(2 * cox$loglik[1]/cox$n)
plogtest <- 1 - pchisq(logtest, df)
pwald.test <- 1 - pchisq(cox$wald.test, df)
pscore <- 1 - pchisq(sctest, df)
if(!is.null(cox$rscore)){
rscore <- cox$rscore
prscore <- 1 - pchisq(cox$rscore, df)
robust <- list(val=rscore, pval=prscore)
}
else robust <- NULL
wald.test <- cox$wald.test
coxsum <- c(list(call=cox$call, coefficients=tmp, intervals=tmpc, df=df,
rsquare=c(Rsquare, maxRsquare),
tests=c(logtest=logtest, wald.test=wald.test, score=sctest),
p.tests = c(logtest=plogtest, wald.test=pwald.test,
score=pscore), robust=robust), coxsum)
class(coxsum) <- "summary.coxph"
coxsum
}
"print.summary.coxph" <-
function (object, table = TRUE, coef = TRUE, conf.int = 0.95,
scale = 1, digits = max(options()$digits - 4, 3), ...)
{
coxsum <- object
if (!is.null(cl <- coxsum$call)) {
cat("Call:\n")
dput(cl)
cat("\n")
}
if (!is.null(coxsum$fail)) {
cat(" Coxreg failed.", coxsum$fail, "\n")
return()
}
savedig <- options(digits = digits)
on.exit(options(savedig))
omit <- coxsum$na.action
if (length(omit))
cat(" n=", coxsum$n, " (", naprint(omit), ")\n", sep = "")
else cat(" n=", coxsum$n, "\n")
if (length(coxsum$icc))
cat(" robust variance based on", coxsum$icc[1], "groups, intra-class correlation =",
format(coxsum$icc[2:3]), "\n")
if (is.null(coxsum$coef)) {
cat(" Null model\n")
return()
}
if (coxsum$invalid)stop("input is invalid")
if (!is.null(coxsum$coef)) {
cat("\n")
coxsum$coef[,"p"] <- signif(coxsum$coef[,"p"],digits-1)
prmatrix(coxsum$coef)
}
if (!is.null(coxsum$intervals)) {
cat("\n")
prmatrix(coxsum$intervals)
}
logtest <- coxsum$tests["logtest"]
sctest <- coxsum$tests["score"]
df <- coxsum$df
cat("\n")
cat("Rsquare=", format(round(coxsum$rsquare[1], 3)),
" (max possible=", format(round(coxsum$rsquare[2],3)), ")\n")
cat("Likelihood ratio test= ", format(round(coxsum$tests["logtest"], 2)),
" on ", df, " df,", " p=", format(coxsum$p.tests["logtest"]),
"\n", sep = "")
cat("Wald test = ", format(coxsum$tests["wald.test.x"]),
" on ", df, " df,", " p=", format(
coxsum$p.tests["wald.test.x"]), "\n", sep = "")
cat("Score (logrank) test = ", format(round(coxsum$tests["score"], 2)),
" on ", df, " df,", " p=", format(coxsum$p.tests["score"]),
sep = "")
if (is.null(coxsum$robust))
cat("\n\n")
else cat(", Robust = ", format(round(coxsum$robust$val, 2)), " p=",
format(coxsum$robust$pval), "\n\n", sep = "")
if (!is.null(coxsum$tests))
cat(" (Note: the likelihood ratio and score tests",
"assume independence of\n observations within a cluster,",
"the Wald and robust score tests do not).\n")
invisible()
}
John Maindonald email: [EMAIL PROTECTED] phone : +61 2 (6125)3473 fax : +61 2(6125)5549 Centre for Bioinformation Science, Room 1194, John Dedman Mathematical Sciences Building (Building 27) Australian National University, Canberra ACT 0200.
______________________________________________ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
