View source: R/summary.formula.s
summary.formula | R Documentation |
summary.formula
summarizes the variables listed in an S formula,
computing descriptive statistics (including ones in a
user-specified function). The summary statistics may be passed to
print
methods, plot
methods for making annotated dot charts, and
latex
methods for typesetting tables using LaTeX.
summary.formula
has three methods for computing descriptive
statistics on univariate or multivariate responses, subsetted by
categories of other variables. The method of summarization is
specified in the parameter method
(see details below). For the
response
and cross
methods, the statistics used to
summarize the data
may be specified in a very flexible way (e.g., the geometric mean,
33rd percentile, Kaplan-Meier 2-year survival estimate, mixtures of
several statistics). The default summary statistic for these methods
is the mean (the proportion of positive responses for a binary
response variable). The cross
method is useful for creating data
frames which contain summary statistics that are passed to trellis
as raw data (to make multi-panel dot charts, for example). The
print
methods use the print.char.matrix
function to print boxed
tables.
The right hand side of formula
may contain mChoice
(“multiple choice”) variables. When test=TRUE
each choice is
tested separately as a binary categorical response.
The plot
method for method="reverse"
creates a temporary
function Key
in frame 0 as is done by the xYplot
and
Ecdf.formula
functions. After plot
runs, you can type
Key()
to put a legend in a default location, or
e.g. Key(locator(1))
to draw a legend where you click the left
mouse button. This key is for categorical variables, so to have the
opportunity to put the key on the graph you will probably want to use
the command plot(object, which="categorical")
. A second function
Key2
is created if continuous variables are being plotted. It is
used the same as Key
. If the which
argument is not
specified to plot
, two pages of plots will be produced. If you
don't define par(mfrow=)
yourself,
plot.summary.formula.reverse
will try to lay out a multi-panel
graph to best fit all the individual dot charts for continuous
variables.
There is a subscripting method for objects created with
method="response"
.
This can be used to print or plot selected variables or summary statistics
where there would otherwise be too many on one page.
cumcategory
is a utility function useful when summarizing an ordinal
response variable. It converts such a variable having k
levels to a
matrix with k-1
columns, where column i
is a vector of zeros and
ones indicating that the categorical response is in level i+1
or
greater. When the left hand side of formula
is cumcategory(y)
,
the default fun
will summarize it by computing all of the relevant
cumulative proportions.
Functions conTestkw
, catTestchisq
, ordTestpo
are
the default statistical test functions for summary.formula
.
These defaults are: Wilcoxon-Kruskal-Wallis test for continuous
variables, Pearson chi-square test for categorical variables, and the
likelihood ratio chi-square test from the proportional odds model for
ordinal variables. These three functions serve also as templates for
the user to create her own testing functions that are self-defining in
terms of how the results are printed or rendered in LaTeX, or plotted.
## S3 method for class 'formula'
summary(formula, data=NULL, subset=NULL,
na.action=NULL, fun = NULL,
method = c("response", "reverse", "cross"),
overall = method == "response" | method == "cross",
continuous = 10, na.rm = TRUE, na.include = method != "reverse",
g = 4, quant = c(0.025, 0.05, 0.125, 0.25, 0.375, 0.5, 0.625,
0.75, 0.875, 0.95, 0.975),
nmin = if (method == "reverse") 100
else 0,
test = FALSE, conTest = conTestkw, catTest = catTestchisq,
ordTest = ordTestpo, ...)
## S3 method for class 'summary.formula.response'
x[i, j, drop=FALSE]
## S3 method for class 'summary.formula.response'
print(x, vnames=c('labels','names'), prUnits=TRUE,
abbreviate.dimnames=FALSE,
prefix.width, min.colwidth, formatArgs=NULL, markdown=FALSE, ...)
## S3 method for class 'summary.formula.response'
plot(x, which = 1, vnames = c('labels','names'), xlim, xlab,
pch = c(16, 1, 2, 17, 15, 3, 4, 5, 0), superposeStrata = TRUE,
dotfont = 1, add = FALSE, reset.par = TRUE, main, subtitles = TRUE,
...)
## S3 method for class 'summary.formula.response'
latex(object, title = first.word(deparse(substitute(object))), caption,
trios, vnames = c('labels', 'names'), prn = TRUE, prUnits = TRUE,
rowlabel = '', cdec = 2, ncaption = TRUE, ...)
## S3 method for class 'summary.formula.reverse'
print(x, digits, prn = any(n != N), pctdig = 0,
what=c('%', 'proportion'),
npct = c('numerator', 'both', 'denominator', 'none'),
exclude1 = TRUE, vnames = c('labels', 'names'), prUnits = TRUE,
sep = '/', abbreviate.dimnames = FALSE,
prefix.width = max(nchar(lab)), min.colwidth, formatArgs=NULL, round=NULL,
prtest = c('P','stat','df','name'), prmsd = FALSE, long = FALSE,
pdig = 3, eps = 0.001, ...)
## S3 method for class 'summary.formula.reverse'
plot(x, vnames = c('labels', 'names'), what = c('proportion', '%'),
which = c('both', 'categorical', 'continuous'),
xlim = if(what == 'proportion') c(0,1)
else c(0,100),
xlab = if(what=='proportion') 'Proportion'
else 'Percentage',
pch = c(16, 1, 2, 17, 15, 3, 4, 5, 0), exclude1 = TRUE,
dotfont = 1, main,
prtest = c('P', 'stat', 'df', 'name'), pdig = 3, eps = 0.001,
conType = c('dot', 'bp', 'raw'), cex.means = 0.5, ...)
## S3 method for class 'summary.formula.reverse'
latex(object, title = first.word(deparse(substitute(object))), digits,
prn = any(n != N), pctdig = 0, what=c('%', 'proportion'),
npct = c("numerator", "both", "denominator", "slash", "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 = 0.001,
auxCol = NULL, dotchart=FALSE, ...)
## S3 method for class 'summary.formula.cross'
print(x, twoway = nvar == 2, prnmiss = any(stats$Missing > 0), prn = TRUE,
abbreviate.dimnames = FALSE, prefix.width = max(nchar(v)),
min.colwidth, formatArgs = NULL, ...)
## S3 method for class 'summary.formula.cross'
latex(object, title = first.word(deparse(substitute(object))),
twoway = nvar == 2, prnmiss = TRUE, prn = TRUE,
caption=attr(object, "heading"), vnames=c("labels", "names"),
rowlabel="", ...)
stratify(..., na.group = FALSE, shortlabel = TRUE)
## S3 method for class 'summary.formula.cross'
formula(x, ...)
cumcategory(y)
conTestkw(group, x)
catTestchisq(tab)
ordTestpo(group, x)
formula |
An R formula with additive effects. For |
x |
an object created by |
y |
a numeric, character, category, or factor vector for |
drop |
logical. If |
data |
name or number of a data frame. Default is the current frame. |
subset |
a logical vector or integer vector of subscripts used to specify the subset of data to use in the analysis. The default is to use all observations in the data frame. |
na.action |
function for handling missing data in the input data. The default is
a function defined here called |
fun |
function for summarizing data in each cell. Default is to take the
mean of each column of the possibly multivariate response variable.
You can specify |
method |
The default is The The |
overall |
For |
continuous |
specifies the threshold for when a variable is considered to be
continuous (when there are at least |
na.rm |
|
na.include |
for |
g |
number of quantile groups to use when variables are automatically
categorized with |
nmin |
if fewer than |
test |
applies if |
conTest |
a function of two arguments (grouping variable and a continuous
variable) that returns a list with components |
catTest |
a function of a frequency table (an integer matrix) that returns a
list with the same components as created by |
ordTest |
a function of a frequency table (an integer matrix) that returns a
list with the same components as created by |
... |
for |
object |
an object created by |
quant |
vector of quantiles to use for summarizing data with
|
vnames |
By default, tables and plots are usually labeled with variable labels
(see the |
pch |
vector of plotting characters to represent different groups, in order
of group levels. For |
superposeStrata |
If |
dotfont |
font for plotting points |
reset.par |
set to |
abbreviate.dimnames |
see |
prefix.width |
see |
min.colwidth |
minimum column width to use for boxes printed with |
formatArgs |
a list containing other arguments to pass to |
markdown |
for |
digits |
number of significant digits to print. Default is to use the current
value of the |
prn |
set to |
prnmiss |
set to |
what |
for |
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. |
npct |
specifies which counts are to be printed to the right of percentages.
The default is to print the frequency (numerator of the percent) in
parentheses. You can specify |
npct.size |
the size for typesetting |
Nsize |
When a second row of column headings is added showing sample sizes,
|
exclude1 |
by default, |
prUnits |
set to |
sep |
character to use to separate quantiles when printing
|
prtest |
a vector of test statistic components to print if |
round |
for |
prmsd |
set to |
msdsize |
defaults to |
long |
set to |
pdig |
number of digits to the right of the decimal place for printing
P-values. Default is |
eps |
P-values less than |
auxCol |
an optional auxiliary column of information, right justified, to add
in front of statistics typeset by
|
twoway |
for |
which |
For |
conType |
For plotting |
cex.means |
character size for means in box-percentile plots; default is .5 |
xlim |
vector of length two specifying x-axis limits. For
|
xlab |
x-axis label |
add |
set to |
main |
a main title. For |
subtitles |
set to |
caption |
character string containing LaTeX table captions. |
title |
name of resulting LaTeX file omitting the |
trios |
If for |
rowlabel |
see |
cdec |
number of decimal places to the right of the decimal point for
|
ncaption |
set to |
i |
a vector of integers, or character strings containing variable names
to subset on. Note that each row subsetted on in an |
j |
a vector of integers representing column numbers |
middle.bold |
set to |
outer.size |
the font size for outer quantiles for |
insert.bottom |
set to |
dcolumn |
see |
na.group |
set to |
shortlabel |
set to |
dotchart |
set to |
group |
for |
tab |
for |
summary.formula
returns a data frame or list depending on
method
. plot.summary.formula.reverse
returns the number
of pages of plots that were made.
plot.summary.formula.reverse
creates a function Key
and
Key2
in frame 0 that will draw legends.
Frank Harrell
Department of Biostatistics
Vanderbilt University
fh@fharrell.com
Harrell FE (2007): Statistical tables and plots using S and LaTeX. Document available from https://hbiostat.org/R/Hmisc/summary.pdf.
mChoice
, smean.sd
, summarize
,
label
, strata
, dotchart2
,
print.char.matrix
, update
,
formula
, cut2
, llist
,
format.default
, latex
,
latexTranslate
bpplt
,
summaryM
, summary
options(digits=3)
set.seed(173)
sex <- factor(sample(c("m","f"), 500, rep=TRUE))
age <- rnorm(500, 50, 5)
treatment <- factor(sample(c("Drug","Placebo"), 500, rep=TRUE))
# Generate a 3-choice variable; each of 3 variables has 5 possible levels
symp <- c('Headache','Stomach Ache','Hangnail',
'Muscle Ache','Depressed')
symptom1 <- sample(symp, 500,TRUE)
symptom2 <- sample(symp, 500,TRUE)
symptom3 <- sample(symp, 500,TRUE)
Symptoms <- mChoice(symptom1, symptom2, symptom3, label='Primary Symptoms')
table(Symptoms)
# Note: In this example, some subjects have the same symptom checked
# multiple times; in practice these redundant selections would be NAs
# mChoice will ignore these redundant selections
#Frequency table sex*treatment, sex*Symptoms
summary(sex ~ treatment + Symptoms, fun=table)
# could also do summary(sex ~ treatment +
# mChoice(symptom1,symptom2,symptom3), fun=table)
#Compute mean age, separately by 3 variables
summary(age ~ sex + treatment + Symptoms)
f <- summary(treatment ~ age + sex + Symptoms, method="reverse", test=TRUE)
f
# trio of numbers represent 25th, 50th, 75th percentile
print(f, long=TRUE)
plot(f)
plot(f, conType='bp', prtest='P')
bpplt() # annotated example showing layout of bp plot
#Compute predicted probability from a logistic regression model
#For different stratifications compute receiver operating
#characteristic curve areas (C-indexes)
predicted <- plogis(.4*(sex=="m")+.15*(age-50))
positive.diagnosis <- ifelse(runif(500)<=predicted, 1, 0)
roc <- function(z) {
x <- z[,1];
y <- z[,2];
n <- length(x);
if(n<2)return(c(ROC=NA));
n1 <- sum(y==1);
c(ROC= (mean(rank(x)[y==1])-(n1+1)/2)/(n-n1) );
}
y <- cbind(predicted, positive.diagnosis)
options(digits=2)
summary(y ~ age + sex, fun=roc)
options(digits=3)
summary(y ~ age + sex, fun=roc, method="cross")
#Use stratify() to produce a table in which time intervals go down the
#page and going across 3 continuous variables are summarized using
#quartiles, and are stratified by two treatments
set.seed(1)
d <- expand.grid(visit=1:5, treat=c('A','B'), reps=1:100)
d$sysbp <- rnorm(100*5*2, 120, 10)
label(d$sysbp) <- 'Systolic BP'
d$diasbp <- rnorm(100*5*2, 80, 7)
d$diasbp[1] <- NA
d$age <- rnorm(100*5*2, 50, 12)
g <- function(y) {
N <- apply(y, 2, function(w) sum(!is.na(w)))
h <- function(x) {
qu <- quantile(x, c(.25,.5,.75), na.rm=TRUE)
names(qu) <- c('Q1','Q2','Q3')
c(N=sum(!is.na(x)), qu)
}
w <- as.vector(apply(y, 2, h))
names(w) <- as.vector( outer(c('N','Q1','Q2','Q3'), dimnames(y)[[2]],
function(x,y) paste(y,x)))
w
}
#Use na.rm=FALSE to count NAs separately by column
s <- summary(cbind(age,sysbp,diasbp) ~ visit + stratify(treat),
na.rm=FALSE, fun=g, data=d)
#The result is very wide. Re-do, putting treatment vertically
x <- with(d, factor(paste('Visit', visit, treat)))
summary(cbind(age,sysbp,diasbp) ~ x, na.rm=FALSE, fun=g, data=d)
#Compose LaTeX code directly
g <- function(y) {
h <- function(x) {
qu <- format(round(quantile(x, c(.25,.5,.75), na.rm=TRUE),1),nsmall=1)
paste('{\\scriptsize(',sum(!is.na(x)),
')} \\hfill{\\scriptsize ', qu[1], '} \\textbf{', qu[2],
'} {\\scriptsize ', qu[3],'}', sep='')
}
apply(y, 2, h)
}
s <- summary(cbind(age,sysbp,diasbp) ~ visit + stratify(treat),
na.rm=FALSE, fun=g, data=d)
# latex(s, prn=FALSE)
## need option in latex to not print n
#Put treatment vertically
s <- summary(cbind(age,sysbp,diasbp) ~ x, fun=g, data=d, na.rm=FALSE)
# latex(s, prn=FALSE)
#Plot estimated mean life length (assuming an exponential distribution)
#separately by levels of 4 other variables. Repeat the analysis
#by levels of a stratification variable, drug. Automatically break
#continuous variables into tertiles.
#We are using the default, method='response'
## Not run:
life.expect <- function(y) c(Years=sum(y[,1])/sum(y[,2]))
attach(pbc)
require(survival)
S <- Surv(follow.up.time, death)
s2 <- summary(S ~ age + albumin + ascites + edema + stratify(drug),
fun=life.expect, g=3)
#Note: You can summarize other response variables using the same
#independent variables using e.g. update(s2, response~.), or you
#can change the list of independent variables using e.g.
#update(s2, response ~.- ascites) or update(s2, .~.-ascites)
#You can also print, typeset, or plot subsets of s2, e.g.
#plot(s2[c('age','albumin'),]) or plot(s2[1:2,])
s2 # invokes print.summary.formula.response
#Plot results as a separate dot chart for each of the 3 strata levels
par(mfrow=c(2,2))
plot(s2, cex.labels=.6, xlim=c(0,40), superposeStrata=FALSE)
#Typeset table, creating s2.tex
w <- latex(s2, cdec=1)
#Typeset table but just print LaTeX code
latex(s2, file="") # useful for Sweave
#Take control of groups used for age. Compute 3 quartiles for
#both cholesterol and bilirubin (excluding observations that are missing
#on EITHER ONE)
age.groups <- cut2(age, c(45,60))
g <- function(y) apply(y, 2, quantile, c(.25,.5,.75))
y <- cbind(Chol=chol,Bili=bili)
label(y) <- 'Cholesterol and Bilirubin'
#You can give new column names that are not legal S names
#by enclosing them in quotes, e.g. 'Chol (mg/dl)'=chol
s <- summary(y ~ age.groups + ascites, fun=g)
par(mfrow=c(1,2), oma=c(3,0,3,0)) # allow outer margins for overall
for(ivar in 1:2) { # title
isub <- (1:3)+(ivar-1)*3 # *3=number of quantiles/var.
plot(s3, which=isub, main='',
xlab=c('Cholesterol','Bilirubin')[ivar],
pch=c(91,16,93)) # [, closed circle, ]
}
mtext(paste('Quartiles of', label(y)), adj=.5, outer=TRUE, cex=1.75)
#Overall (outer) title
prlatex(latex(s3, trios=TRUE))
# trios -> collapse 3 quartiles
#Summarize only bilirubin, but do it with two statistics:
#the mean and the median. Make separate tables for the two randomized
#groups and make plots for the active arm.
g <- function(y) c(Mean=mean(y), Median=median(y))
for(sub in c("D-penicillamine", "placebo")) {
ss <- summary(bili ~ age.groups + ascites + chol, fun=g,
subset=drug==sub)
cat('\n',sub,'\n\n')
print(ss)
if(sub=='D-penicillamine') {
par(mfrow=c(1,1))
plot(s4, which=1:2, dotfont=c(1,-1), subtitles=FALSE, main='')
#1=mean, 2=median -1 font = open circle
title(sub='Closed circle: mean; Open circle: median', adj=0)
title(sub=sub, adj=1)
}
w <- latex(ss, append=TRUE, fi='my.tex',
label=if(sub=='placebo') 's4b' else 's4a',
caption=paste(label(bili),' {\\em (',sub,')}', sep=''))
#Note symbolic labels for tables for two subsets: s4a, s4b
prlatex(w)
}
#Now consider examples in 'reverse' format, where the lone dependent
#variable tells the summary function how to stratify all the
#'independent' variables. This is typically used to make tables
#comparing baseline variables by treatment group, for example.
s5 <- summary(drug ~ bili + albumin + stage + protime + sex +
age + spiders,
method='reverse')
#To summarize all variables, use summary(drug ~., data=pbc)
#To summarize all variables with no stratification, use
#summary(~a+b+c) or summary(~.,data=\dots)
options(digits=1)
print(s5, npct='both')
#npct='both' : print both numerators and denominators
plot(s5, which='categorical')
Key(locator(1)) # draw legend at mouse click
par(oma=c(3,0,0,0)) # leave outer margin at bottom
plot(s5, which='continuous')
Key2() # draw legend at lower left corner of plot
# oma= above makes this default key fit the page better
options(digits=3)
w <- latex(s5, npct='both', here=TRUE)
# creates s5.tex
#Turn to a different dataset and do cross-classifications on possibly
#more than one independent variable. The summary function with
#method='cross' produces a data frame containing the cross-
#classifications. This data frame is suitable for multi-panel
#trellis displays, although `summarize' works better for that.
attach(prostate)
size.quartile <- cut2(sz, g=4)
bone <- factor(bm,labels=c("no mets","bone mets"))
s7 <- summary(ap>1 ~ size.quartile + bone, method='cross')
#In this case, quartiles are the default so could have said sz + bone
options(digits=3)
print(s7, twoway=FALSE)
s7 # same as print(s7)
w <- latex(s7, here=TRUE) # Make s7.tex
library(trellis,TRUE)
invisible(ps.options(reset=TRUE))
trellis.device(postscript, file='demo2.ps')
dotplot(S ~ size.quartile|bone, data=s7, #s7 is name of summary stats
xlab="Fraction ap>1", ylab="Quartile of Tumor Size")
#Can do this more quickly with summarize:
# s7 <- summarize(ap>1, llist(size=cut2(sz, g=4), bone), mean,
# stat.name='Proportion')
# dotplot(Proportion ~ size | bone, data=s7)
summary(age ~ stage, method='cross')
summary(age ~ stage, fun=quantile, method='cross')
summary(age ~ stage, fun=smean.sd, method='cross')
summary(age ~ stage, fun=smedian.hilow, method='cross')
summary(age ~ stage, fun=function(x) c(Mean=mean(x), Median=median(x)),
method='cross')
#The next statements print real two-way tables
summary(cbind(age,ap) ~ stage + bone,
fun=function(y) apply(y, 2, quantile, c(.25,.75)),
method='cross')
options(digits=2)
summary(log(ap) ~ sz + bone,
fun=function(y) c(Mean=mean(y), quantile(y)),
method='cross')
#Summarize an ordered categorical response by all of the needed
#cumulative proportions
summary(cumcategory(disease.severity) ~ age + sex)
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.