Nothing
## $Id: utilities.R 6053 2015-02-22 20:23:45Z bpikouni $
## Functions & objects internal for the cg Library
## and definitely not intended for user-level calls
#####################
stndErr <- function(x) {
##
## PURPOSE: Calculate standard error measure from a group of data.
## NOTE: We use of "n" in the function name
## instead of "stderr" because SV4 "stderr"
## is a system base function that is important in connections.
x <- x[!is.na(x)]
sd(x) / sqrt(length(x))
}
#####################
geoMean <- function(x, na.rm=TRUE) {
##
## PURPOSE: Calculate the geometric mean from a group of data
if(na.rm) {
x <- x[!is.na(x)]
}
if(any(x <= 0)) {
stop(cgMessage("There appear to be some non-positive",
"values in the data. All values",
"must be greater than zero."))
}
exp(mean(log(x)))
}
#####################
pctToRatio <- function(x) { (x/100) + 1 }
#####################
makeZeroScore <- function(x) {
respvec <- unique(stripmiss(x[x > 0]))
return(exp(spline(respvec, log(respvec), method="natural",
xmin=0, n=length(respvec) + 1)$y[1]))
}
#####################
unwind <- function(x) {
as.vector(unlist(x))
}
unwrap <- unwind
#####################
graphStampCG <- function(grid=TRUE) {
##
## PURPOSE: ID Stamp on graphs
##
validBoolean(grid)
## Reads off package object "cgVersion"
dateitems <- unlist(strsplit(date(), split=" "))
## Remove extraneous string if present due to single digit days
dateitems <- dateitems[dateitems!=""]
## sequence: daynumber, month, year, time, weekday
newdate <- paste(dateitems[3]," ",dateitems[2]," ",
dateitems[5]," \n",dateitems[4]," local \n",sep="")
thetext <- paste("cg",
" v", cgVersion, " \n",
newdate,"", sep="")
if(grid) {
grid::grid.text(thetext, x = 0.99, y = 0.01,
default.units = "npc",
just = c("right", "bottom"),
gp=gpar(cex=0.6))
}
else { ## not grid
mtext(text=thetext,
side = 1,
line = -1, outer = TRUE,
adj = 1, cex = 0.6)
}
invisible()
}
#####################
setupAxisTicks <- function(x, ratio=FALSE, percent=FALSE,
difference=FALSE, axis="y",
logscale=TRUE, digits=NULL, remticks=TRUE,
offset=NULL, grid=FALSE, xcex=1, ycex=1) {
##
## PURPOSE: Try to setup good ticks for log spacing
## or original scale spacing:
## Ensure no overlap, and
## Ensure space for later addition of min and max
## NOTE: x is assumed to span the whole axis range of data
##
if(is.null(digits)) digits <- 4
x <- x[!is.na(x)]
if (!is.null(offset)) {
if(logscale==FALSE) {
stop(cgMessage("The offset argument must be set to NULL when",
"logscale=FALSE, as they are not compatible."))
}
x <- (x - offset)
}
if(any(x < 0)) {
## generating ticks on the magnitude of values
x.pos <- abs(x)
}
else {
x.pos <- x
}
log10endpt <- log10(x.pos)
if(ratio==TRUE) {
## to impose that zero change is always represented
log10endpt <- c(0, log10endpt)
}
digitplacement <- floor(max(log10endpt))
log10ticks <- pretty(log10endpt)
tickmarks <- round(10^log10ticks, -digitplacement)
## need to filter extreme values for case where logscaling
## essentially has no effect, i.e. the range is clearly small.
tickmarks <- tickmarks[tickmarks >= round(min(x), digits) &
tickmarks <= round(max(x), digits)]
## Also ensure that powers of 10 from 1 thru range are plotted
powersof10 <- 10^(0:digitplacement)
## Additional adjustments for Ratio/Percent Change
if(ratio==TRUE) {
powersof10 <- (powersof10/10) + 1
}
powersof10 <- powersof10[powersof10 >= round(min(x), digits) &
powersof10 <= round(max(x), digits)]
## and that some interim values that are familiar are there
if(ratio==TRUE) {
interimticks <- as.vector(outer(c(-9.9,-9,-7.5,-5,-2.5,-1,1,2.5,5),
10^(0:digitplacement)))
## Additional adjustments for percent change to ratio
interimticks <- (interimticks/10) + 1
}
else {
interimticks <- as.vector(outer(c(1,2,5),10^(0:digitplacement)))
}
## use only those that fall within range
interimticks <- interimticks[interimticks >= (round(min(x), digits) - 1e-8) &
interimticks <= (round(max(x), digits) + 1e-8) ]
## If Ratio, we also need to include the
## "no-difference" reference which is 1 on the ratio scale
## Also need to make sure extremes of axes are included
if(ratio==TRUE) {
## Cover the two extremes of the axes via percent scaling because we can
## work with integers
ends <- c(floor(100*(min(x.pos) - 1)),
round(100*(min(x.pos) - 1),0),
round(100*(max(x.pos) - 1),0),
ceiling(100*(max(x.pos) - 1)))
## if 0 at lower end replace with 1
ends[ends==0] <- 1
## transform to ratio
ends <- (ends/100) + 1
xticks <- c(1, tickmarks, powersof10, interimticks, ends,
pretty(x.pos))
}
else {
## Cover the two extremes of the axes
ends <- round(range(x), digits)
xticks <- c(tickmarks, powersof10, interimticks, ends,
pretty(x))
## pretty(x.pos))
}
## Remove 0 if it has crept in there if we plan use for log-spacing
if(logscale==TRUE) {
xticks <- xticks[xticks > 0]
}
## Annotation in case the original data had negative values
## Except do not reflect the min and max.
if(any(x < 0)) {
xticks <- c(-xticks, 0, xticks)
xticks <- xticks[!(xticks==round(-min(x), digits) | xticks==round(-max(x),
digits))]
}
## If ratio, then try to ensure symmetry of tickmarks in the range of 0.1 to
## 1.9
if(ratio==TRUE) {
upperbound <- max(min(round(x.pos,1), 1.9), 1.1)
lowerbound <- min(min(round(x.pos,1)), 0.9)
if(upperbound < 1.8) upperbound <- upperbound + 0.2
if(lowerbound > 0.2) lowerbound <- lowerbound - 0.2
xticks <- c(seq(lowerbound, upperbound, by=0.1), xticks)
}
## If difference, then try to ensure symmetry of tickmarks
if(difference==TRUE) {
xticks <- c(0, -xticks, xticks)
}
## If ratio, then try to ensure symmetry of tickmarks for smaller ranges
if(ratio==TRUE) {
xticks <- c(xticks[xticks < max(xticks)] - 1, xticks)
xticks <- xticks[xticks > 0]
}
## If ratio, round off values > 2 to integers
if(ratio==TRUE) {
xticks[xticks > 2] <- round(xticks[xticks > 2], 0)
}
## Exclude anything outside the min or max of original data, Note that 1e-8 is
## added as padding to insure that min and max are included in the comparison
## but that nothing else of consequence will be.
xticks <- xticks[xticks >= (round(min(x), digits) - 1e-8) &
xticks <= (round(max(x), digits) + 1e-8)]
## Take out duplicates and ensure sorted order
xticks <- sort(unique(xticks))
## Empty case handling
if(length(xticks)==0) {
xticks <- 1
names(xticks) <- ""
return(xticks)
}
## Format tick labels
## Try to remove extraneous decimal places for larger numbers
if(ratio==FALSE && axis=="y") {
names(xticks) <- format(chopZeroes(round(xticks, digits)),
small.interval=0, justify="right")
}
else if(ratio==FALSE && axis=="x") {
names(xticks) <- format(chopZeroes(round(xticks, digits)),
small.interval=0, justify="right")
}
else if(percent==TRUE && axis=="x") {
names(xticks) <- fmtRatioToPercent(xticks)
}
else if(percent==TRUE && axis=="y") {
names(xticks) <- fmtRatioToPercent(xticks)
}
else if(difference==TRUE && axis=="x") {
names(xticks) <- fmtDifference(xticks)
}
else if(difference==TRUE && axis=="y") {
names(xticks) <- format(chopZeroes(round(xticks, digits)), justify="right")
}
else if(ratio==TRUE) {
names(xticks) <- fmtRatio(xticks)
}
names(xticks) <- trimWhiteSpace(names(xticks))
## Remove duplicate *labels* (which may have come from rounding digits)
xticks <- xticks[which(!duplicated(names(xticks)))]
if(!is.null(offset)) { xticks <- (xticks + offset) }
if(remticks==TRUE) {
## Remove overlapping ticks
xticks <- rmTicks(xticks, axis=axis, logscale=logscale,
difference=difference, percent=percent,
ratio=ratio,
digits=digits, grid=grid, xcex=xcex, ycex=ycex)
}
return(xticks)
}
#####################
setupLog10AxisTicks <- function(endpt, offset) {
##
## PURPOSE: Create log110 tickmarks intended for
## axis = 4 and that can accomodate log(x + addconstant)
## transformed data
##
prettyticks.log10endpt <- if(is.null(offset)) {
pretty(log10(endpt))
}
else {
pretty(log10(endpt - offset))
}
ticks.log10endpt <- if(is.null(offset)){
prettyticks.log10endpt
}
else {
log10(10^prettyticks.log10endpt + offset)
}
names(ticks.log10endpt) <- prettyticks.log10endpt
return(ticks.log10endpt)
}
#####################
tryAgain <- paste("Please check your argument inputs or/and data, correct as needed, and",
"try again. Typing in traceback() next may also help.")
#####################
seeHelpFile <- function(x) {
return(paste("The help file for", x, "may provide details."))
}
#####################
paragraphWrap <- function (x, width = 80) {
##
## PURPOSE: Format a string for use in stop() or warning().
##
## NOTE: Adapted from R function base:::strwrap
##
indent <- 0
exdent <- 0
prefix <- ""
indentString <- paste(rep(" ", indent), collapse = "")
exdentString <- paste(rep(" ", exdent), collapse = "")
y <- list()
z <- lapply(strsplit(x, "\n[ \n]*\n"), strsplit, "[ \n]")
for (i in seq(along = z)) {
yi <- character(0)
for (j in seq(along = z[[i]])) {
words <- z[[i]][[j]]
nc <- nchar(words)
if (any(nc == 0)) {
zLenInd <- which(nc == 0)
zLenInd <- zLenInd[!(zLenInd %in% (grep("\\.$", words) + 1))]
if (length(zLenInd) > 0) {
words <- words[-zLenInd]
nc <- nc[-zLenInd]
}
}
if (length(words) == 0) {
yi <- c(yi, "", prefix)
next
}
currentIndex <- 0
lowerBlockIndex <- 1
upperBlockIndex <- integer(0)
lens <- cumsum(nc + 1)
first <- TRUE
maxLength <- width - nchar(prefix) - indent
while (length(lens) > 0) {
k <- max(sum(lens < maxLength), 1)
if (first) {
first <- FALSE
maxLength <- maxLength + indent - exdent
}
currentIndex <- currentIndex + k
if (nc[currentIndex] == 0)
upperBlockIndex <- c(upperBlockIndex, currentIndex -
1)
else upperBlockIndex <- c(upperBlockIndex, currentIndex)
if (length(lens) > k) {
if (nc[currentIndex + 1] == 0) {
currentIndex <- currentIndex + 1
k <- k + 1
}
lowerBlockIndex <- c(lowerBlockIndex, currentIndex +
1)
}
if (length(lens) > k)
lens <- lens[-(1:k)] - lens[k]
else lens <- NULL
}
nBlocks <- length(upperBlockIndex)
s <- paste(prefix, c(indentString, rep(exdentString,
nBlocks - 1)), sep = "")
for (k in (1:nBlocks))
s[k] <- paste(s[k],
paste(words[lowerBlockIndex[k]:upperBlockIndex[k]],
collapse = " "), sep = "")
yi <- c(yi, s, prefix)
}
y <- c(y, list(yi[-length(yi)]))
}
## We always want it " simplify'ed " (see base::strwrap)
y <- unlist(y)
## Next to lastly we add the functionality of putting it into a string with
## linebreaks for its intended use in stop() or warning()
y <- paste(y, collapse="\n")
## and add a blank space after periods followed by an upper case characters
y <- gsub('\\.([A-Z])', '\\.\\ \\1', y, perl=TRUE)
## Add a newline so the first line is flush left when shown
## and another at the end.
return(paste("\n", y, "\n",sep=""))
}
#####################
cgMessage <- function(..., warning=FALSE) {
paragraphWrap(paste(..., ifelse(warning==FALSE, tryAgain, "")))
}
#####################
factorInSeq <- function(x, labels) {
##
## PURPOSE: Laziness motivated wrapper for creating
## a factor with the order of its levels determined
## by the data sequence
##
inseqlevels <- unique(x)
if(missing(labels)) {
labels <- as.character(inseqlevels)
}
factor(x, levels=unique(x), labels=labels)
}
#####################
setupGrpNameTicks <- function(x, grplocation=1:length(x),
cexinit=1, cexthreshold=0.80, grid=FALSE) {
##
## PURPOSE: Handle group name labels so there is no overlap
## Either shrink font, or rotate 45 degrees (or both)
##
## NOTE: For traditional graphics,
## Assume this is to be used with a graph that has just been
## plotted, since we need critical values from the device, which
## are given in inches.
##
## x assumed to be a character vector already in preferred order
##
if(!grid) { ## Traditional Graphics
axislength <- par("pin")[1]
axismin <- par("usr")[1]
axisrange <- par("usr")[2] - par("usr")[1]
## Since widths will vary
chrspace <- par("cin")[1]
chrspaces <- chrspace*nchar(x)
}
else { ## grid, work with npc units
usrx <- convertX(unit(0:1, "npc"), "native", TRUE)
usry <- convertY(unit(0:1, "npc"), "native", TRUE)
## grid lattice natural coordinates for panel
axislength <- 1
chrspaces <- convertX(stringWidth(x), "npc", TRUE)
chrspace <- max(convertX(stringWidth(letters), "npc", TRUE))
axismin <- usrx[1]
axisrange <- usrx[2] - usrx[1]
}
## Need to map response scale to device inches scale
n <- length(x)
p <- (axislength/axisrange) * (grplocation - axismin)
## Construct comparisons to determine overlap
ppart <- p[-1] - p[-n]
chrpart <- (chrspaces[-n] + chrspaces[-1])/2
comparisons <- ppart - chrpart
## Initializations
thecex <- cexinit
thesrt <- 0
theadj <- 0.5
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
while(any(overlapchecks==FALSE) && thecex >= cexthreshold) {
thecex <- thecex - 0.05
chrspace <- chrspace*thecex
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
}
if(thecex < cexthreshold) {
thesrt <- 45
theadj <- 1
## Now need to check overlap in the rotated scale, which
## can be done assuming 90 degree rotation, since
## character height now serves as the magic overlap msmt.
chrspace <-
if(!grid) {
par("cin")[2]
}
else {
max(convertX(stringHeight(letters), "inches", TRUE))
}
## we only need to worry about 1 character since
## they labels are virtually perpendicular to the axis
chrspaces <- rep(chrspace, length(x))
chrpart <- (chrspaces[-n] + chrspaces[-1])/2
comparisons <- ppart - chrpart
thecex <- 1
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
while(any(overlapchecks==FALSE) && thecex >= cexthreshold) {
thecex <- thecex - 0.05
chrspace <- chrspace*thecex
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
}
}
if(!grid) {
## If the labels are to be rotated, then we need to make sure they are
## shrunk within the margin available.
## Only worry about this for traditional graphics
if(thesrt!=0) {
themaxlength <- max(strwidth(x,
units="inches"))/sqrt(2)
thecex <- min(thecex, (par("mai")[1]/themaxlength)*thecex)
}
}
return(list(cex=thecex, srt=thesrt, adj=theadj))
}
xTicksCex <- function(x, thenames=names(x),
cexinit=1, cexthreshold=0.70,
grid=FALSE) {
##
##
thenames <- names(x)
if(!grid) { ## Traditional Graphics
axislength <- par("pin")[1]
axismin <- par("usr")[1]
axisrange <- par("usr")[2] - par("usr")[1]
## Since widths will vary
chrspace <- par("cin")[1]
chrspaces <- chrspace*nchar(thenames)
}
else { ## grid lattice natural coordinates for panel
usrx <- convertX(unit(0:1, "npc"), "native", TRUE)
axislength <- 1
axismin <- usrx[1]
axisrange <- usrx[2] - usrx[1]
chrspaces <- convertX(stringWidth(thenames), "npc", TRUE)
chrspace <- max(convertX(stringWidth(letters), "npc", TRUE))
}
## Need to map response scale to device inches scale for traditional,
## npc for grid
n <- length(x)
p <- (axislength/axisrange) * (x - axismin)
## Construct comparisons to determine overlap
ppart <- p[-1] - p[-n]
chrpart <- (chrspaces[-n] + chrspaces[-1])/2
comparisons <- ppart - chrpart
## Initializations
thecex <- cexinit
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
while(any(overlapchecks==FALSE) && thecex > cexthreshold) {
thecex <- thecex - 0.05
chrspace <- chrspace*thecex
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
}
return(thecex)
}
yTicksCex <- function(y,
cexinit=1, cexthreshold=0.70,
grid=FALSE) {
##
##
thenames <- names(y)
if(!grid) { ## Traditional Graphics
axislength <- par("pin")[2]
axismin <- par("usr")[2]
axisrange <- par("usr")[4] - par("usr")[3]
chrspace <- par("cin")[2]
chrspaces <- chrspace*nchar(thenames)
}
else {
usry <- convertY(unit(0:1, "npc"), "native", TRUE)
## grid lattice natural coordinates for panel
axislength <- 1
axismin <- usry[1]
axisrange <- usry[2] - usry[1]
chrspaces <- convertY(stringHeight(thenames), "npc", TRUE)
chrspace <- max(convertY(stringHeight(letters), "npc", TRUE))
}
## Need to map response scale to device inches scale
n <- length(y)
p <- (axislength/axisrange) * (y - axismin)
## Construct comparisons to determine overlap
ppart <- p[-1] - p[-n]
chrpart <- (chrspaces[-n] + chrspaces[-1])/2
comparisons <- ppart - chrpart
## Initializations
thecex <- cexinit
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
while(any(overlapchecks==FALSE) && thecex > cexthreshold) {
thecex <- thecex - 0.05
chrspace <- chrspace*thecex
overlapchecks <- c(TRUE, comparisons > chrspace/1.2)
}
return(thecex)
}
#####################
rmTicks <- function(x, axis="y", logscale=TRUE, ratio=FALSE, percent=FALSE,
difference=FALSE, digits=NULL, grid=FALSE, xcex=1, ycex=1) {
##
## PURPOSE: Try to prevent overlapping tickmarks
##
## NOTE: Assume this is to be used with a graph that has just been
## plotted, since we need critical values from the device, which
## are given in inches.
##
x.orig <- x
x <- sort(x)
if(logscale) {
x <- log10(x)
}
if(!grid) { ## Traditional Graphics
if(axis=="x") {
chrspace <- par("cin")[1]
axislength <- par("pin")[1]
axismin <- par("usr")[1]
axisrange <- par("usr")[2] - par("usr")[1]
chrspaces <- chrspace*nchar(names(x))
thecex <- xcex
}
else { ## axis="y"
chrspace <- par("cin")[2]
axislength <- par("pin")[2]
axismin <- par("usr")[3]
axisrange <- par("usr")[4] - par("usr")[3]
## in the y-axis case we only need to worry about 1 character since
## the y labels are perpendicular to the axis
chrspaces <- rep(chrspace, length(x))
thecex <- ycex
}
}
else if(grid) {
usrx <- convertX(unit(0:1, "npc"), "native", TRUE)
usry <- convertY(unit(0:1, "npc"), "native", TRUE)
## grid lattice natural coordinates for panel
axislength <- 1
if(axis=="x") {
chrspaces <- convertX(stringWidth(names(x)), "npc", TRUE)
axismin <- usrx[1]
axisrange <- usrx[2] - usrx[1]
chrspace <- convertX(stringWidth("0"), "npc", TRUE)
thecex <- xcex
}
else { ## axis="y"
chrspaces <- convertY(stringHeight(names(x)), "npc", TRUE)
axismin <- usry[1]
axisrange <- usry[2] - usry[1]
chrspace <- convertY(stringHeight("0"), "npc", TRUE)
thecex <- ycex
}
}
## Need to map response scale to device inches scale for traditional,
## npc for grid
p <- (axislength/axisrange) * (x - axismin)
n <- length(p)
## Stay in the p-space for conducting comparisons
## TRUE means to keep (no overlap)
## Keep min and max, and order statistics 2 and n - 1
if(n > 4) {
overlapchecks <- rep(TRUE, n)
lastOKtickindx <- 2
for(i in 3:(n-2)) {
if( ( 1.0*(p[i] - p[lastOKtickindx]) -
(thecex*(chrspaces[i] + chrspaces[lastOKtickindx])/2) ) > 0) {
overlapchecks[i] <- TRUE
lastOKtickindx <- i
}
else {
overlapchecks[i] <- FALSE
}
}
## Make sure the order n - 1 statistic is included
if(overlapchecks[n-1]==FALSE) {
overlapchecks[n-2] <- FALSE
overlapchecks[n-1] <- TRUE
}
## Purposely trim the minimum and maximum
overlapchecks[c(1, n)] <- FALSE
keep <- p[overlapchecks]
## Go back to the x space
out <- x[is.element(p, keep)]
}
else out <- x
## Ensure that zero/one no-difference tickmark is always included
if(difference==TRUE || percent==TRUE) {
if(!any(out==0)) {
out <- c(0, out)
if(ratio==TRUE && percent==FALSE) {
names(out)[1] <- "1"
}
else names(out)[1] <- "0"
out <- sort(out)
}
}
if(logscale) {
out <- 10^(out)
}
## If no entries
if(all(is.na(out))) {
out <- mean(range(x.orig))
names(out) <- " "
}
return(out)
}
#####################
minmaxTicks <- function(x, theaxis="y", logscale=TRUE,
ratio=FALSE, percent=FALSE, difference=FALSE,
digits=NULL, offset=NULL, zeroscore=NULL,
maxonly=FALSE, grid=FALSE, cex=0.7,...) {
##
## PURPOSE: Add min and max to axis tick labels
##
if(!is.null(zeroscore) && !is.null(offset)) {
stop(cgMessage("The zeroscore and offset arguments",
"cannot be specified at the same time."))
}
axisrange <- range(x)
if(is.numeric(offset)) { x <- x - offset }
labelrange <- range(x)
if(is.null(digits)) digits <- 4
if(logscale) {
if(percent==TRUE && theaxis=="x") {
thelabels <- fmtRatioToPercent(axisrange)
}
else if(percent==TRUE && theaxis=="y") {
thelabels <- fmtRatioToPercent(axisrange)
}
else {
thelabels <- trimWhiteSpace(chopZeroes(fround(labelrange, digits)))
}
if(!is.null(zeroscore)) thelabels[1] <- "0"
}
else {
if(!is.null(offset) || !is.null(zeroscore)) {
stop(cgMessage("The offset and zeroscore arguments",
"must be set to NULL when",
"logscale=FALSE, as they are not compatible."))
}
if(difference==TRUE && (theaxis=="x" || theaxis=="y")) {
thelabels <- fmtDifference(axisrange, digits=digits)
}
else {
thelabels <- trimWhiteSpace(chopZeroes(fround(labelrange, digits)))
}
}
if(maxonly) {
thelabels <- c("", thelabels[2])
}
if(!grid) {
if(theaxis=="y") {
text(x=rep(par("usr")[1], 2), y=if(logscale) log10(axisrange) else axisrange,
labels=paste("", thelabels), col="blue", adj=0, cex=cex)
}
else if(theaxis=="x") {
text(y=rep(par("usr")[3], 2), x=if(logscale) log10(axisrange) else axisrange,
labels=paste(thelabels,"\n"), col="blue", adj=0.5, cex=cex, xpd=NA)
}
}
else { ## grid==TRUE
if(theaxis=="y") {
panel.text(x=rep(current.panel.limits()$xlim[1], 2),
y=if(logscale) log10(axisrange) else axisrange,
labels=paste("", thelabels), col="blue", adj=0, cex=cex)
}
else if(theaxis=="x") {
panel.text(y=rep(current.panel.limits()$ylim[1], 2),
x=if(logscale) log10(axisrange) else axisrange,
labels=paste(thelabels, "\n"), col="blue", adj=0.5, cex=cex)
}
}
invisible()
}
#####################
plotGrpNameTicks <- function(grpnames, settings, grid=FALSE) {
##
## PURPOSE: Add group labels to x-axis
##
padEnd <- function(x, length=2) {
paste(x, paste(rep(" ", length), collapse=""), sep="")
}
numberofgrps <- length(grpnames)
if(!grid) {
if(settings$srt==0) {
axis(1, at=1:numberofgrps, labels=grpnames,
cex=settings$cex, adj=settings$adj)
}
else if(settings$srt!=0) {
axis(1, at=1:numberofgrps, labels=FALSE)
text(1:numberofgrps, par("usr")[3], srt = settings$srt,
adj = settings$adj, cex=settings$cex,
labels = padEnd(grpnames), ## to prevent axis overlap
xpd = TRUE)
}
}
else { ## grid, lattice
panel.axis(side="bottom",
at=1:length(grpnames),
labels=if(settings$srt!=0) {
padEnd(grpnames)
}
else {
grpnames
},
tck=0, text.cex=settings$cex,
rot=settings$srt,
outside=TRUE)
}
invisible()
}
#####################
boxplotStamp <- function() {
##
## PURPOSE: Add explanatory message to boxplot rendering
##
u <- par("usr")
if(FALSE) {
text(u[2] - 0.15 *(u[2] - u[1]),
u[4] + 0.15 * (u[4] - u[3]),
labels=paste("Open circles \"o\" ",
"inside boxes are medians;",
"Plus signs \"+\" are means", sep="\n"),
cex=0.6, adj=0, xpd=TRUE)
}
mtext(text=paste("Open circles \"o\" ",
"inside boxes are medians;",
"Plus signs \"+\" are means", sep=" "),
side=3,
cex=0.6, col="gray", line=0.1, adj=0.5,
xpd=NA)
invisible()
}
#####################
errorBarGraphStamp <- function(mcadjust=FALSE, alphapercent, grid=FALSE) {
##
## PURPOSE: Add explanatory message to errorBarGraph rendering
##
firstline <-
if(!grid) {
paste("Non-overlapping error bars indicate statistically\n",
"significant difference(s)",
sep=" ")
}
else {
paste("Non-overlapping error bars within a panel indicate statistically\n",
"significant difference(s) ",
sep=" ")
}
secondline <- paste("at", alphapercent,
"% significance level",
if(mcadjust) {
"(Multiplicity adjusted)"
} else "", sep=" ")
themessage <- paste(firstline,
secondline,
collapse="")
if(!grid) {
mtext(text=themessage, side=3, line=0.75, col="red", cex=0.7,
xpd=NA, outer=FALSE, adj=0.5, font=2)
}
else {
seekViewport(trellis.vpname("main"))
grid.text(themessage,
y=unit(-1, "lines"),
gp=gpar(col="red", cex=0.7, fontface="bold"),
just=c("centre","centre"))
}
invisible()
}
#####################
comparisonsGraphStamp <- function(mcadjust=FALSE, alphapercent, grid=FALSE, desc="") {
##
## PURPOSE: Add explanatory message to comparisonsGraph rendering
##
themessage <- paste(
paste(" Error bars that do not cross the zero line indicate statistically\n",
"significant difference(s) at ",
sep=" "),
alphapercent,
"% significance level",
if(mcadjust) {
" (Multiplicity adjusted)"
} else "",
sep="")
if(!grid) {
mtext(text=themessage, side=3, line=0.75, col="red", cex=0.7,
xpd=NA, outer=FALSE, adj=0, font=2)
if(desc!="") {
mtext(text=desc, side=3,
cex=0.6, col="gray", line=0.1, adj=0.5,
xpd=NA)
}
}
else {
seekViewport(trellis.vpname("main"))
grid.text(themessage,
y=unit(-1, "lines"),
gp=gpar(col="red", cex=0.7, fontface="bold"),
just=c("centre","centre"))
if(desc!="") {
grid.text(desc,
x=unit(0.5, "npc"),
y=unit(-2.75, "lines"),
gp=gpar(col="gray", cex=0.6))
}
}
invisible()
}
#####################
errorBarGraphApproximateStamp <- function(grid=FALSE, msg=NULL) {
##
## PURPOSE: Add warning message about unequal sample sizes
## and/or standard error estimates
##
if(is.null(msg)) msg <- "Error bars are approximate."
if(!grid) {
if(msg!="") {
mtext(text=msg, side=3,
cex=0.6, col="gray", line=0.1, adj=0.5,
xpd=NA)
}
}
else {
seekViewport(trellis.vpname("main"))
if(msg!="") {
grid.text(msg,
x=unit(0.5, "npc"),
y=unit(-2.75, "lines"),
gp=gpar(col="gray", cex=0.6))
}
}
invisible()
}
#####################
trimWhiteSpace <- function(x) {
##
## PURPOSE: Laziness motivated wrapper to
## trim leading and trailing spaces in a character string
## NOTE: x can be a vector of strings, and empty contents
## will be converted to an empty string.
## EXAMPLE USAGE: trimWhiteSpace(c(" 37 ","abc ","jfk", 4576))
##
if(length(x)==0) return("")
y <- gsub("^\\s+","", x, perl=TRUE) ## leading
return(gsub("\\s+$","", y, perl=TRUE)) ## trailing
}
#####################
chopZeroes <- function(x) {
##
## PURPOSE: For tickmark labels, delete any superfluous trailing zeroes after a
## decimal point. In effect, convert the character item to look like an integer.
## EXAMPLE USAGE: chopZeroes(c("30","12.0","33.000","478.020","65.20","0.1205"))
##
if(length(x)==0) return("")
gsub("\\.0*0+$", "", as.character(x), perl=TRUE)
}
#####################
fmtRatioToPercent <- function(x, digitsbelow10=1) {
##
## PURPOSE: Map ratios to percent units for tick labels
##
y <- 100 * (x - 1)
outfmt <- trimWhiteSpace(ifelse(abs(y) < 9.5 & y!=0, formatC(round(y, digitsbelow10),
format="f", digits=digitsbelow10),
chopZeroes(round(y, 0))))
outfmt[as.numeric(outfmt)==0] <- "0"
return(outfmt)
}
#####################
fmtDifference <- function(x, digits) {
##
## PURPOSE: Prettify Differences for tick labels
##
trimWhiteSpace(ifelse(abs(x) < 10,
chopZeroes(formatC(round(x, digits=digits),
format="f", digits=digits)),
chopZeroes(round(x, 0))))
}
#####################
fmtRatio <- function(x) {
##
## PURPOSE: Prettify Ratios for tick labels
##
trimWhiteSpace(ifelse(abs(x) < 0.1, formatC(round(x, 3), format="f", digits=3),
chopZeroes(round(x, 2))))
}
#####################
fmtPercent <- function(x, decimalthreshold=10) {
out <- fround(x, 0)
for(i in seq(along=out)) {
if(out[i]=="-100") {
out[i] <- "< -99.9"
}
## else if(out[i]=="0") {
else if(abs(x[i]) < decimalthreshold) {
out[i] <- fround(x[i], 1)
}
}
return(out)
}
#####################
fmtPvalue <- function(x) {
out <- fround(x, 3)
for(i in seq(along=x)) {
if(x[i] < 0.0005) {
out[i] <- "<0.001"
}
}
return(out)
}
#####################
cgDevice <- function(cgtheme=TRUE, new=FALSE, ...) {
##
## PURPOSE: Ensure trellis device with limited color scheme
if(cgtheme) {
thetheme <- list(background = list(col = "white"),
strip.shingle = list(col = "white"),
strip.background = list(col = "white"))
}
else {
thetheme <- trellis.par.get()
}
if(.Device=="null device") {
new <- TRUE
}
trellis.device(theme=thetheme, new=new, ...)
invisible()
}
#####################
contrastMatrix <- function(listoflevels,
type = c("allgroupstocontrol",
"pairwisereflect",
"pairwise",
"successivetimes")
) {
##
## PURPOSE: Utility function to construct comparisons specifications
## NOTE: The contrast matrix generation code is borrowed and adapted from the
## contrMat function of the R 'multcomp' package
##
if(mode(listoflevels) != "list") {
listoflevels <- list(listoflevels)
}
## start of R multcomp::contrMat code adaptation
type <- match.arg(type) ## ensure default is 'allgroupstocontrol'
## Note that "grpnames" is just a generic term for
## the different levels of a factor
genAtOneLevel <- function(grpnames, type=type) {
CM <- c()
rnames <- c()
k <- length(grpnames)
kindx <- 1:k
switch(type,
allgroupstocontrol = {
for (i in 2:k) CM <- rbind(CM, as.numeric(kindx == i) -
as.numeric(kindx == 1))
rnames <- paste(grpnames[2:k], " vs. ", grpnames[1], sep = "")
},
pairwisereflect = {
for (i in 1:k) {
for (j in (1:k)[-i]) {
CM <- rbind(CM, as.numeric(kindx == j) - as.numeric(kindx ==
i))
rnames <- c(rnames, paste(grpnames[j], " vs. ", grpnames[i],
sep = ""))
}
}
},
pairwise = {
for (i in 1:(k - 1)) {
for (j in (i + 1):k) {
CM <- rbind(CM, as.numeric(kindx == j) - as.numeric(kindx ==
i))
rnames <- c(rnames, paste(grpnames[j], " vs. ", grpnames[i],
sep = ""))
}
}
},
successivetimes = {
for(i in 2:k) {
CM <- rbind(CM, as.numeric(kindx == i) -
as.numeric(kindx == (i-1)))
rnames <- c(rnames, paste(grpnames[i], " vs. ", grpnames[i-1], sep = ""))
}
}
)
dimnames(CM) <- list(rnames, grpnames)
CM
}
## end of R multcomp::contrMat code adaptation
listLs <- lapply(listoflevels, genAtOneLevel, type=type)
if(length(listLs)==1) {
Lcomparisonmat <- listLs[[1]]
}
else {
Lcomparisonmat <- blockDiag(listLs)
}
return(Lcomparisonmat)
}
#####################
blockDiag <- function(listmat) {
##
## PURPOSE: Create a block diagonal matrix from list of matrices
## NOTE: adapted from Ben Bolker contribution on R-help
## http://www.r-project.org/nocvs/mail/r-help/2002/1533.html
##
nr <- sapply(listmat,nrow)
nc <- sapply(listmat,ncol)
cumnc <- cumsum(nc)
NR <- sum(nr)
NC <- sum(nc)
rowfun <- function(m,zbefore,zafter) {
cbind(matrix(0,ncol=zbefore,nrow=nrow(m)),m,
matrix(0,ncol=zafter,nrow=nrow(m)))
}
ret <- rowfun(listmat[[1]],0,NC-ncol(listmat[[1]]))
for (i in 2:length(listmat)) {
ret <- rbind(ret,rowfun(listmat[[i]],cumnc[i-1],NC-cumnc[i]))
}
dimnames(ret)[[2]] <- as.vector(unlist(lapply(listmat, function(x)
dimnames(x)[[2]])))
attr(ret, "colsperblock") <- nc
attr(ret, "rowsperblock") <- nr
ret
}
#####################
rangeExtend <- function(x, pct=4) {
##
## PURPOSE: Utility function to extend axis ranges, as sometimes trellis
## or lattice/grid functions in R seem to come a bit short, and points or
## labels spill over the plot region boundaries
##
if(length(pct)==1) {
pct <- list(minside=4, maxside=4)
}
rangelength <- diff(range(x, na.rm=TRUE))
minx <- max(x, na.rm=TRUE) - (1+(pct$minside/100))*rangelength
maxx <- min(x, na.rm=TRUE) + (1+(pct$maxside/100))*rangelength
c(minx, maxx)
}
#####################
getNumDigits <- function(dfr, max=4) {
##
## PURPOSE: Utility function to determine the maximum number of digits
## used in values of an input set of data
## argument dfr can be a data frame, list, matrix, or atomic vector
##
thedata <- as.character(as.vector(unlist(dfr)))
maxdigits <- max(nchar(gsub("(\\d*)(\\.*)(\\d*)", "\\3", thedata,
perl=TRUE)))
if(maxdigits > 4) maxdigits <- max
return(maxdigits)
}
#####################
makeCensored <- function(x, boundary, direction=">") {
fmtx <- x
for(i in seq(along=x)) {
if(x[i] >= boundary) {
fmtx[i] <- paste(direction, x[i], sep="")
}
}
return(fmtx)
}
#####################
cgLineColors <- c("black", "blue", "green", "red", "orange", "brown", "yellow",
"darkblue", "darkgreen", "darkgrey",
"black", "blue", "green", "red", "orange", "darkgrey")
#####################
multcompInform <- function() {
cat(cgMessage("Some time may be needed as the",
"critical point from the multcomp::summary.glht",
"function call is calculated. Please wait...\n\n",
warning=TRUE))
flush.console()
invisible()
}
multcompDone <- function(model) {
cat(cgMessage(paste("...Done. Critical point from", model,"fit",
"is calculated.\n\n"),
warning=TRUE))
flush.console()
invisible()
}
#####################
isAllEqual <- function(x) {
## x is assumed to be an atomic vector
return(length(unique(x))==1)
}
#####################
makeEndptLabel <- function(endptname, endptunits) {
## set up for input into expression
e1c <- if(is.expression(endptname)) {
paste(endptname) }
else {
paste("plain('",endptname,"')", sep="")
}
e2c <-
if(is.expression(endptunits)) {
paste(endptunits) }
else {
paste("plain('",endptunits,"')", sep="")
}
noe2c <- (e2c=="plain('')")
expr <- parse(text=paste(
"paste(",
e1c,
if(!noe2c) " , ",
if(!noe2c) "' (',",
if(!noe2c) e2c,
if(!noe2c) " , ",
if(!noe2c) "')'",
")",
sep=" "))
return(expr)
}
#####################
catCharExpr <- function(char, expr, rev=FALSE) {
if(!rev) {
out <- parse(text=paste("paste(",
paste("plain('",char,"')", sep=""),
" , ",
paste(expr),
")",
sep=" "))
}
else {
out <- parse(text=paste("paste(",
paste(expr),
" , ",
paste("plain('",char,"')", sep=""),
")",
sep=" "))
}
return(out)
}
#####################
residualgrptrend.helper <- function(tukeyresids, fitteds, tukeyresids2=NULL,
status=NULL, grpf, desc="") {
##
## PURPOSE: Determine a trend for residuals (y) against
## fitted values (x), with likely option of rescaling the fitted values
## to integer values afterwards: 1 to number of groups.
##
trend.index <- order(fitteds)
thegrpf <- factorInSeq(grpf[trend.index])
tukeyresids <- tukeyresids[trend.index]
fitteds <- fitteds[trend.index]
if(!is.null(status)) {
## Assume proper censored data format (survival:::Surv())
## no checking
has.censored <- TRUE
tukeyresids2 <- tukeyresids2[trend.index]
status <- status[trend.index]
thesurvobject <- survival::Surv(time=tukeyresids, time2=tukeyresids2,
event=status, type="interval")
if(nlevels(thegrpf) < 7) {
## Fit "group means" nominal AFT-normal model
thefit <- try(survreg(thesurvobject ~ -1 + thegrpf,
dist="gaussian", maxiter=100))
}
else {
## Need different data format for VGAM:::vgam()
forvgam.dfr <- data.frame(x=fitteds, y=tukeyresids)
extra <- list(leftcensored = (status==2),
rightcensored = (status==0))
thefit <- try(VGAM::vgam(y ~ s(x, 4), VGAM::cennormal(zero=2),
data=forvgam.dfr, trace=FALSE,
extra=extra))
if(inherits(thefit, "try-error")) {
thefit <- try(VGAM::vgam(y ~ s(x, 3), VGAM::cennormal(zero=2),
data=forvgam.dfr, trace=FALSE,
extra=extra))
}
if(inherits(thefit, "try-error")) {
thefit <- try(VGAM::vgam(y ~ s(x, 2), VGAM::cennormal(zero=2),
data=forvgam.dfr, trace=FALSE,
extra=extra))
}
else if(inherits(thefit, "try-error")) {
## Fit "group means" nominal AFT-normal model if smoothing fails
thefit <- try(survreg(thesurvobject ~ -1 + fitteds,
dist="gaussian", maxiter=100))
}
}
trendtukeyresids <- if(class(thefit)[1]=="survreg") {
predict(thefit) } else if(class(thefit)[1]=="vgam") { predict(thefit)[, 1] }
}
else {
has.censored <- FALSE
trendtukeyresids <- lowess(fitteds, tukeyresids)$y
}
dfr <- data.frame(grpc=I(as.character(thegrpf)),
grpn=as.numeric(thegrpf),
tukeyresids=tukeyresids,
trendtukeyresids=trendtukeyresids,
type=rep(desc, length(trend.index)))
if(has.censored) {
dfr$status <- status
}
return(dfr)
}
#####################
fround <- function(x, digits=0) {
format(round(x, digits=digits), nsmall=digits)
}
######################
fround.charcens <- function(x, digits=0) {
## generalization when x has "<" and ">" representations
for(i in seq(along=x)) {
if(is.na(x[i]) || x[i]=="<NA>") next
else if(regexpr("<|>", x[i]) > 0) {
pieces <- strsplit(x[i], "<|>")
x[i] <- paste(substr(x[i], 1, 1), fround(as.numeric(pieces[[1]][2]), digits),
sep="")
}
else {
x[i] <- fround(as.numeric(x[i]), digits)
}
}
return(x)
}
#####################
chop.matrix <- function(x) {
n <- nrow(x)
return(x[1:(n-1), 1:(n-1)])
}
#####################
prepare <- function(type, ... ) {
##
## PURPOSE: Read-in the input data frame, and create
## an object with it
## and add metadata for cg methods.
## Essentially supplies alias to specific prepare* methods
type <- validArgMatch(type, c("onefactor","unpairedgroups",
"pairedgroups","paireddifference",
"repmes","repeatedmeasures"))
if(is.element(type, c("onefactor","unpairedgroups"))) {
return(switch(type,
"onefactor"={ prepareCGOneFactorData(...) },
"unpairedgroups"={ prepareCGOneFactorData(...) }
))
}
else if(is.element(type, c("paireddifference","pairedgroups"))) {
return(switch(type,
"paireddifference"={ prepareCGPairedDifferenceData(...) },
"pairedgroups"={ prepareCGPairedDifferenceData(...) }
))
}
invisible()
}
#####################
stripmiss <- function(x) x[!is.na(x)]
#####################
makeTickMarks <- function(ticklabels, tickmarks, percent=FALSE,
offset=NULL) {
validList(ticklabels,
names=c("mod","marks"),
argname="ticklabels")
validBoolean(percent)
mod <- eval(parsePartialName(names(ticklabels), "mod",
prefix="ticklabels$"))
mod <- validArgMatch(mod, choices=c("add", "replace"))
marks <- eval(parsePartialName(names(ticklabels), "marks",
prefix="ticklabels$"))
if(mod=="replace") {
tickmarks <- as.numeric(marks)
if(percent) { tickmarks <- pctToRatio(tickmarks) }
if(!is.null(offset)) tickmarks <- (tickmarks + offset)
names(tickmarks) <- marks
}
else if(mod=="add") {
addedtickmarks <- as.numeric(marks)
if(percent) { addedtickmarks <- pctToRatio(addedtickmarks) }
if(!is.null(offset)) addedtickmarks <- (addedtickmarks + offset)
names(addedtickmarks) <- marks
tickmarks <- sort(c(tickmarks, addedtickmarks))
}
return(tickmarks)
}
#####################
scaleVar <- function(x, endptscale, logbase=10, percent=FALSE) {
##
## PURPOSE: Simple function to help work with variables
## that may or may not be transformed
endptscale <- validArgMatch(endptscale, c("log","original"))
validBoolean(percent)
if(endptscale=="log") {
if(percent) { x <- pctToRatio(x) }
return(log(x, logbase))
}
else return(x)
}
#####################
makeContrastVec <- function(posplaces=NULL, negplaces=NULL, length,
posvalues=rep(1, length(posplaces)),
negvalues=rep(-1, length(negplaces))
) {
##
## PURPOSE: Simple function to create vectors in
## a contrast matrix for custom comparisons
## WARNING: Argument checking needs to be implemented!
vec <- numeric(length) ## initialization
vec[posplaces] <- posvalues
vec[negplaces] <- negvalues
return(vec)
}
#####################
orderPairedGroups <- function(x, refgrp) {
## PURPOSE: For the PairedDifference
## methods, order the 2 group names
## such that the refgrp is first and
## derive the difference between then
if(x[1]==refgrp) {
return(x)
}
else {
c(x[2], x[1])
}
}
###########################
qminmin <- function(y, x, q) {
## PURPOSE: compute quantiles from survival curve estimate
## to aid in boxplot graphing
## NOTE: adapted from minmin helper function within
## survmean helper function from the survival package.
## (The minmin function used to be in the print.survfit function
## in earlier survival package versions).
p <- 1 - q
keep <- (!is.na(y) & y <= p)
if (!any(keep))
NA
else {
x <- x[keep]
y <- y[keep]
if (y[1] == p && any(y < p))
(x[1] + x[min(which(y < p))])/2
else x[1]
}
}
######################
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.