As promised here is a try at a function that is like Trellis' xyplot function
but that allows for a "multivariate" response. The most common application is
adding confidence bars to group means. xYplot is now in the Hmisc library
(which will be updated on Statlib tonight or you can get it from the web page
listed below). xYplot allows for "error bars", alternating high- and low error
bars, and "bands". For uses wanting not having the Hmisc library who want
to give xYplot a try, set the default for the label.curve argument in
panel.xYplot
to F instead of T.
Example: xYplot(Cbind(ymean, lower,upper) ~ month | year, groups=treatment,
method='alt bars')
xYplot(Cbind(ymean,1.96*se) ~ month)
Thanks to all who gave me ideas. Bug reports/suggestions welcome.
---------------------------------------------------------------------------
Frank E Harrell Jr
Professor of Biostatistics and Statistics
Director, Division of Biostatistics and Epidemiology
Dept of Health Evaluation Sciences
University of Virginia School of Medicine
http://www.med.virginia.edu/medicine/clinical/hes/biostat.htm
Cbind <- function(...) { # See llist function with Hmisc label function
dotlist <- list(...)
lname <- names(dotlist)
name <- vname <- as.character(sys.call()[-1])
for(i in 1:length(dotlist)) {
vname[i] <- lname[i]
if(length(vname[i])==0 || vname[i]=='') vname[i] <- name[i]
}
lab <- attr(y <- dotlist[[1]],'label')
if(!length(lab)) lab <- vname[1]
other <- as.matrix(as.data.frame(dotlist))[,-1,drop=F]
dimnames(other)[[2]] <- vname[-1]
structure(y, class='Cbind', label=lab, other=other)
}
'[.Cbind' <- function(x, ...) {
structure(unclass(x)[...], class='Cbind',
label=attr(x,'label'),
other=attr(x,'other')[...,,drop=F])
}
prepanel.xYplot <- function(x, y, ...) {
xlim <- range(x)
ylim <- range(y, attr(y,'other'), na.rm=T)
list(xlim=xlim, ylim=ylim, dx=diff(xlim), dy=diff(ylim))
}
panel.xYplot <- function(x, y, subscripts, groups=NULL,
type='b',
method=c('bars','bands','upper bars','lower bars','alt bars'),
label.curves=T, cap=.015, lty.bar=1,
lwd = plot.line$lwd,
lty = plot.line$lty,
pch = plot.symbol$pch,
cex = plot.symbol$cex,
font= plot.symbol$font,
col = NULL, ...) {
method <- match.arg(method)
groups <- as.factor(groups)
other <- attr(y, 'other')
if(ncol(other)==1) {
lower <- y - other
upper <- y + other
} else {
lower <- other[,1]
upper <- other[,2]
}
y <- unclass(y)
g <- unclass(groups)[subscripts]
ng <- if(length(groups)) max(g) else 1
plot.symbol <- trellis.par.get(if(ng>1)"superpose.symbol" else "plot.symbol")
plot.line <- trellis.par.get(if(ng>1)"superpose.line" else "plot.line")
lty <- rep(lty, length = ng)
lwd <- rep(lwd, length = ng)
pch <- rep(pch, length = ng)
cex <- rep(cex, length = ng)
font <- rep(font,length = ng)
if(!length(col)) col <- if(type=='p') plot.symbol$col else plot.line$col
col <- rep(col, length = ng)
if(ng > 1) {
levnum <- sort(unique(g))
panel.superpose(x, y, subscripts, groups,
lwd=lwd, lty=lty, pch=pch, cex=cex,
font=font, col=col, type=type)
## Note potential problem - panel.superpose does not handle type='b'
if(type!='p' && !(is.logical(label.curves) && !label.curves)) {
lc <- if(is.logical(label.curves)) list(lwd=lwd, cex=cex[1]) else
c(list(lwd=lwd, cex=cex[1]), label.curves)
curves <- vector('list',length(levnum)); names(curves) <- levels(groups)
i <- 0
for(gg in levnum) {
i <- i+1
s <- g==gg
curves[[i]] <- list(x[s], y[s])
}
labcurve(curves, lty=lty[levnum], lwd=lwd[levnum], col=col[levnum], opts=lc,
...)
}
} else panel.xyplot(x, y,
lwd=lwd, lty=lty, pch=pch, cex=cex,
font=font, col=col, type=type)
if(method=='bands') {
for(j in 1:ncol(other)) {
if(ng==1) panel.xyplot(x, other[,j],
lwd=lwd, lty=lty, pch=pch, cex=cex,
font=font, col=col, type='l') else
panel.superpose(x, other[,j], subscripts, groups,
lwd=lwd, lty=lty, pch=pch, cex=cex,
font=font, col=col, type='l')
}
} else {
errbr <- function(x, y, lower, upper, cap, lty, lwd, col, connect) {
smidge <- .5*cap*diff(par('usr')[1:2])
switch(connect,
all={segments(x, lower, x, upper, lty = lty, lwd = lwd, col=col)
segments(x - smidge, lower, x + smidge, lower, lwd = lwd, lty=1, col=col)
segments(x - smidge, upper, x + smidge, upper, lwd = lwd, lty=1, col=col)},
upper={segments(x, y, x, upper, lty=lty, lwd=lwd, col=col)
segments(x - smidge, upper, x + smidge, upper, lwd=lwd, lty=1, col=col)},
lower={segments(x, y, x, lower, lty=lty, lwd=lwd, col=col)
segments(x - smidge, lower, x + smidge, lower, lwd=lwd, lty=1, col=col)})
}
if(ng==1) errbr(x, y, lower, upper, cap, lty.bar, lwd, col,
switch(method,
bars='all','upper bars'='upper','lower bars'='lower','alt bars'='lower'))
else {
if(method=='alt bars') medy <- median(y, na.rm=T)
for(gg in levnum) {
s <- g==gg
connect <- switch(method,
bars='all','upper bars'='upper','lower bars'='lower',
'alt bars'=if(median(y[s],na.rm=T) > medy)'upper' else 'lower')
errbr(x[s], y[s], lower=lower[s], upper=upper[s], cap,
lty.bar, lwd[gg], col[gg], connect)
}
}
}
if(type!='l' && ng>1) { ##set up for key() if points plotted
Key <- function(x, y, lev, cex, col, font, pch, ...) {
if(!missing(x)) {
if(is.list(x)) {y <- x$y; x <- x$x}
key(x=x, y=y, text=list(lev, col=col),
points=list(cex=cex,col=col,font=font,pch=pch),
transparent=T, ...) } else
key(text=list(lev, col=col),
points=list(cex=cex,col=col,font=font,pch=pch),transparent=T, ...)
invisible()
}
Key$lev <- levels(groups)
Key$cex <- cex; Key$col <- col; Key$font <- font; Key$pch <- pch
assign('Key', Key, frame=0)
}
}
xYplot <- function(formula, data = sys.parent(1),
groups = NULL,
prepanel=prepanel.xYplot, panel=panel.xYplot, ...,
ylab=attr(eval(formula[[2]], data),'label'), subset=T) {
subset <- eval(substitute(subset), data)
setup.2d.trellis(formula, data = data,
prepanel=prepanel, panel=panel,
groups = eval(substitute(groups), data)[subset], ...,
ylab=ylab, subset = subset)
}
DESCRIPTION
A utility function Cbind returns the first argument as a vector and combines
all other arguments into a matrix stored as an
attribute called "other". The arguments can be named (e.g.,
Cbind(pressure=y,ylow,yhigh)) or a label attribute may be pre-attached
to the first argument. In either case, the name or label of the first argument
is stored as an attribute "label" of the object
returned by Cbind. Storing other vectors as a matrix attribute facilitates
plotting error bars, etc., as trellis really wants the
response variable to be a vector, not a matrix. A subscript method for Cbind
objects subscripts the other matrix along with the
main y vector.
The xYplot function is a substitute for xyplot that allows for simulated
multi-column y. It uses by default the panel.xYplot and
prepanel.xYplot functions to do the actual work. The method argument passed to
panel.xYplot from xYplot allows you to make error
bars, the upper-only or lower-only portions of error bars, alternating
lower-only and upper-only bars, or bands. panel.xYplot
decides how to alternate upper and lower bars according to whether the median y
value of the current main data line is above the
median y for all groups of lines or not. If the median is above the overall
median, only the upper bar is drawn. For bands, any
number of other columns of y will be drawn as lines having the same thickness,
color, and type as the main data line. If plotting
bars and only one additional column is specified for the response variable,
that column is taken as the half width of a precision
interval for y, and the lower and upper values are computed automatically as y
plus or minus the value of the additional column
variable.
When a groups variable is present, panel.xYplot will create a function in frame
0 called Key that when invoked will draw a key
describing the groups labels, point symbols, and colors. By default, the key is
outside the graph. If Key(locator(1)) is specified,
the key will appear so that its upper left corner is at the coordinates of the
mouse click.
Unlike xyplot, xYplot senses the presence of a groups variable and
automatically invokes panel.superpose instead of panel.xyplot.
USAGE
Cbind(y, ...)
xYplot(formula, data = sys.parent(1), groups = NULL,
prepanel=prepanel.xYplot, panel=panel.xYplot, ...,
ylab=attr(eval(formula[[2]], data),'label'), subset=T)
panel.xYplot(x, y, subscripts, groups=NULL, type="b",
method=c("bars", "bands", "upper bars", "lower bars", "alt bars"),
label.curves=T, cap=0.015, lty.bar=1,
lwd=plot.line$lwd, lty=plot.line$lty, pch=plot.symbol$pch,
cex=plot.symbol$cex, font=plot.symbol$font, col=NULL, ...)
prepanel.xYplot(x, y, ...)
REQUIRED ARGUMENTS
y a vector representing the main variable to plot, i.e., the variable used to
draw the main lines
... for Cbind ... can be any number of additional numeric vectors. Unless you
are using method="bands", vectors after the first two
are ignored. If drawing bars and only one extra variable is given in ...,
upper and lower values are computed as described above.
formula a trellis formula consistent with xyplot
x, y x and y-axis variables
OPTIONAL ARGUMENTS
x x-axis variable
y an object created by Cbind
subscripts, groups, type see trellis.args
method defaults to "bars" to draw error-bar type plots. See meaning of other
values above.
label.curves set to F to suppress invocation of labcurve to label primary
curves where they are most separated or to draw a legend
in an empty spot on the panel. You can also set label.curves to a list of
options to pass to labcurve. These options can also be
passed as ... to xYplot. See the examples below.
cap the half-width of horizontal end pieces for error bars, as a fraction of
the length of the x-axis
lty.bar line type for bars
lwd, lty, pch, cex, font, col see trellis.args. These are vectors when groups
is present.
... other arguments to pass to labcurve or Key
VALUE
Cbind returns a matrix with attributes. Other functions return standard
trellis results.
SIDE EFFECTS
plots, and panel.xYplot creates the Key function in the session frame.
AUTHOR
Frank Harrell
Division of Biostatistics and Epidemiology
Department of Health Evaluation Sciences
University of Virginia Health Sciences Center
fharrell@virginia.edu
REFERENCES
SEE ALSO
xyplot, panel.xyplot, trellis, label, labcurve, error.bar
EXAMPLES
dfr <- expand.grid(month=1:12, continent=c('Europe','USA'),
sex=c('female','male'))
attach(dfr)
set.seed(13)
y <- month/10 + 1*(sex=='female') + 2*(continent=='Europe') +
runif(48,-.15,.15)
lower <- y - runif(48,.05,.15)
upper <- y + runif(48,.05,.15)
xYplot(Cbind(y,lower,upper) ~ month,subset=sex=='male'&continent=='USA')
# add ,label.curves=F to suppress use of labcurve to label curves where
farthest apart
xYplot(Cbind(y,lower,upper) ~ month|continent,subset=sex=='male')
xYplot(Cbind(y,lower,upper) ~ month|continent,groups=sex); Key()
xYplot(Cbind(y,lower,upper) ~ month,groups=sex,subset=continent=='Europe')
xYplot(Cbind(y,lower,upper) ~
month,groups=sex,subset=continent=='Europe',keys='lines')
# keys='lines' causes labcurve to draw a legend where the panel is most empty
xYplot(Cbind(y,lower,upper) ~
month,groups=sex,subset=continent=='Europe',method='bands')
xYplot(Cbind(y,lower,upper) ~
month,groups=sex,subset=continent=='Europe',method='upper')
label(y) <- 'Quality of Life Score'
# label is in Hmisc library = attr(y,'label') <- 'Quality...'; will be y-axis
label
# can also specify Cbind('Quality of Life Score'=y,lower,upper)
xYplot(Cbind(y,lower,upper) ~
month,groups=sex,subset=continent=='Europe',method='alt bars',
offset=.4) # offset passed to labcurve to label .4 y units away from curve
-----------------------------------------------------------------------
This message was distributed by s-news@wubios.wustl.edu. To unsubscribe
send e-mail to s-news-request@wubios.wustl.edu with the BODY of the
message: unsubscribe s-news
|