[R] Package Hmisc, functions summary.formula() and latex(), options pdig, pctdig, eps and prmsd

Frank E Harrell Jr f.harrell at vanderbilt.edu
Fri Jul 25 14:25:05 CEST 2008


David Hajage wrote:
> Hello R users,
> 
> I have several problems with the functions summary.formula and latex in
> package Hmisc. Here an example :
> 
> ####
> library(Hmisc)
> sex <- factor(sample(c("m","f", "?"), 50, rep=TRUE))
> age <- rnorm(50, 50, 5)
> treatment <- factor(sample(c("Drug","Placebo"), 50, rep=TRUE))
> symp <- c('Headache','Stomach Ache','Hangnail',
>           'Muscle Ache','Depressed')
> symptom1 <- sample(symp, 50,TRUE)
> symptom2 <- sample(symp, 50,TRUE)
> symptom3 <- sample(symp, 50,TRUE)
> Symptoms <- mChoice(symptom1, symptom2, symptom3, label='Primary Symptoms')
> 
> f <- summary(treatment ~ age + sex + Symptoms, method="reverse", test=TRUE)
> print(f, digits = 5, pdig = 2, pctdig = 3, eps = 0.5, prmsd = F)
> print(f, digits = 5, pdig = 2, pctdig = 3, eps = 0.5, prmsd = T)
> latex(f, long = T, pdig = 2, pctdig = 3, eps = 0.5, prmsd = F, file = "")
> ###
> 
> Here the problems :
>  - The first print(f, ...) doesn't replace all p-value <0.5 by "P<0.5" :
> 
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |                            |Drug
> |Placebo                  |  Test                       |
> |                            |(N=31)
> |(N=19)                   |Statistic                    |
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |age                         |     45.926/48.750/54.019|
> 47.344/50.728/53.696|    F=0.9 d.f.=1,48 P<0.5    |
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |sex : ?                     |         25.806% ( 8)    |         26.316% (
> 5)    | Chi-square=0.3 d.f.=2 P=0.86|
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |    f                       |         38.710% (12)    |         31.579% (
> 6)    |                             |
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |    m                       |         35.484% (11)    |         42.105% (
> 8)    |                             |
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |Primary Symptoms : Depressed|         41.935% (13)    |         63.158%
> (12)    | Chi-square=2.12 d.f.=1 P<0.5|
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |    Hangnail                |         48.387% (15)    |         42.105% (
> 8)    |Chi-square=0.19 d.f.=1 P=0.67|
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |    Stomach Ache            |         45.161% (14)    |         68.421%
> (13)    | Chi-square=2.57 d.f.=1 P<0.5|
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |    Muscle Ache             |         54.839% (17)    |         26.316% (
> 5)    | Chi-square=3.89 d.f.=1 P<0.5|
> +----------------------------+-------------------------+-------------------------+-----------------------------+
> |    Headache                |         51.613% (16)    |         36.842% (
> 7)    | Chi-square=1.03 d.f.=1 P<0.5|
> +----------------------------+-------------------------+-------------------------+-----------------------------+

I did not see an instance of P < 0.5 that was not replaced by "P<0.5". 
This is an unusual cutoff, and you are printing many more digits of 
precision than offered by the data.

> 
>  - The second print(f, ..., prmsd = T) has the same problem for p-value.
> There is also 4 decimals in the first line instead of 3 :

Try modifying the digits argument.

> 
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |                            |Drug
> |Placebo                                   |
> Test                              |
> |                            |(N=31)
> |(N=19)
> |Statistic                           |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |age                         |45.9263/48.7501/54.0194  49.1817+/-
> 5.2751|47.3436/50.7276/53.6963  51.1313+/- 5.5868|           F=0.9 d.f.=1,48
> P<0.5    |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |sex : ?                     |               25.806% ( 8)
> |               26.316% ( 5)               |       Chi-square=0.3 d.f.=2
> P=0.86 |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |    f                       |               38.710% (12)
> |               31.579% ( 6)
> |                                    |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |    m                       |               35.484% (11)
> |               42.105% ( 8)
> |                                    |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |Primary Symptoms : Depressed|               41.935% (13)
> |               63.158% (12)               |       Chi-square=2.12 d.f.=1
> P<0.5 |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |    Hangnail                |               48.387% (15)
> |               42.105% ( 8)               |       Chi-square=0.19 d.f.=1
> P=0.67|
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |    Stomach Ache            |               45.161% (14)
> |               68.421% (13)               |       Chi-square=2.57 d.f.=1
> P<0.5 |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |    Muscle Ache             |               54.839% (17)
> |               26.316% ( 5)               |       Chi-square=3.89 d.f.=1
> P<0.5 |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> |    Headache                |               51.613% (16)
> |               36.842% ( 7)               |       Chi-square=1.03 d.f.=1
> P<0.5 |
> +----------------------------+------------------------------------------+------------------------------------------+------------------------------------+
> 
>  - In the latex(f, ...), only the first-line p-value is replaced by "< 0.5",
> and the quantiles for age have too much decimal :

