The R Project SVN R-packages

#SCCS 03/25/97 @(#)summary.coxph.s  4.8
summarycoxph <-
 function(object, table = TRUE, coef = TRUE, conf.int = 0.95, scale = 1,
          digits=getOption("digits"),...)
    {
        cox<-object
    if (!is.null(cl<- cox$call)) {
    cat("Call:\n")
    dput(cl)
    cat("\n")
    }
    if (!is.null(cox$fail)) {
    cat(" Coxreg failed.", cox$fail, "\n")
    return()
    }
    savedig <- options(digits = digits)
    on.exit(options(savedig))

    omit <- cox$na.action
    if (length(omit))
    cat("  n=", cox$n, " (", naprint(omit), ")\n", sep="")
    else cat("  n=", cox$n, "\n")
    if (length(cox$icc))
    cat("  robust variance based on", cox$icc[1],
        "groups, intra-class correlation =", format(cox$icc[2:3]), "\n")
    if (is.null(cox$coefficients)) {   # Null model
    cat ("   Null model\n")
    return()
    }

    beta <- cox$coefficients
    nabeta <- !(is.na(beta))          #non-missing coefs
    beta2 <- beta[nabeta]
    if(is.null(beta) | is.null(cox$var))
        stop("Input is not valid")
    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,
           signif(1 - pchisq((beta/ se)^2, 1), digits -1))
        dimnames(tmp) <- list(names(beta), c("coef", "exp(coef)",
        "se(coef)", "z", "p"))
        }
    else {
        tmp <- cbind(beta, exp(beta), nse, se, beta/se,
           signif(1 - pchisq((beta/ se)^2, 1), digits -1))
        dimnames(tmp) <- list(names(beta), c("coef", "exp(coef)",
        "se(coef)", "robust se", "z", "p"))
        }
        cat("\n")
        prmatrix(tmp)
        }
    if(conf.int) {
        z <- qnorm((1 + conf.int)/2, 0, 1)
        beta <- beta * scale
        se <- se * scale
        tmp <- cbind(exp(beta), exp(-beta), exp(beta - z * se),
            exp(beta + z * se))
        dimnames(tmp) <- list(names(beta), c("exp(coef)", "exp(-coef)",
            paste("lower .", round(100 * conf.int, 2), sep = ""),
            paste("upper .", round(100 * conf.int, 2), sep = "")))
        cat("\n")
        prmatrix(tmp)
        }
    logtest <- -2 * (cox$loglik[1] - cox$loglik[2])
    sctest <- cox$score
    df <- length(beta2)
    cat("\n")
    cat("Rsquare=", format(round(1-exp(-logtest/cox$n),3)),
    "  (max possible=", format(round(1-exp(2*cox$loglik[1]/cox$n),3)),
    ")\n" )
    cat("Likelihood ratio test= ", format(round(logtest, 2)), "  on ",
    df, " df,", "   p=", format(1 - pchisq(logtest, df)),
    "\n", sep = "")
    cat("Wald test            = ", format(round(cox$wald.test, 2)), "  on ",
    df, " df,", "   p=", format(1 - pchisq(cox$wald.test, df)),
    "\n", sep = "")
    cat("Score (logrank) test = ", format(round(sctest, 2)), "  on ", df,
        " df,", "   p=", format(1 - pchisq(sctest, df)), sep ="") 
    if (is.null(cox$rscore)) cat("\n\n")
    else cat(",   Robust = ", format(round(cox$rscore, 2)), 
       "  p=", format(1 - pchisq(cox$rscore, df)), "\n\n", sep="")   

    if (!is.null(cox$naive.var))
    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()
    }

