R/export.R

In yapomif: Yet another package of misc functions

## Funzioni utili per l'esportazione tabellare di output
## da oggetti di varia natura, principalmente da modelli
export <- function(x, ...) {
    UseMethod("export")
}


export.coxph <- function(
    x, ## coxph object
    file = NULL,
    coef=FALSE,            ##print beta?
    se.coef=FALSE,        ##print se(beta)?
    z=FALSE,              ##print z
    p=TRUE,               ##print p
    exp.coef=TRUE,        ##print exp(beta)
    exp.coef.lci=TRUE,    ##print confint(beta) low
    exp.coef.uci=TRUE,     ##print confint(beta) low
    print.formula=TRUE,
    print.n=TRUE,
    print.tests=TRUE,
    print.ph.test=TRUE,
    ...

)  {
    if (is.null(file)) {
        stop("file must be specified")
    }
    
    xs <- summary(x)
    
    coeffs <- xs$coefficients
    ## qui ci vuole qualcosa che stabilisca se l'analisi è univariata
    ## o multivariata (poi serve anche alla 
    
    nvars <- nrow(coeffs)
    
    ci <- xs$conf.int[,3:4] ## lower e upper ci solamente
    ## ci: se la stima è univariata restituisce un vettore, 
    ## che va trattato opportunamente
    
    if(nvars==1) {
        nomi <- names(ci)
        dim(ci) <- c(1,2)
        colnames(ci) <- nomi
    }
    results <- cbind(coeffs, ci)

    ## Variabili da tenere (per ora ignora altre tipo robust se e
    ## altre possibili colonne
    results <- results[,colnames(results) %in%
                       c("coef",
                         "exp(coef)",
                         "se(coef)",         
                         "z",    
                         "Pr(>|z|)", 
                         "lower .95", 
                         "upper .95") ]
    
    ## ordine attuale
    ## --------------
    ## 1) coef 
    ## 2) exp(coef)    
    ## 3) se(coef)         
    ## 4) z    
    ## 5) Pr(>|z|) 
    ## 6) lower .95 
    ## 7) upper .95
    
    ## ordine da impostare
    ## c(1,3,4,2,6:7,5)

    ## 1) coef 
    ## 3) se(coef)         
    ## 4) z    
    ## 2) exp(coef)    
    ## 6) lower .95 
    ## 7) upper .95
    ## 5) Pr(>|z|) 
    
    my.order <- c(1,3,4,2,6:7,5)

    if(nvars > 1) {
        results <- results[, my.order, drop=F]
    } else if (nvars==1) {
        results <- results[my.order]
    }

    stat.select <- c(coef,
                     se.coef,
                     z,
                     exp.coef,
                     exp.coef.lci,
                     exp.coef.uci,
                     p)
    
    ## Fase finale di consolidamento output
    var.name <- row.names(coeffs)
    out <- data.frame("var"=var.name)
    if (nvars>1) {
        out <- cbind(out, results[, stat.select] )
    } else if (nvars==1) {
        out <- cbind(out, t(results[stat.select]))
    }
    row.names(out) <- NULL

    ## Pulisco eventuali precedenti file
    cat("\n", file = file)
    
    ## Intestazione: chiamata, numero e numero eventi
    if (print.formula) {
        cat(Reduce(paste, deparse(x$formula)),
            sep=" ", file=file, append=TRUE)
        cat("\n", file=file, append=TRUE)
    }
    if (print.n) {
        cat("\n\t","n\t", xs$n, "\t", "events\t", xs$nevent,"\n\n",
            file=file, append=TRUE, sep="") 
    }
    
    ## Statistiche principali
    write.delim(out, file=file, row.names = FALSE, append = TRUE)
    cat("\n", file=file, append=TRUE)

    ## Test
    if (print.tests) {
        ## header
        cat("STAT","test", "df","pvalue", sep="\t", file=file, append=TRUE)
        cat("\n", file=file, append=TRUE)
        ## lrt
        cat("LR test", xs$logtest, sep="\t",file=file, append=TRUE)
        cat("\n", file=file, append=TRUE)
        ## wald
        cat("Wald test", xs$waldtest, sep="\t",file=file, append=TRUE)
        cat("\n", file=file, append=TRUE)
        ## Logrank
        cat("Score (logrank) test", xs$sctest, sep="\t",file=file, append=TRUE)
        cat("\n", file=file, append=TRUE)
        
    }

    ## phtest
    if (print.ph.test) {
        ph <- as.data.frame(cox.zph(x)$table)
        ph <- cbind("var"=row.names(ph), ph)
        row.names(ph) <- NULL
        cat("\nProportional Hazard Tests\n",
            file = file, append = TRUE)
        write.delim(ph, file=file, append=TRUE)
    }
    
    
    ## example(summary.coxph)
    ## sfit <- summary(fit)
    ## export(sfit)
}



## export.cox.zph <- function(
##     x, ## cox.zph(coxph) object
##     file=NULL,
##     ...
##     )  {

##     ## Scrivi una intestazione: se è in appending scrivila in append
##     ## se no scrivila sopra
    
##     dots <- as.list(substitute(list(...)))[-1]
    
##     if ("append" %in% names(dots)) {
##         if (dots$append == TRUE) {
##             append.to.file <- TRUE
##         } else {
##             append.to.file <- FALSE
##         }
##     } else { ## append non specificato a linea di comando
##         append.to.file <- FALSE
##     }
    
##     ## Questo di sotto deve comunque essere in appending rispetto
##     ## all'intestazione 
##     cat("\n\nProportional Hazard Tests\n",
##         file = file, append = append.to.file)
##     write.delim(x, file=file, append=TRUE)
    
##     ## example(summary.coxph)
##     ## sfit <- summary(fit)
##     ## export(sfit)
## }