This is a bug.  An override source is below, to use until the next release.

> 
> % latex.default(cstats, title = title, caption = caption, rowlabel =
> rowlabel,      col.just = col.just, numeric.dollar = FALSE, insert.bottom =
> legend,      rowname = lab, dcolumn = dcolumn, extracolheads =
> extracolheads,      extracolsize = Nsize, ...)
> %
> \begin{table}[!tbp]
>  \caption{Descriptive Statistics by treatment\label{f}}
>  \begin{center}
>  \begin{tabular}{lccc}\hline\hline
> \multicolumn{1}{l}{}&
> \multicolumn{1}{c}{Drug}&
> \multicolumn{1}{c}{Placebo}&
> \multicolumn{1}{c}{Test Statistic}
> \\   &\multicolumn{1}{c}{{\scriptsize
> $N=31$}}&\multicolumn{1}{c}{{\scriptsize $N=19$}}&\\ \hline
> age&{\scriptsize 45.92627049~}{48.75014585 }{\scriptsize 54.01942542}
> &{\scriptsize 47.34358486~}{50.72757132 }{\scriptsize 53.69634575} &$
> F_{1,48}=0.9 ,~ P<0.5 ^{1} $\\
> sex&&&$ \chi^{2}_{2}=0.3 ,~ P=0.859 ^{2} $\\
> ~~~~?&25.806\%~{\scriptsize~(~8)}&26.316\%~{\scriptsize~(~5)}&\\
> ~~~~f&38.710\%~{\scriptsize~(12)}&31.579\%~{\scriptsize~(~6)}&\\
> ~~~~m&35.484\%~{\scriptsize~(11)}&42.105\%~{\scriptsize~(~8)}&\\
> Primary~Symptoms&&&\\
> ~~~~Depressed&41.935\%~{\scriptsize~(13)}&63.158\%~{\scriptsize~(12)}&$
> \chi^{2}_{1}=2.12 ,~ P=0.145 ^{2} $\\
> ~~~~Hangnail&48.387\%~{\scriptsize~(15)}&42.105\%~{\scriptsize~(~8)}&$
> \chi^{2}_{1}=0.19 ,~ P=0.665 ^{2} $\\
> ~~~~Stomach~Ache&45.161\%~{\scriptsize~(14)}&68.421\%~{\scriptsize~(13)}&$
> \chi^{2}_{1}=2.57 ,~ P=0.109 ^{2} $\\
> ~~~~Muscle~Ache&54.839\%~{\scriptsize~(17)}&26.316\%~{\scriptsize~(~5)}&$
> \chi^{2}_{1}=3.89 ,~ P=0.049 ^{2} $\\
> ~~~~Headache&51.613\%~{\scriptsize~(16)}&36.842\%~{\scriptsize~(~7)}&$
> \chi^{2}_{1}=1.03 ,~ P=0.309 ^{2} $\\
> \hline
> \end{tabular}
> 
> \end{center}
> 
> \noindent {\scriptsize $a$\ }{$b$\ }{\scriptsize $c$\ } represent the lower
> quartile $a$, the median $b$, and the upper quartile $c$\ for continuous
> variables.\\
> Numbers after percents are frequencies.\\
>  \indent Tests used: $^{1}$Wilcoxon test; $^{2}$Pearson test
> \end{table}
> 
> I would like to understand :
>  - why does option 'pctdig' affect both percentages and quantiles in the
> 'print' command (pctdig : 'number of digits to the right of the decimal
> place for *printing percentages*. The default is zero, so percents will be
> rounded to the nearest percent');
>  - how can I change the number of decimal for continuous variables in the
> 'latex' command ?

I hope to address these 2 questions soon.

Frank

>  - why doesn't option 'eps' affect all p-value ?
> 
> I've just discover these wonderfull functions. Thank you to Frank Harrell
> for this package.
> 
> DH
> 

latex.summary.formula.reverse <-
   function(object, title=first.word(deparse(substitute(object))),
            digits, prn = any(n!=N), pctdig=0,
            npct=c('numerator','both','denominator','none'),
            npct.size='scriptsize', Nsize='scriptsize',
            exclude1=TRUE,  vnames=c("labels","names"), prUnits=TRUE,
            middle.bold=FALSE, outer.size="scriptsize",
            caption, rowlabel="",
            insert.bottom=TRUE, dcolumn=FALSE, formatArgs=NULL, round=NULL,
            prtest=c('P','stat','df','name'), prmsd=FALSE, msdsize=NULL,
            long=dotchart, pdig=3, eps=.001, auxCol=NULL, 
dotchart=FALSE, ...)
{
   x      <- object
   npct   <- match.arg(npct)
   vnames <- match.arg(vnames)
   if(is.logical(prtest) && !prtest)
     prtest <- 'none'

   stats  <- x$stats
   nv     <- length(stats)
   cstats <- lab <- character(0)
   nn     <- integer(0)
   type   <- x$type
   n      <- x$n
   N      <- x$N
   nams   <- names(stats)
   labels <- x$labels
   Units  <- x$units
   nw     <- if(lg <- length(x$group.freq)) lg
             else 1  #23Nov98
   gnames <- names(x$group.freq)
   test   <- x$testresults
   if(!length(test))
     prtest <- 'none'

   gt1.test <-
     if(all(prtest=='none'))
       FALSE
     else
       length(unique(sapply(test,function(a)a$testname))) > 1

   if(!missing(digits)) {   #.Options$digits <- digits 6Aug00
     oldopt <- options(digits=digits)
     on.exit(options(oldopt))
   }

   if(missing(caption))
     caption <- paste("Descriptive Statistics",
                      if(length(x$group.label))
                        paste(" by",x$group.label)
                      else
                        paste("  $(N=",x$N,")$",sep=""), sep="")

   bld <- if(middle.bold) '\\bf '
          else ''

   cstats <- NULL
   testUsed <- auxc <- character(0)

   for(i in 1:nv) {
     if(length(auxCol))
       auxc <- c(auxc, auxCol[[1]][i])

     nn <- c(nn, n[i])   ## 12aug02
     nam <- if(vnames=="names") nams[i]
            else labels[i]

     if(prUnits && nchar(Units[i]) > 0)
       nam <- paste(nam, '~\\hfill\\tiny{',translate(Units[i],'*',' 
'),'}',sep='')

     tr  <- if(length(test) && all(prtest!='none')) test[[nams[i]]]
            else NULL

     if(length(test) && all(prtest!='none'))
       testUsed <- unique(c(testUsed, tr$testname))

     if(type[i]==1 || type[i]==3) {
       cs <- formatCats(stats[[i]], nam, tr, type[i],
                        if(length(x$group.freq)) x$group.freq else x$n[i],
                        npct, pctdig, exclude1, long, prtest,
                        latex=TRUE, testUsed=testUsed,
                        npct.size=npct.size,
                        pdig=pdig, eps=eps,
                        footnoteTest=gt1.test, dotchart=dotchart)
       nn <- c(nn, rep(NA, nrow(cs)-1))
     } else cs <- formatCons(stats[[i]], nam, tr, x$group.freq, prmsd,
                             prtest=prtest, formatArgs=formatArgs, 
round=round,
                             latex=TRUE, testUsed=testUsed,
                             middle.bold=middle.bold,
                             outer.size=outer.size, msdsize=msdsize,
                             pdig=pdig, eps=eps, footnoteTest=gt1.test)

     cstats <- rbind(cstats, cs)
     if(length(auxc) && nrow(cstats) > 1)
       auxc <- c(auxc, rep(NA, nrow(cs)-1))
   }

   lab <- dimnames(cstats)[[1]]
   gl <- names(x$group.freq)
   ##gl <- if(length(gl)) paste(gl, " $(N=",x$group.freq,")$",sep="") 
else " "
   ## Thanks: Eran Bellin <ebellin at montefiore.org>   3Aug01
   if(!length(gl))
     gl <- " "

   lab <- sedit(lab,c(" ","&"),c("~","\\&"))  #was format(lab) 21Jan99
   lab <- latexTranslate(lab, greek=.R.)
   gl  <- latexTranslate(gl, greek=.R.)
   ## if(any(gl != " ")) gl <- paste(gl, " 
$(N=",x$group.freq,")$",sep="") # 3Aug01
   ## Added any( ) 26Mar02  21jan03
   extracolheads <-
     if(any(gl != " "))
       c(if(prn)'', paste('$N=',x$group.freq,'$',sep=''))
     else NULL # 21jan03

   if(length(test) && !all(prtest=='none')) {
     gl <- c(gl,
             if(length(prtest)==1 && prtest!='stat')
               if(prtest=='P') 'P-value'
               else prtest
             else 'Test Statistic')

     if(length(extracolheads)) extracolheads <- c(extracolheads,'') # 
21jan03
   }

   dimnames(cstats) <- list(NULL,gl)
   ## was dimnames(cstats) <- list(lab, gl) 12aug02
   cstats <- data.frame(cstats, check.names=FALSE)

   ## Added row.names=lab below 10jul02 - S+ was dropping dimnames[[1]]
   ##attr(cstats,'row.names') <- lab  12aug02
   col.just <- rep("c",length(gl))
   if(dcolumn && all(prtest!='none') &&
      gl[length(gl)] %in% c('P-value','Test Statistic'))
     col.just[length(col.just)] <- '.'

   if(prn) {
     cstats <- data.frame(N=nn, cstats, check.names=FALSE)
     col.just <- c("r",col.just)
   }

   if(!insert.bottom)
     legend <- NULL
   else {
     legend <- paste(if(any(type==2)) {
                       paste("\\noindent {\\",outer.size," $a$\\ 
}{",bld,"$b$\\ }{\\",
                             outer.size," $c$\\ } represent the lower 
quartile $a$, the median $b$, and the upper quartile $c$\\ for 
continuous variables.",
                             if(prmsd) '~~$x\\pm s$ represents 
$\\bar{X}\\pm 1$ SD.'
                             else '',
                             '\\\\', sep="")
                     } else NULL,
                     if(prn) '$N$\\ is the number of non--missing 
values.\\\\',
                     if(any(type==1) && npct=='numerator')
                       'Numbers after percents are frequencies.\\\\',
                     sep="\n")
     legend <- NULL
     if(any(type==2)) {
       legend <- paste("\\noindent {\\", outer.size, " $a$\\ }{", bld,
                       "$b$\\ }{\\", outer.size,
                       " $c$\\ } represent the lower quartile $a$, the 
median $b$, and the upper quartile $c$\\ for continuous variables.",
                       if(prmsd) '~~$x\\pm s$ represents $\\bar{X}\\pm 
1$ SD.'
                       else '',
                       '\\\\\n', sep="")
     }

     if(prn) {
       legend <- paste(legend,
                       '$N$\\ is the number of non--missing values.\\\\\n',
                       sep='')
     }

     if(any(type==1) && npct=='numerator') {
       legend <- paste(legend,
                       'Numbers after percents are frequencies.\\\\\n',
                       sep='')
     }

     if(length(testUsed))
       legend <-paste(legend,
                      if(length(testUsed)==1)'\\noindent Test used:'
                      else '\\indent Tests used:',
                      if(length(testUsed)==1) paste(testUsed,'test')
                      else
                        paste(paste('$^{',1:length(testUsed),'}$',testUsed,
                                    ' test',sep=''),collapse='; '))

     ## added rowname=lab 12aug02  added '\n\n' 4mar03 for ctable=T
   }

   if(length(auxc)) {
     if(length(auxc) != nrow(cstats))
       stop(paste('length of auxCol (',length(auxCol[[1]]),
                  ') is not equal to number or variables in table (',
                  nv,').', sep=''))
     auxcc <- format(auxc)
     auxcc[is.na(auxc)] <- ''
     cstats <- cbind(auxcc, cstats)
     nax <- names(auxCol)
     heads <- get2rowHeads(nax)
     names(cstats)[1] <- heads[[1]]
     if(length(col.just)) col.just <- c('r', col.just)
     if(length(extracolheads)) extracolheads <- c(heads[2], extracolheads)
   }
   resp <- latex.default(cstats, title=title, caption=caption, 
rowlabel=rowlabel,
                         col.just=col.just, numeric.dollar=FALSE,
                         insert.bottom=legend,  rowname=lab, 
dcolumn=dcolumn,
                         extracolheads=extracolheads, extracolsize=Nsize,
                         ...)

   if(dotchart)
     resp$style <- unique(c(resp$style, 'calc', 'epic', 'color'))

   resp
}


-- 
Frank E Harrell Jr   Professor and Chair           School of Medicine
                      Department of Biostatistics   Vanderbilt University



More information about the R-help mailing list