Nothing
if(!exists("NROW", mode='function')) {
NROW <- function(x)
if (is.array(x) || is.data.frame(x)) nrow(x) else length(x)
}
if(!exists("NCOL", mode='function')) {
NCOL <- function(x)
if (is.array(x) && length(dim(x)) > 1 || is.data.frame(x)) ncol(x) else as.integer(1)
}
prn <- function(x, txt, file='')
{
calltext <- as.character(sys.call())[2]
if(file != '') sink(file, append=TRUE)
if(!missing(txt)) {
if(nchar(txt) + nchar(calltext) +3 > .Options$width)
calltext <- paste('\n\n ',calltext,sep='')
else
txt <- paste(txt, ' ', sep='')
cat('\n', txt, calltext, '\n\n', sep='')
}
else cat('\n',calltext,'\n\n',sep='')
print(x)
if(file != '') sink()
invisible()
}
formatSep <- function(x, digits, ...)
{
y <- character(length(x))
for(i in 1:length(x))
y[i] <- if(missing(digits)) format(x[i], ...)
else format(x[i],digits=digits, ...)
names(y) <- names(x)
y
}
nomiss <- function(x)
{
if(is.data.frame(x)) na.exclude(x)
else if(is.matrix(x))
x[!is.na(x %*% rep(1,ncol(x))),]
else x[!is.na(x)]
}
fillin <- function(v, p)
{
v.f <- ifelse(is.na(v),p,v)
if(length(p)==1)
label(v.f) <- paste(label(v),"with",sum(is.na(v)),
"NAs replaced with",format(p))
else
label(v.f) <- paste(label(v),"with",sum(is.na(v)),"NAs replaced")
v.f
}
spearman <- function(x, y)
{
x <- as.numeric(x)
y <- as.numeric(y) ## 17Jul97
notna <- !is.na(x+y) ##exclude NAs
if(sum(notna) < 3)
c(rho=NA)
else
c(rho=cor(rank(x[notna]), rank(y[notna])))
}
plotCorrPrecision <- function(rho=c(0,0.5), n=seq(10,400,length.out=100),
conf.int=0.95, offset=.025, ...)
{
## Thanks to Xin Wang for computations
curves <- vector('list', length(rho))
names(curves) <- paste('r',format(rho),sep='=')
zcrit <- qnorm(1-(1-conf.int)/2)
for(i in 1:length(rho)) {
r <- rho[i]
z <- .5*log((1+r)/(1-r))
lo <- z - zcrit/sqrt(n-3)
hi <- z + zcrit/sqrt(n-3)
rlo <- (exp(2*lo)-1)/(exp(2*lo)+1)
rhi <- (exp(2*hi)-1)/(exp(2*hi)+1)
precision <- pmax(rhi-r, r-rlo)
curves[[i]] <- list(N=n, Precision=precision)
}
labcurve(curves, pl=TRUE, xrestrict=quantile(n,c(.25,1)), offset=offset, ...)
invisible()
}
trap.rule <- function(x,y) sum(diff(x)*(y[-1]+y[-length(y)]))/2
uncbind <- function(x, prefix="", suffix="")
{
nn <- dimnames(x)[[2]]
warning("You are using uncbind. That was a really bad idea. If you had any variables in the global environment named ", paste(prefix, nn, suffix, sep="", collapse=", "), " they are now over writen.\n\nYou are now warned.", immediate. = TRUE, )
for(i in 1:ncol(x))
assign(paste(prefix,nn[i],suffix,sep=""), x[,i], pos=parent.env())
invisible()
}
## Function to pick off ordinates of a step-function at user-chosen abscissas
stepfun.eval <- function(x, y, xout, type=c("left","right"))
{
s <- !is.na(x+y)
type <- match.arg(type)
approx(x[s], y[s], xout=xout, method="constant", f=if(type=="left")0 else 1)$y
}
km.quick <- function(S, times, q)
{
sRequire('survival')
S <- S[!is.na(S),]
n <- nrow(S)
stratvar <- factor(rep(1,nrow(S)))
f <- survival::survfitKM(stratvar, S, se.fit=FALSE, conf.type='none')
tt <- c(0, f$time)
ss <- c(1, f$surv)
if(missing(times))
min(tt[ss <= q])
else
approx(tt, ss, xout=times, method='constant', f=0)$y
}
oPar <- function()
{
## Saves existing state of par() and makes changes suitable
## for restoring at the end of a high-level graphics functions
oldpar <- par()
oldpar$fin <- NULL
oldpar$new <- FALSE
invisible(oldpar)
}
setParNro <- function(pars)
{
## Sets non-read-only par parameters from the input list
i <- names(pars) %nin%
c('cin','cra','csi','cxy','din','xlog','ylog','gamma','page')
invisible(par(pars[i]))
}
mgp.axis.labels <- function(value, type=c('xy','x','y','x and y'))
{
type <- match.arg(type)
if(missing(value)) {
value <- .Options$mgp.axis.labels
pr <- par(c('mgp', 'las'))
mgp <- pr$mgp
if(! length(value)) value <- c(.7, .7)
return(switch(type,
xy = value,
x = c(mgp[1], value[1], mgp[3]),
y = c(mgp[1], value[2], mgp[3]),
'x and y' = list(x = c(mgp[1], value[1], mgp[3]),
y = c(mgp[1], value[2], mgp[3]))))
}
if(value[1]=='default')
value <- c(.7,.7)
##c(.6, if(par('las')==1) 1.3 else .6)
options(mgp.axis.labels=value, TEMPORARY=FALSE)
invisible()
}
mgp.axis <-
function(side, at=NULL, ...,
mgp=mgp.axis.labels(type=if(side==1 | side==3)'x' else 'y'),
axistitle=NULL, cex.axis=par('cex.axis'), cex.lab=par('cex.lab'))
{
## Version of axis() that uses appropriate mgp from mgp.axis.labels and
## gets around bug in axis(2, ...) that causes it to assume las=1
mfrow <- par('mfrow')
tcl <- max(par('tcl'), -0.4)
nr <- mfrow[1]; nc <- mfrow[2]
w <- list(side=side)
w <- c(w, list(...))
w$cex.axis <- cex.axis
if(length(at)) w$at <- at
if(side == 1 || side == 3) {
w$mgp <- mgp / nr
w$tcl <- tcl / nr
if(side==1 && length(axistitle))
title(xlab=axistitle, mgp = mgp / min(2.25, nr), cex.lab=cex.lab)
} else {
w$mgp <- mgp / nc
w$tcl <- tcl / nc
las <- par('las')
w$srt <- 90 * (las == 0)
w$adj <- if(las == 0) 0.5
else 1
if(side == 2 && length(axistitle))
title(ylab=axistitle, mgp=mgp / min(2.25, nc), cex.lab=cex.lab)
}
do.call('axis', w)
invisible()
}
trellis.strip.blank <- function()
{
s.b <- lattice::trellis.par.get("strip.background")
s.b$col <- 0
lattice::trellis.par.set("strip.background", s.b)
s.s <- lattice::trellis.par.get("strip.shingle")
s.s$col <- 0
lattice::trellis.par.set("strip.shingle", s.s)
invisible()
}
lm.fit.qr.bare <- function(x, y,
tolerance = NULL,
intercept=TRUE, xpxi=FALSE,
singzero=FALSE)
{
if(!length(tolerance)) tolerance <- 1e-7
if(intercept)
x <- cbind(Intercept=1, x)
else x <- as.matrix(x)
z <- lm.fit(x, y, tol=tolerance)
coef <- z$coefficients
if(singzero && any(isna <- is.na(coef))) coef[isna] <- 0.
res <- z$residuals
sse <- sum(res^2)
sst <- sum((y - mean(y))^2)
res <- list(coefficients = coef, residuals = res,
rsquared = 1 - sse / sst,
fitted.values = z$fitted.values)
if(xpxi) {
p <- 1L : z$rank
res$xpxi <- chol2inv(z$qr$qr[p, p, drop=FALSE])
}
res
}
all.is.numeric <- function(x, what=c('test','vector','nonnum'),
extras=c('.','NA'))
{
what <- match.arg(what)
x <- sub('[[:space:]]+$', '', x)
x <- sub('^[[:space:]]+', '', x)
xs <- x[x %nin% c('', extras)]
if(! length(xs) || all(is.na(x)))
return(switch(what, test = FALSE, vector=x, nonnum=x[0]))
isnon <- suppressWarnings(! is.na(xs) & is.na(as.numeric(xs)))
isnum <- ! any(isnon)
# suppressWarnings below handles extras present in x
switch(what,
test = isnum,
vector = if(isnum) suppressWarnings(as.numeric(x)) else x,
nonnum = xs[isnon])
}
Lag <- function(x, shift=1)
{
## Lags vector x shift observations, padding with NAs or blank strings
## preserving attributes of x
xLen <- length(x)
if(shift == 0) return(x)
# Create base vector use character to generate "" for mode "character"
# Coerce base vector to be type of x
ret <- as.vector(character(xLen), mode=storage.mode(x))
# set resp attributes equal to x attributes
attrib <- attributes(x)
if(length(attrib$label))
attrib$label <- paste(attrib$label, 'lagged', shift, 'observations')
if(abs(shift) < xLen)
{
if(shift > 0) ret[-(1:shift)] <- x[1:(xLen - shift)]
else ret[1:(xLen+shift)] <- x[(1-shift):xLen]
}
attributes(ret) <- attrib
return(ret)
}
xySortNoDupNoNA <- function(x, y)
{
if(is.list(x)) {
y <- x[[2]]; x <- x[[1]]
}
s <- !is.na(x + y)
if(any(s)) {
x <- x[s]; y <- y[s]
}
i <- order(x)
x <- x[i]
y <- y[i]
i <- !duplicated(x)
list(x=x[i], y=y[i])
}
outerText <-
function(string, y, cex=par('cex'), ...) {
usr <- par('usr'); plt <- par('plt')
pos <- usr[2] + (usr[2] - usr[1])/(plt[2] - plt[1]) * (1 - plt[2])
axis(2, at=y, labels=string, tick=FALSE, las=1,
pos=pos, cex.axis=cex, xpd=NA)
}
## if(missing(space)) space <- max(nchar(string))*.5
## mtext(string, side=side, las=1, at=y, adj=adj, cex=cex, line=space)
# This method does not survive shrinking the graphics window
# Right justifies (if adj=1) a vector of strings against the right margin
# (side=4) or against the y-axis (side=2)
#outerText <-
# function(string, y, side=4, cex=par('cex'), adj=1, ...) {
# if(side %nin% c(2,4)) stop('only works for side=2 or 4')
# x <- if(side==4) grconvertX(1, from='nfc', to='user') else
# par('usr')[1]
# text(x, y, paste(string,''), cex=cex, adj=adj, xpd=NA)
#}
## Old method [dropped because does not scale upon resizing device]
## Use text() to put test strings in left or right margins
## Temporarily sets par(xpd=NA) if using R
## For adj=1 side=4, setAside is a character string used to determine
## the space to set aside for all strings
## space is the number of extra characters to leave to the left of
## the string(s) (adj=0) or to the right (adj=1)
if(FALSE) outerText <- function(string, y, setAside=string[1], side=4, space=1,
adj=1, cex=par('cex'))
{
usr <- par('usr')
xpd <- par('xpd')
if(!is.na(xpd)) {
on.exit(par(xpd=xpd))
par(xpd=NA)
}
ie <- is.expression(string) ## 1sep02
if(ie)
adj <- 0 ## adj=1 not work well for expressions in R
if(side!=4)
stop('only side=4 implemented')
if(adj==0)
text(usr[2], y,
if(ie)
string
else
paste(space,string,sep=''),
adj=0)
else {
usr.space.needed <- strwidth(setAside, units='user', cex=cex)
text(usr[2]+0.5*strwidth(space, units='user', cex=cex)+usr.space.needed,
y, string, adj=1, cex=cex) # was usr[2]- 18jul02;added 0* 25jul02
## was 0*strwidth(space,...) 31jan03
}
invisible()
}
if(FALSE) {
expandUsrCoord <- function()
{
## Expands usr coordinates of current plot to entire figure region
## so that out of range plots may be plotted
pr <- par()
usr <- pr$usr
p <- pr$plt
invisible(pr)
}
}
## Author: Patrick Connolly <P.Connolly@hortresearch.co.nz>
## HortResearch
## Mt Albert
## Auckland, New Zealand
print.char.matrix <-
function (x, file = "",
col.name.align = "cen", col.txt.align = "right",
cell.align = "cen", hsep = "|", vsep = "-", csep = "+",
row.names = TRUE, col.names = FALSE,
append = FALSE, top.border = TRUE, left.border = TRUE, ...)
{
### To print a data frame or matrix to a text file or screen
### and having names line up with stacked cells
###
### First, add row names as first column (might be removed later)
ndimn <- names(dimnames(x)) ## FEH
rownames <- dimnames(x)[[1]]
x <- cbind(rownames, x)
names(dimnames(x)) <- ndimn ## FEH
cnam <- dimnames(x)[[2]] ## FEH
if(length(ndimn))
cnam[1] <- ndimn[1] ## FEH
##dimnames(x)[[1]] <- seq(nrow(x)) 25Mar02 for R FEH
dimnames(x) <- list(as.character(seq(nrow(x))), cnam)
names(dimnames(x)) <- ndimn ## 26Mar02 FEH
### Set up some padding functions:
###
pad.left <- function(z, pads)
{
## Pads spaces to left of text
padding <- paste(rep(" ", pads), collapse = "")
paste(padding, z, sep = "")
}
pad.mid <- function(z, pads)
{
## Centres text in available space
padding.right <- paste(rep(" ", pads%/%2), collapse = "")
padding.left <- paste(rep(" ", pads - pads%/%2), collapse = "")
paste(padding.left, z, padding.right, sep = "")
}
pad.right <- function(z, pads) {
## Pads spaces to right of text
padding <- paste(rep(" ", pads), collapse = "")
paste(z, padding, sep = "")
}
## (Padding happens on the opposite side to alignment)
pad.types <- c("left", "mid", "right")
names(pad.types) <- c("right", "cen", "left")
pad.name <- pad.types[col.name.align]
pad.txt <- pad.types[col.txt.align]
pad.cell <- pad.types[cell.align]
## Padding character columns
## Need columns with uniform number of characters
pad.char.col.right <- function(y)
{
## For aligning text to LHS of column
col.width <- nchar(y)
biggest <- max(col.width)
smallest <- min(col.width)
padding <- biggest - col.width
out <- NULL
for (i in seq(y))
out[i] <- pad.right(y[i], pads = padding[i])
out
}
pad.char.col.left <- function(y)
{
## For aligning text to RHS of column
col.width <- nchar(y)
biggest <- max(col.width)
smallest <- min(col.width)
padding <- biggest - col.width
out <- NULL
for (i in seq(y))
out[i] <- pad.left(y[i], pads = padding[i])
out
}
pad.char.col.mid <- function(y) {
## For aligning text to centre of column
col.width <- nchar(y)
biggest <- max(col.width)
smallest <- min(col.width)
padding <- biggest - col.width
out <- NULL
for (i in seq(y))
out[i] <- pad.mid(y[i], pads = padding[i])
out
}
## which functions to use this time.
pad.name.fn <- get(paste("pad.", pad.name, sep = ""))
pad.txt.fn <- get(paste("pad.char.col.", pad.txt, sep = ""))
pad.cell.fn <- get(paste("pad.", pad.cell, sep = ""))
## Remove troublesome factors
x <- as.data.frame(x)
fac.col <- names(x)[sapply(x, is.factor)]
for (i in fac.col)
x[, i] <- I(as.character(x[, i]))
## ARE ANY LINE BREAKS IN ANY COLUMNS?
break.list <- list()
for (i in seq(nrow(x))) {
x.i <- unlist(x[i, ])
rows.i <- sapply(strsplit(unlist(x[i, ]), "\n"), length)
rows.i[rows.i < 1] <- 1
break.list[[i]] <- rows.i
}
break.row <- sapply(break.list, function(x) any(x > 1))
names(break.row) <- seq(nrow(x))
xx <- x
if (any(break.row)) {
## add in extra row/s
xx <- NULL
reprow <- lapply(break.list, unique)
for (k in seq(nrow(x))) {
x.k <- unlist(x[k, ])
x.k[x.k == ""] <- " "
if (break.row[k]) {
l.k <- strsplit(x.k, "\n")
add.blanks <- max(break.list[[k]]) - break.list[[k]]
names(l.k) <- names(add.blanks) <- seq(length(l.k))
if (any(add.blanks > 0)) {
for (kk in names(add.blanks[add.blanks > 0]))
l.k[[kk]] <- c(l.k[[kk]], rep(" ", add.blanks[kk]))
}
l.k.df <- as.data.frame(l.k)
names(l.k.df) <- names(x)
xx <- rbind(xx, as.matrix(l.k.df))
}
else xx <- rbind(xx, x.k)
}
row.names(xx) <- paste(rep(row.names(x), sapply(reprow,
max)),
unlist(reprow), sep = ".")
## Make an index for the rows to be printed
rn <- row.names(xx)
rnb <- strsplit(rn, "\\.")
rpref <- as.numeric(factor(sapply(rnb, function(z) z[1])))
## was codes( ) 10oct03
}
else
rpref <- seq(nrow(x))
x <- as.data.frame(xx)
## Character columns need different treatment from numeric columns
char.cols <- sapply(x, is.character)
if (any(char.cols))
x[char.cols] <- sapply(x[char.cols], pad.txt.fn)
## Change numeric columns into character
if (any(!char.cols))
x[!char.cols] <- sapply(x[!char.cols], format)
## now all character columns each of which is uniform element width
##
## Lining up names with their columns
## Sometimes the names of columns are wider than the columns they name,
## sometimes vice versa.
names.width <- nchar(names(x))
if (!col.names)
names.width <- rep(0, length(names.width))
cell.width <- sapply(x, function(y) max(nchar(as.character(y))))
## (the width of the characters in the cells as distinct
## from their names)
name.pads <- cell.width - names.width
cell.pads <- -name.pads
name.pads[name.pads < 0] <- 0
cell.pads[cell.pads < 0] <- 0
pad.names <- name.pads > 0
pad.cells <- cell.pads > 0
## Pad out the column names if necessary:
if (any(pad.names)) {
stretch.names <- names(x)[pad.names]
for (i in stretch.names) {
names(x)[names(x) == i] <- pad.name.fn(i, name.pads[i])
}
}
## likewise for the cells and columns
if (any(pad.cells)) {
stretch.cells <- names(x)[pad.cells]
for (j in stretch.cells) x[, j] <- pad.cell.fn(x[, j],
cell.pads[j])
}
## Remove row names if not required
if (!row.names)
x <- x[-1]
## Put the column names on top of matrix
if (col.names)
mat2 <- rbind(names(x), as.matrix(x))
else
mat2 <- as.matrix(x)
mat.names.width <- nchar(mat2[1, ])
## character string to separate rows
space.h <- ""
for (k in seq(along=mat.names.width)) { ## added along= FEH 26Mar02
space.h <- c(space.h, rep(vsep, mat.names.width[k]), csep)
}
line.sep <- paste(c(ifelse(left.border, csep, ""), space.h),
collapse = "")
if (col.names)
rpref <- c(0, rpref, 0)
else
rpref <- c(rpref, 0)
## print to screen or file
if(top.border && line.sep !='') {
write(line.sep, file = file, append = append)
append <- TRUE
}
for (i in 1:nrow(mat2)) {
if (left.border)
write(paste(paste(c("", mat2[i, ]), collapse = hsep),
hsep, sep = ""), file = file, append = append)
else
write(paste(paste(mat2[i, ], collapse = hsep), hsep,
sep = ""), file = file, append = append)
append <- TRUE
## print separator if row prefix is not same as next one
if (rpref[i] != rpref[i + 1] && line.sep != '')
write(line.sep, file = file, append = TRUE)
}
}
unPaste <- function(str, sep='/')
{
w <- strsplit(str, sep)
w <- matrix(unlist(w), ncol=length(str))
nr <- nrow(w)
ans <- vector('list', nr)
for(j in 1:nr)
ans[[j]] <- w[j,]
ans
}
get2rowHeads <- function(str) {
w <- strsplit(str, '\n')
## strsplit returns character(0) when element="" 23may03
list(sapply(w, function(x)if(length(x)) x[[1]] else ''),
sapply(w, function(x)if(length(x) > 1)x[[2]] else ''))
}
## Note: can't say f[vector of names] <- list(...) to update args
## In R you have to put ALL arguments in list(...) so sometimes we set
## unneeded ones to NULL. Ignore this assignment in S<
## Two lists of functions, one for primitives for S+ or R (either Trellis
## or low-level), one for R grid
## Note: rect is only defined in R, not S+
ordGridFun <- function(grid)
{
if(!grid)
list(lines = function(...) lines(...),
points = function(..., size=NULL)
{
if(length(size))
warning('size not implemented yet')
points(...)
},
text = function(...) text(...),
segments = function(...) segments(...),
arrows = function(..., open, size) arrows(..., length=size*.8),
rect = function(...) rect(...),
polygon = function(x, y=NULL, ..., type=c('l','s'))
{
type <- match.arg(type)
if(!length(y))
{
y <- x$y
x <- x$x
}
j <- !is.na(x+y)
x <- x[j]
y <- y[j]
if(type=='s') polygon(makeSteps(x, y), ..., border=NA)
else polygon(x, y, ..., border=NA)
},
abline = function(...) abline(...),
unit = function(x, units='native')
{
if(units!='native')
stop('units="native" is only units implemented outside of grid')
x
},
axis = function(...) axis(...))
else {
sRequire('lattice')
list(lines = function(x, y, ...)
{
if(is.list(x)) {
y <- x[[2]]; x <- x[[1]]
}
lattice::llines(if(is.unit(x))
convertX(x, 'native', valueOnly=TRUE)
else x,
if(is.unit(y))
convertY(y, 'native', valueOnly=TRUE)
else y,
...)
},
points = function(x, y, ...)
{
if(is.list(x)) {
y <- x[[2]]; x <- x[[1]]
}
lattice::lpoints(if(is.unit(x))
convertX(x, 'native', valueOnly=TRUE)
else x,
if(is.unit(y))
convertY(y, 'native', valueOnly=TRUE)
else y,
...)
},
text = function(x, y, ...)
{
if(is.list(x)) {
y <- x[[2]]; x <- x[[1]]
}
lattice::ltext(if(is.unit(x))
convertX(x, 'native', valueOnly=TRUE)
else x,
if(is.unit(y))
convertY(y, 'native', valueOnly=TRUE)
else y,
...)
},
segments = function(x0, y0, x1, y1, ...)
{
grid.segments(x0, y0, x1, y1, default.units='native',
gp=gpar(...))
},
arrows = function(...) lattice::larrows(...),
rect = function(xleft, ybottom, xright, ytop, density, angle,
border, xpd, ...)
{
grid.rect(xleft, ybottom, width=xright-xleft,
height=ytop-ybottom, just='left',
default.units='native', gp=gpar(...))
},
polygon = function(x, y=NULL, col=par('col'), type=c('l','s'), ...)
{
type <- match.arg(type)
if(!length(y))
{
y <- x$y
x <- x$x
}
j <- !is.na(x+y)
x <- x[j]
y <- y[j]
if(type=='s') grid.polygon(makeSteps(x, y),
default.units='native',
gp=gpar(fill=col, col='transparent', ...))
else grid.polygon(x, y, default.units='native',
gp=gpar(fill=col,col='transparent',...))
},
abline=function(...) lattice::panel.abline(...),
unit = function(x, units='native', ...) unit(x, units=units, ...),
axis = function(side=1, at=NULL, labels, ticks=TRUE,
distn, line, pos, outer, ...)
{
if(!length(at))stop('not implemented for at= unspecified')
if(side > 2) stop('not implemented for side=3 or 4')
## ticks=ticks removed from grid.?axis FEH 30Aug09
if(side==1) grid.xaxis(at=at, label=labels, gp=gpar(...))
if(side==2) grid.yaxis(at=at, label=labels, gp=gpar(...))
})
}
}
parGrid <- function(grid=FALSE)
{
pr <- par()
cin <- pr$cin
cex <- pr$cex
lwd <- pr$lwd
if(grid) {
## cvp <- current.viewport()
## usr <- c(cvp$xscale, cvp$yscale)
usr <- c(convertX(unit(0:1, "npc"), "native", valueOnly=TRUE),
convertY(unit(0:1, "npc"), "native", valueOnly=TRUE))
pin <-
c(convertWidth(unit(1, "npc"), "inches", valueOnly=TRUE),
convertHeight(unit(1, "npc"), "inches", valueOnly=TRUE))
uin <-
c(convertWidth(unit(1, "native"), "inches", valueOnly=TRUE),
convertHeight(unit(1, "native"), "inches", valueOnly=TRUE))
}
else {
usr <- pr$usr
pin <- pr$pin
uin <- c(pin[1]/(usr[2]-usr[1]), pin[2]/(usr[4]-usr[3]))
## 22Mar01 - R does not have par(uin)
}
list(usr=usr, pin=pin, uin=uin, cin=cin, cex=cex, lwd=lwd)
}
## Replaces R's xinch, yinch, extending them to grid
## These convert inches to data units
xInch <- function(x=1, warn.log=!grid, grid=FALSE)
{
if (warn.log && par("xlog"))
warning("x log scale: xInch() is nonsense")
pr <- parGrid(grid)
x * diff(pr$usr[1:2])/pr$pin[1]
}
yInch <- function (y = 1, warn.log=!grid, grid=FALSE)
{
if (warn.log && par("ylog"))
warning("y log scale: yInch is nonsense")
pr <- parGrid(grid)
y * diff(pr$usr[3:4])/pr$pin[2]
}
na.include <- function(obj) {
if(inherits(obj,'data.frame'))
for(i in seq(along=obj))
obj[[i]] <- na.include(obj[[i]])
else {
if(length(levels(obj)) && any(is.na(obj)))
obj <- factor(obj,exclude=NULL)
}
obj
}
if(FALSE) {
whichClosest <- function(x, w)
{
## x: vector of reference values
## w: vector of values to find closest matches in x
## Returns: subscripts in x corresponding to w
i <- order(x)
x <- x[i]
n <- length(x)
br <- c(-1e30, x[-n]+diff(x)/2,1e30)
m <- length(w)
i[.C("bincode", as.double(w), m, as.double(br),
length(br), code = integer(m), right = TRUE,
include = FALSE, NAOK = TRUE, DUP = FALSE,
PACKAGE = "base")$code]
}
NULL
}
## Just as good, ties shuffled to end
## function(x, w) round(approx(x,1:length(x),xout=w,rule=2,ties='ordered')$y)
## Remove ties= for S-Plus. Note: does not work when 2nd arg to
## approx is not uniformly spaced
## NO! ties='ordered' bombs in x not ordered
## Try
## approx(c(1,3,5,2,4,2,4),1:7,xout=c(1,3,5,2,4,2,4),rule=2,ties=function(x)x[1])
## NO: only works in general if both x and y are already ordered
## The following runs the same speed as the previous S version (in R anyway)
whichClosest <- function(x, w)
{
## x: vector of reference values
## w: vector of values for which to lookup closest matches in x
## Returns: subscripts in x corresponding to w
## Assumes no NAs in x or w
.Fortran(F_wclosest,as.double(w),as.double(x),
length(w),length(x),
j=integer(length(w)))$j
}
whichClosePW <- function(x, w, f=0.2) {
lx <- length(x)
lw <- length(w)
.Fortran(F_wclosepw,as.double(w),as.double(x),
as.double(runif(lw)),as.double(f),
lw, lx, double(lx), j=integer(lw))$j
}
whichClosek <- function(x, w, k) {
## x: vector of reference values
## w: vector of values for which to lookup close matches in x
## Returns: subscripts in x corresponding to w
## Assumes no NAs in x or w
## First jitters x so there are no ties
## Finds the k closest matches and takes a single random pick of these k
y <- diff(sort(x))
mindif <- if(all(y == 0)) 1 else min(y[y > 0])
x <- x + runif(length(x), -mindif/100, mindif/100)
z <- abs(outer(w, x, "-"))
s <- apply(z, 1, function(u) order(u)[1:k])
if(k == 1) return(s)
apply(s, 2, function(u) sample(u, 1))
}
if(FALSE) {
sampWtdDist <- function(x, w)
{
## x: vector of reference values
## w: vector of values to find closest matches in x
## Returns: subscripts in x corresponding to w
## 25% slower but simpler method:
## z <- abs(outer(w, x, "-"))
## s <- apply(z, 1, max)
## z <- (1 - sweep(z, 1, s, FUN='/')^3)^3
## sums <- apply(z, 1, sum)
## z <- sweep(z, 1, sums, FUN='/')
lx <- length(x)
lw <- length(w)
z <- matrix(abs( rep( x , lw ) - rep( w, each = lx ) ),
nrow=lw, ncol=lx, byrow=TRUE) ## Thanks: Chuck Berry
## s <- pmax( abs( w - min(x) ), abs( w - max(x) ) ) # to use max dist
s <- rowSums(z)/lx/3 # use 1/3 mean dist for each row
tricube <- function(u) (1 - pmin(u,1)^3)^3
## z <- (1 - (z/rep(s,length.out=lx*lw))^3)^3 # Thanks: Tim Hesterberg
z <- tricube(z/s) # Thanks: Tim Hesterberg
sums <- rowSums(z)
z <- z/sums
as.vector(rMultinom(z, 1))
}
NULL
}
approxExtrap <- function(x, y, xout, method='linear', n=50, rule=2,
f=0, ties='ordered', na.rm=FALSE)
{
## Linear interpolation using approx, with linear extrapolation
## beyond the data
if(is.list(x)) {
y <- x[[2]]; x <- x[[1]]
}
## remove duplicates and order so can do linear extrapolation
if(na.rm) {
d <- ! is.na(x + y)
x <- x[d]; y <- y[d]
}
x <- as.numeric(x) # handles dates etc.
y <- as.numeric(y)
d <- ! duplicated(x)
x <- x[d]
y <- y[d]
d <- order(x)
x <- x[d]
y <- y[d]
w <- approx(x, y, xout=xout, method=method, n=n,
rule=2, f=f, ties=ties)$y
r <- range(x)
d <- xout < r[1]
if(any(is.na(d)))
stop('NAs not allowed in xout')
if(any(d))
w[d] <- (y[2]-y[1])/(x[2]-x[1])*(xout[d]-x[1])+y[1]
d <- xout > r[2]
n <- length(y)
if(any(d))
w[d] <- (y[n]-y[n-1])/(x[n]-x[n-1])*(xout[d]-x[n-1])+y[n-1]
list(x=xout, y=w)
}
inverseFunction <- function(x, y) {
d <- diff(y)
xd <- x[-1]
dl <- c(NA, d[-length(d)])
ic <- which(d>=0 & dl<0 | d>0 & dl<=0 | d<=0 & dl>0 | d<0 & dl>=0)
nt <- length(ic)
k <- nt + 1
if(k==1) {
h <- function(y, xx, yy, turns, what, coef)
approx(yy, xx, xout=y, rule=2)$y
formals(h) <- list(y=numeric(0), xx=x, yy=y, turns=numeric(0),
what=character(0), coef=numeric(0))
return(h)
}
turns <- x[ic]
turnse <- c(-Inf, turns, Inf)
xrange <- yrange <- matrix(NA, nrow=k, ncol=2)
for(j in 1:k) {
l <- which(x >= turnse[j] & x <= turnse[j+1])
xrange[j,] <- x[l[c(1,length(l))]]
yrange[j,] <- y[l[c(1,length(l))]]
}
for(j in 1:length(ic)) {
l <- (ic[j]-1):(ic[j]+1)
turns[j] <- approxExtrap(d[l], xd[l], xout=0, na.rm=TRUE)$y
}
hh <- function(y, xx, yy, turns, xrange, yrange, what, coef) {
what <- match.arg(what)
## Find number of monotonic intervals containing a given y value
ylo <- pmin(yrange[,1],yrange[,2])
yhi <- pmax(yrange[,1],yrange[,2])
n <- outer(y, ylo, function(a,b) a >= b) &
outer(y, yhi, function(a,b) a <= b)
## Columns of n indicate whether or not y interval applies
ni <- nrow(yrange)
fi <- matrix(NA, nrow=length(y), ncol=ni)
turnse <- c(-Inf, turns, Inf)
for(i in 1:ni) {
w <- n[,i]
if(any(w)) {
l <- xx >= turnse[i] & xx <= turnse[i+1]
fi[w,i] <- approx(yy[l], xx[l], xout=y[w])$y
}
}
noint <- !apply(n, 1, any)
if(any(noint)) {
## Determine if y is closer to yy at extreme left or extreme right
## of an interval
m <- length(yy)
yl <- as.vector(yrange); xl <- as.vector(xrange)
fi[noint,1] <- xl[whichClosest(yl, y[noint])]
}
if(what=='sample')
apply(fi, 1, function(x) {
z <- x[!is.na(x)]
if(length(z)==1) z else if(length(z)==0) NA else sample(z, size=1)
}) else fi
}
formals(hh) <- list(y=numeric(0), xx=x, yy=y, turns=turns,
xrange=xrange, yrange=yrange,
what=c('all', 'sample'), coef=numeric(0))
## coef is there for compatibility with areg use
hh
}
Names2names <- function(x)
{
if(is.list(x)) {
}
else {
n <- names(attributes(x))
if(any(n=='.Names'))
names(attributes(x)) <- ifelse(n=='.Names','names',n)
}
x
}
##xedit <- function(file, header, title, delete.file=FALSE) {
## In R, use e.g. options(pager=xedit); page(x,'p')
## sys(paste('xedit -title "', title, '" ', file, ' &',
## sep=''))
## invisible()
##}
if(FALSE) {
gless <- function(x, ...)
{
## Usage: gless(x) - uses print method for x, puts in window with
## gless using name of x as file name prefixed by ~, leaves window open
nam <- substring(deparse(substitute(x)), 1, 40)
file <- paste('/tmp/',nam,sep='~') #tempfile('Rpage.')
sink(file)
## cat(nam,'\n' )
## if(length(attr(x,'label')) && !inherits(x,'labelled'))
## cat(attr(x,'label'),'\n')
## cat('\n')
print(x, ...)
sink()
sys(paste('gless --geometry=600x400 "',file,'" &',sep=''))
## gless does not have a title option
invisible()
}
NULL
}
xless <-
function(x, ..., title=substring(deparse(substitute(x)),1,40))
{
## Usage: xless(x) - uses print method for x, puts in persistent window with
## xless using name of x as title (unless title= is specified)
## If running under MacOS, use the system open command instead of xless
file <- tempfile()
sink(file)
print(x, ...)
sink()
cmd <- if(Sys.info()['sysname'] == 'Darwin') paste('open -a TextEdit', file) else
paste('xless -title "',title,'" -geometry "90x40" "',
file,'" &',sep='')
system(cmd)
invisible()
}
pasteFit <- function(x, sep=',', width=.Options$width)
{
## pastes as many elements of character vector x as will fit in a line
## of width 'width', starting new lines when needed
## result is the lines of pasted text
m <- nchar(x)
out <- character(0)
cur <- ''
n <- 0
for(i in 1:length(x)) {
if(cur=='' | (m[i] + nchar(cur) <= width))
cur <- paste(cur, x[i],
sep=if(cur=='')''
else sep)
else {
out <- c(out, cur)
cur <- x[i]
}
}
if(cur != '') out <- c(out, cur)
out
}
## Determine if variable is a date, time, or date/time variable in R.
## The following 2 functions are used by describe.vector
## timeUsed assumes is date/time combination variable and has no NAs
testDateTime <- function(x, what=c('either','both','timeVaries'))
{
what <- match.arg(what)
cl <- class(x)
if(!length(cl))
return(FALSE)
dc <- c('Date', 'POSIXt','POSIXct','dates','times','chron')
dtc <- c('POSIXt','POSIXct','chron')
switch(what,
either = any(cl %in% dc),
both = any(cl %in% dtc),
timeVaries = {
if('chron' %in% cl || 'Date' %in% cl) {
## chron or S+ timeDate
y <- as.numeric(x)
length(unique(round(y - floor(y),13))) > 1
}
else length(unique(format(x,'%H%M%S'))) > 1
})
}
## Format date/time variable from either R or S+
## x = a numeric summary of the original variable (e.g., mean)
## at = attributes of original variable
formatDateTime <- function(x, at, roundDay=FALSE)
{
cl <- at$class
w <- if(any(cl %in% c('chron','dates','times'))){
attributes(x) <- at
fmt <- at$format
if(roundDay) {
if (!requireNamespace("chron", quietly = TRUE))
stop("'roundDay = TRUE' requires the 'chron' package.")
if(length(fmt)==2 && is.character(fmt))
format(chron::dates(x), fmt[1])
else
format(chron::dates(x))
}
else x
} else {
attributes(x) <- at
if(roundDay && 'Date' %nin% at$class)
as.POSIXct(round(x, 'days'))
else x
}
format(w)
}
## Try to guess whether a factor or character variable is a date, and
## allow for partial dates of the form YYYY and mm/YYYY, the former
## only used if at least one observation has a month in it. If a minority
## fraction of observations fracnn or less is not convertable to a date,
## set those observations to NA
## Allows for mixing of common date forms across observations
##
## Example:
##
## x <- convertPdate(c(rep('2019-03-04',7), '2018', '03/2018', 'junk', '3/11/17','2017-01-01','2017-01-01', NA, ''))
## x
## describe(x)
convertPdate <- function(x, fracnn=0.3, considerNA=NULL) {
xo <- x
if(is.factor(x)) x <- as.character(x)
if(! is.character(x)) return(xo)
x <- trimws(x)
if(all(is.na(x)) || all(x =='')) return(xo)
ymd <- grepl('^[0-9]{4}-[0-9]{1,2}-[0-9]{1,2}$', x)
mdy <- grepl('^[0-9]{1,2}/[0-9]{1,2}/[0-9]{4}$', x)
y <- grepl('^[0-9]{4}$', x)
my <- grepl('^[0-9]{1,2}/[0-9]{4}$', x)
m <- sum(x %nin% c('', considerNA) & ! is.na(x)) * (1 - fracnn)
dny <- ymd | mdy | my # date other than just YYYY
d <- dny | (y & any(dny))
## The variable is a date if at least m values are dates
if(sum(d) < m) return(xo)
special <- obs <- NULL
ndnm <- ! d & ! is.na(x) # not date and not missing
if(any(ndnm)) {
special <- setdiff(x[ndnm], c('', ' '))
obs <- x %in% special
x[ndnm] <- NA # values such as text comments
}
x <- ifelse(y, paste0(x, '-07-03'),
ifelse(my, gsub('^([0-9]{1,2})/([0-9]{4})$', '\\2-\\1-15', x),
ifelse(ymd, x,
ifelse(mdy, gsub('^([0-9]{1,2})/([0-9]{1,2})/([0-9]{4})$',
'\\3-\\1-\\2', x), NA))))
x <- as.Date(x)
if(length(special))
attr(x, 'special.miss') <- list(codes=special, obs=obs)
if(any(y | my)) attr(x, 'imputed') <- which(y | my)
x
}
getHdata <-
function(file, what=c('data','contents','description','all'),
where='https://hbiostat.org/data/repo') {
what <- match.arg(what)
fn <- as.character(substitute(file))
localrepo <- .Options$localHfiles
localrepo <- length(localrepo) && is.logical(localrepo) && localrepo
if(localrepo) where <- '~/web/data/repo'
ads <- readLines(paste0(where, '/Rcontents.txt'))
a <- unlist(strsplit(ads,'.sav|.rda'))
if(missing(file)) return(a)
wds <- paste(substitute(file), c('rda','sav'), sep='.')
if(!any(wds %in% ads))
stop(paste(paste(wds, collapse=','),
'are not on the web site.\nAvailable datasets:\n',
paste(a, collapse=' ')))
wds <- wds[wds %in% ads]
if(what %in% c('contents','all')) {
w <- paste(if(fn=='nhgh')'' else 'C',fn,'.html',sep='')
browseURL(paste(where, w, sep='/'))
}
if(what %in% c('description','all')) {
ades <- scan(paste(where,'Dcontents.txt',sep='/'),list(''),
quiet=TRUE)[[1]]
i <- grep(paste(fn,'\\.',sep=''),ades)
if(!length(i))
warning(paste('No description file available for',fn))
else {
w <- ades[i[1]]
browseURL(paste(where, w, sep='/'))
}
}
if(what %nin% c('data','all'))
return(invisible())
f <- paste(where, wds, sep='/')
if(length(f) > 1)
warning(paste('More than one file matched; using the first:',
paste(f, collapse=', ')))
if(localrepo) return(invisible(load(f, .GlobalEnv)))
tf <- tempfile()
download.file(f, tf, mode='wb', quiet=TRUE)
load(tf, .GlobalEnv)
invisible()
}
hdquantile <- function(x, probs=seq(0, 1, 0.25), se=FALSE,
na.rm=FALSE, names=TRUE, weights=FALSE)
{
if(na.rm) {
na <- is.na(x)
if(any(na))
x <- x[!na]
}
x <- sort(x, na.last=TRUE)
n <- length(x)
if(n < 2)
return(rep(NA, length(probs)))
m <- n + 1
ps <- probs[probs > 0 & probs < 1]
qs <- 1 - ps
a <- outer((0:n)/n, ps,
function(x,p,m) pbeta(x, p*m, (1-p)*m), m=m)
w <- a[-1,,drop=FALSE] - a[-m,,drop=FALSE]
r <- drop(x %*% w)
rp <- range(probs)
pp <- ps
if(rp[1]==0) {
r <- c(x[1], r); pp <- c(0,pp)
}
if(rp[2]==1) {
r <- c(r, x[n]); pp <- c(pp,1)
}
r <- r[match(pp, probs)]
if(names) names(r) <- format(probs)
if(weights)
attr(r,'weights') <- structure(w, dimnames=list(NULL,format(ps)))
if(!se)
return(r)
if(n < 3)
stop('must have n >= 3 to get standard errors')
l <- n - 1
a <- outer((0:l)/l, ps,
function(x,p,m) pbeta(x, p*m, (1-p)*m), m=m)
w <- a[-1,,drop=FALSE] - a[-n,,drop=FALSE]
storage.mode(x) <- 'double'
storage.mode(w) <- 'double'
nq <- length(ps)
## Get all n leave-out-one quantile estimates
S <- matrix(.Fortran(F_jacklins, x, w, as.integer(n), as.integer(nq),
res=double(n*nq))$res, ncol=nq)
se <- l * sqrt(diag(var(S))/n)
if(rp[1]==0)
se <- c(NA, se)
if(rp[2]==1)
se <- c(se, NA)
se <- se[match(pp,probs)]
if(names)
names(se) <- names(r)
attr(r, 'se') <- se
r
}
sepUnitsTrans <- function(x,
conversion=c(day=1, month=365.25/12, year=365.25, week=7),
round=FALSE, digits=0)
{
if(!any(is.present(x)))
return(x)
target <- names(conversion[conversion==1])
if(!length(target))
stop('must specify a target unit with conversion factor=1')
lab <- attr(x,'label')
x <- ifelse(is.present(x),casefold(as.character(x)),'')
for(w in names(conversion)) {
i <- grep(w, x)
if(length(i)) x[i] <-
as.character(as.numeric(gsub(paste(w,'s*',sep=''), '', x[i]))*
conversion[w])
}
i <- grep('[a-z]', x)
if(any(i))
warning(paste('variable contains units of measurement not in',
paste(names(conversion), collapse=','),':',
paste(unique(x[i]),collapse=' ')))
x <- as.numeric(x)
if(round)
x <- round(x, digits)
units(x) <- target
if(length(lab))
label(x) <- lab
x
}
makeNames <- function(names, unique=FALSE, allow=NULL)
{
## Runs make.names with exceptions in vector allow
## By default, R 1.9 make.names is overridden to convert _ to . as
## with S-Plus and previous versions of R. Specify allow='_' otherwise.
n <- make.names(names, unique)
if(!length(allow))
n <- gsub('_', '.', n)
n
}
Load <- function(object)
{
nam <- deparse(substitute(object))
path <- .Options$LoadPath
if(length(path))
path <- paste(path,'/',sep='')
file <- paste(path, nam, '.rda', sep='')
load(file, .GlobalEnv)
}
Save <- function(object, name=deparse(substitute(object)), compress=TRUE)
{
path <- .Options$LoadPath
if(length(path))
path <- paste(path, '/', sep='')
.FileName <- paste(path, name, '.rda', sep='')
assign(name, object)
if(is.logical(compress) && compress) compress <- 'gzip'
eval(parse(text=paste('save(', name, ', file="',
.FileName, '", compress="', compress, '")', sep='')))
}
getZip <- function(url, password=NULL) {
## Allows downloading and reading a .zip file containing one file
## File may be password protected. Password will be requested unless given.
## Password is 'foo'
## url may also be a local file
## Note: to make password-protected zip file z.zip, do zip -e z myfile
if(grepl("^https?://", tolower(url))) {
f <- tempfile()
download.file(url, f)
} else f <- url
cmd <- if(length(password))
paste('unzip -p -P', password) else 'unzip -p'
pipe(paste(cmd, f))
}
getLatestSource <- function(x=NULL, package='Hmisc',
recent=NULL, avail=FALSE) {
urlf <- paste0('https://hbiostat.org/R/', package, '/dir.txt')
fs <- scan(urlf, what=list('', ''), sep=' ', quiet=TRUE)
dates <- fs[[1]]
files <- fs[[2]]
url <- if(length(recent))
paste0('https://github.com/harrelfe/', package, '/commits/master/R')
else
paste0('https://github.com/harrelfe/', package, '/tree/master/R/')
if(avail) return(data.frame(file=files, date=as.Date(dates)))
if(length(recent)) x <- files[1:recent]
if(length(x)==1 && x=='all') x <- files
for(fun in x) {
i <- which(files==fun)
if(!length(i)) stop(paste('no file ', fun,' in ',package, sep=''))
cat('Fetching', fun, dates[i],'\n')
url <- paste0('https://raw.githubusercontent.com/harrelfe/', package,
'/master/R/', fun)
source(url)
}
}
clowess <- function(x, y=NULL, iter=3, ...) {
## to get around bug in lowess with occasional wild values with iter>0
r <- range(if(length(y)) y else x$y)
f <- lowess(x, y, iter=iter, ...)
if(iter != 0 && any(f$y < r[1] | f$y > r[2]))
f <- lowess(x, y, iter=0)
f
}
prselect <- function(x, start=NULL, stop=NULL, i=0, j=0, pr=TRUE)
{
f <- function(pattern, x)
{
y <- grep(pattern, x)
if(length(y) > 1) y <- y[1]
y
}
lx <- length(x)
k <- if(length(start)) f(start, x) else 1
if(length(k))
{
k <- k + i
m <- if(length(stop))
{
w <- f(stop, x[k:lx])
if(length(w)) w + k - 1 + j else -1
}
else lx
if(m > 0) x <- if(k==1) (if(m==lx) '...' else c('...', x[-(k:m)]))
else
{
if(m==lx) c(x[-(k:m)], '...')
else c(x[1:(k-1)], '...', x[(m+1):lx])
}
}
else # no start specified; keep lines after stop
{
m <- f(stop, x)
if(length(m) > 0)
{
m <- if(length(m)) m + j - 1 else lx
x <- if(m==lx) '...' else c('...', x[-(1:m)])
}
}
if(pr) cat(x, sep='\n')
invisible(x)
}
## The following is taken from survival:::plot.survfit internal dostep function
## Remove code to remove duplicates in y
makeSteps <- function(x, y)
{
if (is.na(x[1] + y[1]))
{
x <- x[-1]
y <- y[-1]
}
n <- length(x)
if (n > 2)
{
xrep <- rep(x, c(1, rep(2, n - 1)))
yrep <- rep(y, c(rep(2, n - 1), 1))
list(x = xrep, y = yrep)
}
else if (n == 1)
list(x = x, y = y)
else list(x = x[c(1, 2, 2)], y = y[c(1, 1, 2)])
}
latexBuild <- function(..., insert=NULL, sep='') {
w <- list(...)
l <- length(w)
if(l %% 2 != 0) stop('# arguments must be multiple of 2')
k <- l / 2
j <- 1
txt <- op <- character(0)
for(i in 1 : k) {
a <- w[[j]]
if(length(a)) {
txt <- c(txt, a)
if(w[[j + 1]] != '') op <- c(op, w[[j + 1]])
}
j <- j + 2
}
txt <- paste(txt, collapse=sep)
w <- character(0)
close <- if(length(op)) {
for(y in rev(op)) {
if(length(insert))
for(ins in insert)
if(length(ins) &&
ins[[1]] == y && ins[[2]] == 'before')
w <- c(w, '\n', ins[[3]])
w <- c(w,
if(y == '(') ')'
else if(y == '{') '}'
else if(y == '[') ']'
else sprintf('\\end{%s}', y))
if(length(insert))
for(ins in insert)
if(length(ins) &&
ins[[1]] == y && ins[[2]] == 'after')
w <- c(w, '\n', ins[[3]])
}
paste(w, collapse=sep)
}
structure(txt, close=close)
}
getRs <- function(file=NULL,
guser='harrelfe', grepo='rscripts',
gdir='raw/master', dir=NULL,
browse=c('local', 'browser'), cats=FALSE,
put=c('source', 'rstudio')) {
browse <- match.arg(browse)
put <- match.arg(put)
localrepo <- .Options$localHfiles
localrepo <- length(localrepo) && is.logical(localrepo) && localrepo
if(localrepo) where <- '~/R/rscripts'
else {
where <- paste('https://github.com', guser, grepo, gdir, sep='/')
if(length(dir)) where <- paste(where, dir, sep='/')
}
trim <- function(x) sub('^[[:space:]]+','',sub('[[:space:]]+$','', x))
pc <- function(s) {
wr <- function(x) {
n <- length(x)
z <- character(n)
for(i in 1 : n) z[i] <- paste(strwrap(x[i], width=15), collapse='\n')
z
}
s <- with(s, cbind(Major = wr(Major),
Minor = wr(Minor),
File = wr(File),
Type = wr(Type),
Description = wr(Description)))
print.char.matrix(s, col.names=TRUE)
}
read.table.HTTPS <- function(url) {
res <- tryCatch(read.table(url,
sep='|', quote='', header=TRUE, as.is=TRUE),
error=function(e) e)
if(inherits(res, "simpleError")) {
if(res$message == "https:// URLs are not supported") {
res$message <- paste(res$message, "Try installing R version >= 3.2.0", sep="\n\n")
}
stop(res)
}
res
}
download.file.HTTPS <- function(url, file, method='libcurl',
quiet=TRUE, extra='--no-check-certificate') {
res <- tryCatch(download.file(url, file, method, quiet=quiet, extra=extra),
error=function(e) e)
if(inherits(res, "simpleError")) {
if(res$message == "download.file(method = \"libcurl\") is not supported on this platform") {
warning(paste(res$message, "Try installing R version >= 3.2.0", "Attempting method=\"wget\"", sep="\n\n"))
return(download.file.HTTPS(url, file, method='wget'))
}
if(res$message == "https:// URLs are not supported") {
res$message <- paste(res$message, "Try installing R version >= 3.2.0", sep="\n\n")
}
stop(res)
}
invisible(res)
}
if(! length(file)) {
s <- read.table.HTTPS(paste(where, 'contents.md', sep='/'))
s <- s[-1,]
names(s) <- c('Major', 'Minor', 'File', 'Type', 'Description')
sd <- s; n <- nrow(s) # sd = s with dittoed items duplicated
for(x in c('Major', 'Minor')) {
u <- v <- gsub('\\*\\*', '', trim(s[[x]]))
for(i in 2 : n) if(u[i] == '"') u[i] <- u[i - 1]
v <- gsub('"', '', v)
s[[x]] <- v; sd[[x]] <- u
}
s$File <- trim(gsub('\\[(.*)\\].*', '\\1', s$File))
d <- trim(gsub('\\[.*\\]\\(.*\\)', '', s$Description))
s$Description <- gsub('\\[report\\].*', '', d)
if(is.logical(cats)) {
if(cats) {
## List all major and minor categories
maj <- sort(unique(sd$Major))
min <- setdiff(sort(unique(sd$Minor)), '')
cat('\nMajor categories:\n', maj,
'\nMinor categories:\n', min, '', sep='\n')
return(invisible(list(Major=maj, Minor=min)))
}
} else { ## list all scripts whose "first hit" major category contains cats
i <- grepl(tolower(cats), tolower(sd$Major))
if(! any(i)) cat('No scripts with', cats, 'in major category\n')
else pc(s[i, ])
return(invisible(s[i, ]))
}
if(browse == 'local') pc(s)
else
browseURL(if(localrepo) '~/R/rscripts/contents.md'
else
'https://github.com/harrelfe/rscripts/blob/master/contents.md')
return(invisible(s))
}
if(put == 'source')
return(invisible(source(paste(where, file, sep='/'))))
if(localrepo) file.copy(paste(where, file, sel='/'), file)
else download.file.HTTPS(paste(where, file, sep='/'), file)
if(requireNamespace('rstudioapi', quietly=TRUE) &&
rstudioapi::isAvailable()) rstudioapi::navigateToFile(file)
else file.edit(file)
invisible()
}
knitrSet <-
function(basename = NULL,
w=if(! bd) 4,
h=if(! bd) 3,
wo=NULL, ho=NULL,
fig.path = if(length(basename)) basename else '',
fig.align = if(! bd) 'center',
fig.show = 'hold',
fig.pos = if(! bd) 'htbp',
fig.lp = if(! bd) paste('fig', basename, sep=':'),
dev = switch(lang,
latex='pdf', markdown='png',
blogdown=NULL, quarto=NULL),
tidy=FALSE, error=FALSE,
messages=c('messages.txt', 'console'),
width=61, decinline=5, size=NULL, cache=FALSE,
echo=TRUE, results='markup', capfile=NULL,
lang=c('latex','markdown','blogdown','quarto')) {
if(! requireNamespace('knitr')) stop('knitr package not available')
messages <- match.arg(messages)
lang <- match.arg(lang)
options(knitrSet.lang = lang)
bd <- lang %in% c('blogdown', 'quarto')
## Specify e.g. dev=c('pdf','png') or dev=c('pdf','postscript')
## to produce two graphics files for each plot
## But: dev='CairoPNG' is preferred for png
if(length(basename)) basename <- paste(basename, '-', sep='')
## Default width fills Sweavel boxes when font size is \small and svmono.cls
## is in effect (use 65 without svmono)
if(lang == 'latex') knitr::render_listings()
if(messages != 'console') {
unlink(messages) # Start fresh with each run
hook_log = function(x, options) cat(x, file=messages, append=TRUE)
knitr::knit_hooks$set(warning = hook_log, message = hook_log)
}
else
knitr::opts_chunk$set(message=FALSE, warning=FALSE)
if(length(size)) knitr::opts_chunk$set(size = size)
## For htmlcap see http://stackoverflow.com/questions/15010732
## Causes collisions in html and plotly output; Original (no better)
## enclosed in <p class="caption"> ... </p>
# if(lang == 'markdown')
# knitr::knit_hooks$set(htmlcap = function(before, options, envir) {
# if(! before) options$htmlcap
# htmltools::HTML(paste0('<br><div style="font-size: 75%;">',
# options$htmlcap, "</div><br>"))
# })
if(length(decinline)) {
rnd <- function(x, dec) if(!is.numeric(x)) x else round(x, dec)
formals(rnd) <- list(x=NULL, dec=decinline)
knitr::knit_hooks$set(inline = rnd)
}
knitr::knit_hooks$set(par=function(before, options, envir)
if(before && options$fig.show != 'none') {
p <- c('bty','mfrow','ps','bot','top','left','rt','lwd',
'mgp','las','tcl','axes','xpd')
pars <- knitr::opts_current$get(p)
pars <- pars[! is.na(names(pars))]
## knitr 1.6 started returning NULLs for unspecified pars
i <- sapply(pars, function(x) length(x) > 0)
.spar. <-
function(mar=if(!axes)
c(2.25+bot-.45*multi,2*(las==1)+2+left,.5+top+.25*multi,
.5+rt) else
c(3.25+bot-.45*multi,2*(las==1)+3.5+left,.5+top+.25*multi,
.5+rt),
lwd = if(multi)1 else 1.75,
mgp = if(!axes) mgp=c(.75, .1, 0) else
if(multi) c(1.5, .365, 0) else c(2.4-.4, 0.475, 0),
tcl = if(multi)-0.25 else -0.4, xpd=FALSE, las=1,
bot=0, left=0, top=0, rt=0, ps=if(multi) 14 else 12,
mfrow=NULL, axes=TRUE, cex.lab=1.15, cex.axis=1,
...) {
multi <- length(mfrow) > 0
par(mar=mar, lwd=lwd, mgp=mgp, tcl=tcl, ps=ps, xpd=xpd,
cex.lab=cex.lab, cex.axis=cex.axis, las=las, ...)
if(multi) par(mfrow=mfrow)
}
if(any(i)) do.call(.spar., pars[i]) else .spar.()
})
knitr::opts_knit$set(width=width)
if(length(capfile)) {
options(FigCapFile=capfile)
cf <- function(before, options, envir) {
if(before) return()
lang <- getOption('knitrSet.lang')
label <- knitr::opts_current$get('label')
prefx <- if(lang == 'quarto') '' else options$fig.lp
figname <- paste0(prefx, label)
## Quarto uses a chunk figure label convention fig-...
## and figures are referenced by @fig-...
figref <- if(grepl('^fig-', figname))
paste0('@', figname) else paste0('\\@ref(', figname, ')')
cap <- options$fig.cap
scap <- options$fig.scap
if(length(cap) && is.call(cap)) cap <- eval(cap)
if(length(scap) && is.call(scap)) scap <- eval(scap)
if( ! length(scap) || scap == '') scap <- cap
if(length(scap) && scap != '')
cat(label, figref, paste0('"', scap, '"\n'), sep=',',
append=TRUE, file=getOption('FigCapFile'))
}
knitr::knit_hooks$set(capfileFun=cf)
}
## May want to see https://stackoverflow.com/questions/37116632/r-markdown-html-number-figures
## aliases=c(h='fig.height', w='fig.width', cap='fig.cap', scap='fig.scap'))
## eval.after = c('fig.cap','fig.scap'),
## error=error) #, keep.source=keep.source (TRUE))
## See if need to remove dev=dev from below because of plotly graphics
w <- list(fig.path=fig.path, fig.align=fig.align,
fig.width=w, fig.height=h,
out.width=wo,out.height=ho,
fig.show=fig.show, fig.lp=fig.lp, fig.pos=fig.pos,
dev=dev, par=TRUE, capfileFun=length(capfile) > 0,
tidy=tidy, cache=cache,
echo=echo, error=error, comment='', results=results)
if(bd) w$fig.path <- NULL
w <- w[sapply(w, function(x) length(x) > 0)]
## knitr doesn't like null fig.align etc.
do.call(knitr::opts_chunk$set, w)
if(lang %in% c('markdown', 'blogdown'))
knitr::knit_hooks$set(uncover=markupSpecs$html$uncover)
hook_chunk = knitr::knit_hooks$get('chunk')
## centering will not allow too-wide figures to go into left margin
if(lang == 'latex') knitr::knit_hooks$set(chunk = function(x, options) {
res = hook_chunk(x, options)
if (options$fig.align != 'center') return(res)
gsub('\\{\\\\centering (\\\\includegraphics.+)\n\n\\}',
'\\\\centerline{\\1}', res)
})
knitr::set_alias(w = 'fig.width', h = 'fig.height',
wo = 'out.width', ho = 'out.height',
cap = 'fig.cap', scap ='fig.scap')
}
## see http://yihui.name/knitr/options#package_options
## Use caption package options to control caption font size
grType <- function() {
if(! length(find.package('plotly', quiet=TRUE))) return('base')
if(length(g <- .Options$grType) && g == 'plotly') 'plotly' else 'base'
}
prType <- function() {
g <- .Options$prType
if(! length(g)) 'plain' else g
}
htmlSpecialType <- function() {
if(length(g <- .Options$htmlSpecialType) && g == '&') '&' else 'unicode'
}
## Save a plotly graphic with name foo.png where foo is the name of the
## current chunk
## http://stackoverflow.com/questions/33959635/exporting-png-files-from-plotly-in-r
plotlySave <- function(x, ...) {
if (!requireNamespace("plotly"))
stop("This function requires the 'plotly' package.")
chunkname <- knitr::opts_current$get("label")
path <- knitr::opts_chunk$get('fig.path')
if(is.list(x) & ! inherits(x, 'plotly_hash')) {
for(w in names(x)) {
file <- paste0(path, chunkname, '-', w, '.png')
plotly::plotly_IMAGE(x[[w]], format='png', out_file=file, ...)
}
}
else {
file <- paste0(path, chunkname, '.png')
plotly::plotly_IMAGE(x, format='png', out_file=file, ...)
}
invisible()
}
## Miscellaneous functions helpful for plotly specifications
plotlyParm = list(
## Needed height in pixels for a plotly dot chart given the number of
## rows in the chart
heightDotchart = function(rows, per=25, low=200, high=800)
min(high, max(low, per * rows)),
## Given a vector of row labels that appear to the left on a dot chart,
## compute the needed chart height taking label line breaks into account
## Since plotly devotes the same vertical space to each category,
## just need to find the maximum number of breaks present
heightDotchartb = function(x, per=40,
low=c(200, 200, 250, 300, 375)[min(nx, 5)],
high=1700) {
x <- if(is.factor(x)) levels(x) else sort(as.character(x))
nx <- length(x)
m <- sapply(strsplit(x, '<br>'), length)
# If no two categories in a row are at the max # lines,
# reduce max by 1
mx <- max(m)
lm <- length(m)
mlag <- if(lm == 1) 0 else c(0, m[1:(lm - 1)])
if(! any(m == mx & mlag == mx)) mx <- mx - 1
z <- 1 + (if(mx > 1) 0.5 * (mx - 1) else 0)
min(high, max(low, per * length(x) * z))
},
## Colors for unordered categories
colUnorder = function(n=5, col=colorspace::rainbow_hcl) {
if(! is.function(col)) rep(col, length.out=n)
else col(n)
},
## Colors for ordered levels
colOrdered = function(n=5, col=viridis::viridis) {
if(! is.function(col)) rep(col, length.out=n)
else col(n)
},
## Margin to leave enough room for long labels on left or right as
## in dotcharts
lrmargin = function(x, wmax=190, mult=7) {
if(is.character(x)) x <- max(nchar(x))
min(wmax, max(70, x * mult))
}
)
## Function written by Dirk Eddelbuettel:
tobase64image <- function (file, Rd = FALSE, alt = "image") {
input <- normalizePath(file, mustWork = TRUE)
buf <- readBin(input, raw(), file.info(input)$size)
base64 <- base64enc::base64encode(buf)
sprintf("%s<img src=\"data:image/png;base64,\n%s\" alt=\"%s\" />%s",
if (Rd)
"\\out{"
else "", base64, alt, if (Rd)
"}"
else "")
}
plotp <- function(data, ...) UseMethod("plotp")
keepHattrib <- function(obj) {
g <- function(x) {
a <- attributes(x)
i <- intersect(names(a), c('label', 'units'))
if(length(i)) a[i]
}
if(! is.list(obj)) list(.single.variable.=g(obj)) else sapply(obj, g)
}
restoreHattrib <- function(obj, attribs) {
nam <- names(obj)
for(n in names(attribs)) {
a <- attribs[[n]]
if(length(a)) {
sv <- n == '.single.variable.'
if(sv || n %in% nam) {
x <- if(sv) obj else obj[[n]]
if(length(a$label)) label(x) <- a$label
if(length(a$units)) units(x) <- a$units
if(sv) return(x)
obj[[n]] <- x
}
}
}
obj
}
if(FALSE) {
Hglossary <-
list(Gmd=list('Gini\'s mean difference', 'a measure of dispersion defined as the mean absolute difference over all possible pairs of different observations. It is more robust than the standard deviation.', 'https://www.researchgate.net/publication/5182211_Gini\'s_Mean_Difference_A_Superior_Measure_of_Variability_for_Non-Normal_Distributions'),
Info=list('Information index', 'a measure of the information content in a numeric variable relative to the information in a continuous numeric variable with no ties. The lowest value of Info occurs in a very imbalanced binary variable. Info comes from the approximate formula for the variance of a log odds ratio for a proportional odds model/Wilcoxon test, due to Whitehead (1993). Info is the ratio of the variance if there no ties in the data to the variance for the frequency distribution of observed values.', 'https://hbiostat.org/bib/r2.html')
)
rHglossary <- function(x, html=TRUE, collapse=TRUE) {
nams <- names(Hglossary)
i <- which(tolower(nams) == tolower(x))
if(! length(i))
stop(paste(x, 'is not defined in Hglossary'))
w <- Hglossary[[i]]
sname <- nams[i]
lname <- w[[1]]
def <- w[[2]]
href <- w[[3]]
if(html) {
lname <- paste0('<a href="', href[1], '">', lname, '</a>')
if(length(href) > 1)
def <- paste0(def, ' <a href="', href[2], '">More information</a>')
}
if(collapse) paste0('<details><summary>', sname, '</summary>',
lname, ': ', def, '</details>')
else paste0(sname, ': ', lname, ', ', def)
}
}
## Function to render HTML
## Converts argument to one character string with \n delimiters
## If knitr is currently running, runs this string through
## knitr::asis_output
## Otherwise, makes it browsable HTML using htmltools so that
## an RStudio Viewer or a new browser window will display the result
## See https://github.com/quarto-dev/quarto-cli/discussions/4248 which
## explains that you meed to enclose the text to keep from fooling
## Pandoc's reader
rendHTML <- function(x, html=TRUE) {
x <- paste(x, collapse='\n')
raw <- getOption('rawmarkup', FALSE)
if(raw) {
cat(x, '\n')
return(invisible())
}
if(length(getOption('knitr.in.progress'))) {
if(html)
return(htmltools::knit_print.html(x)) # includes htmlPreserve
## if(html) x <- paste0('```{=html}\n\n', x, '\n```\n')
return(knitr::asis_output(x))
}
if(! html) { # Convert from RMarkdown to html
tf <- tempfile(fileext='.Rmd')
o <- tempfile(fileext='.html')
cat('---\ntitle: ""\npagetitle: x\noutput: html_document\n---\n',
x, '\n', sep='', file=tf)
rmarkdown::render(tf, output_file=o, quiet=TRUE)
x <- readLines(o)
}
## The following has prettier output for model fits than the kableExtra method
print(htmltools::browsable(htmltools::HTML(x)))
}
## See ~/r/rmarkdown/html/render.qmd
## kableExtra print method did not render regression fit output
## as nicely as print(htmltools::browsable(...)) method above
## See kableExtra:::print.kableExtra
# class(x) <- 'kableExtra'
# dep <- list(rmarkdown::html_dependency_jquery(),
# rmarkdown::html_dependency_bootstrap(theme = "cosmo"),
# kableExtra::html_dependency_kePrint(),
# kableExtra::html_dependency_lightable())
# ht <- htmltools::browsable(htmltools::HTML(as.character(x),
# "<script type=\"text/x-mathjax-config\">MathJax.Hub.Config({tex2jax: {inlineMath: [[\"$\",\"$\"]]}})</script><script async src=\"https://mathjax.rstudio.com/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML\"></script>"))
# htmltools::htmlDependencies(ht) <- dep
# htmltools::html_print(ht)
sRequire <- function(package) {
if(! requireNamespace(package, quietly=TRUE))
stop(paste('package', package, 'is required but not installed'))
invisible()
}
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.