summary.coxph <-
 function(object,  coef = TRUE, conf.int = 0.95, scale = 1,
          digits = max(options()$digits - 4, 3),...)
{
     cox<-object
     beta <- cox$coefficients
     nabeta <- !(is.na(beta))          #non-missing coefs
     beta2 <- beta[nabeta]
     if(is.null(beta) | is.null(cox$var))
         stop("Input is not valid")
     se <- sqrt(diag(cox$var))
     if (!is.null(cox$naive.var)) nse <- sqrt(diag(cox$naive.var))

     rval<-list(call=cox$call,fail=cox$fail, na.action=cox$na.action,
                n=cox$n,icc=cox$icc)
     if(coef) {
         if (is.null(cox$naive.var)) {
        tmp <- cbind(beta, exp(beta), se, beta/se,
           signif(1 - pchisq((beta/ se)^2, 1), digits -1))
        dimnames(tmp) <- list(names(beta), c("coef", "exp(coef)",
        "se(coef)", "z", "p"))
        }
    else {
        tmp <- cbind(beta, exp(beta), nse, se, beta/se,
           signif(1 - pchisq((beta/ se)^2, 1), digits -1))
        dimnames(tmp) <- list(names(beta), c("coef", "exp(coef)",
        "se(coef)", "robust se", "z", "p"))
        }
        rval$coef<-tmp
        }
    if(conf.int) {
        z <- qnorm((1 + conf.int)/2, 0, 1)
        beta <- beta * scale
        se <- se * scale
        tmp <- cbind(exp(beta), exp(-beta), exp(beta - z * se),
            exp(beta + z * se))
        dimnames(tmp) <- list(names(beta), c("exp(coef)", "exp(-coef)",
            paste("lower .", round(100 * conf.int, 2), sep = ""),
            paste("upper .", round(100 * conf.int, 2), sep = "")))
        rval$conf.int<-tmp
    }
    df <- length(beta2)
     logtest <- -2 * (cox$loglik[1] - cox$loglik[2])
    rval$logtest <- c(test=logtest,
                      df=df,
                      pvalue=1 - pchisq(logtest, df))
    rval$sctest <- c(test=cox$score,
                     df=df,
                     pvalue=1 - pchisq(cox$score, df))
     rval$rsq<-c(rsq=1-exp(-logtest/cox$n),
                 maxrsq=1-exp(2*cox$loglik[1]/cox$n))
     rval$waldtest<-c(test=as.vector(round(cox$wald.test, 2)),
                      df=df,
                      pvalue=1 - pchisq(as.vector(cox$wald.test), df))
     if (!is.null(cox$rscore))
         rval$robscore<-c(test=cox$rscore,
                          df=df,
                          pvalue=1 - pchisq(cox$rscore, df))
     rval$used.robust<-!is.null(cox$naive.var)
     class(rval)<-"summary.coxph"
     rval
 }


print.summary.coxph <-
 function(x, digits = max(options()$digits - 4, 3), ...)
    {
    if (!is.null(x$call)) {
    cat("Call:\n")
    dput(x$call)
    cat("\n")
    }
    if (!is.null(x$fail)) {
    cat(" Coxreg failed.", x$fail, "\n")
    return()
    }
    savedig <- options(digits = digits)
    on.exit(options(savedig))

    omit <- x$na.action
    if (length(omit))
    cat("  n=", x$n, " (", naprint(omit), ")\n", sep="")
    else cat("  n=", x$n, "\n")
    if (length(x$icc))
    cat("  robust variance based on", x$icc[1],
        "groups\n")
    if (nrow(x$coef)==0) {   # Null model
    cat ("   Null model\n")
    return()
    }


    if(!is.null(x$coef)) {
        prmatrix(x$coef)
        }
    if(!is.null(x$conf.int)) {
        cat("\n")
        prmatrix(x$conf.int)
        }
    cat("\n")
    cat("Rsquare=", format(round(x$rsq["rsq"],3)),
    "  (max possible=", format(round(x$rsq["maxrsq"],3)),
        ")\n" )
    cat("Likelihood ratio test= ", format(round(x$logtest["test"], 2)), "  on ",
    x$logtest["df"], " df,", "   p=", format(x$logtest["pvalue"]),
        "\n", sep = "")
    cat("Wald test            = ", format(round(x$waldtest["test"], 2)), "  on ",
    x$waldtest["df"], " df,", "   p=", format(x$waldtest["pvalue"]),
    "\n", sep = "")
    cat("Score (logrank) test = ", format(round(x$sctest["test"], 2)), "  on ",
        x$sctest["df"]," df,", "   p=", format(x$sctest["pvalue"]), sep ="")
    if (is.null(x$robscore))
        cat("\n\n")
    else cat(",   Robust = ", format(round(x$robscore["test"], 2)), 
             "  p=", format(x$robscore["pvalue"]), "\n\n", sep="")   

    if (x$used.robust)
    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()
}