termdeps.R
# gives formula terms as a list of symbols (interactions as sub-list)
# [note that it does not expand formulas]
# termlist(terms(~a * b+ c, simplify = TRUE))
## termlist(~a+b+a:b) --> list(a, b, list(a, b))
termlist <- function(x) {
is.plus <- function(x) is.call(x) && x[[1L]] == "+"
## parses interaction expression into list: a:b:c --> list(a,b,c)
intr <- function(x) {
# is it an expression for interaction? (e.g. a:b:c)
is.intr <- function(x) is.call(x) && x[[1L]] == ":"
if(is.intr(x)) {
res <- list()
repeat {
res <- c(x[[3L]], res)
x <- x[[2L]]
if(!is.intr(x)) break
}
list(c(x, res))
} else x
}
if(x[[1L]] == "~") x <- x[[length(x)]]
res <- list()
while(is.plus(x)) {
res <- c(intr(x[[3L]]), res)
x <- x[[2L]]
}
res <- c(intr(x), res)
res
}
# calculates all lower order term names:
# expandintr(1:3) --> c("1", "2", "1:2", "3", "1:3", "2:3", "1:2:3")
expandintr <- function(x) {
asstr <- function(x) asChar(x, backtick = TRUE)
if(!is.language(x)) {
a <- sapply(x, asstr)
k <- length(a)
vapply(seq.int(2L^k - 1L), function(y) paste0(a[as.logical(intToBits(y)[1L:k])],
collapse = ":"), "")
} else asstr(x)
}
# given a formula, 'term dependency matrix', i.e. dependency of higer
# order terms on other lower order terms
termdepmat <- function(f) {
trm <- terms(f, simplify = TRUE)
tl <- termlist(trm)
v <- attr(trm, "term.labels")
n <- length(v)
mat <- matrix(FALSE, n, n, dimnames = list(v, v))
for(i in seq.int(n)) mat[match(expandintr(tl[[i]]), v), i] <- TRUE
mat
}
# alternative to 'termdepmat', gives matrix dimension names as numbers
# so a,b,a:b become 1,2,1:2
termdepmat2 <- function(f) {
filist <- formula2idx(f, asCall = FALSE)
n <- length(filist)
v <- vapply(filist, paste0, "", collapse = ":")
mat <- matrix(FALSE, n, n, dimnames = list(v, v))
for(i in seq.int(n)) mat[match(expandintr(filist[[i]]), v), i] <- TRUE
mat
}
## combines term-dependency-matrices
#termdepmat_list <- function(fl)
# termdepmat_combine(lapply(fl, termdepmat))
termdepmat_combine <- function(x) {
dm <- sum(vapply(x, nrow, 1L))
mat <- matrix(FALSE, dm, dm)
j <- 1L
for(i in seq_along(x)) {
n <- nrow(x[[i]])
k <- seq.int(j, length.out = n)
mat[k, k] <- x[[i]]
j <- j + n
}
dn <- unlist(lapply(x, rownames))
dimnames(mat) <- list(dn, dn)
mat
}
# converts formula to a(n unevaluated) list of numeric indices
# e.g. a*b --> list(1,2,1:2)
formula2idx <- function(x, asCall = TRUE) {
if(!is.call(x) || !inherits(x, "formula")) stop("'x' is not a formula")
fac <- attr(delete.response(terms(x)), "factors")
dimnames(fac) <- NULL
ret <- apply(fac > 0L, 2L, which)
if(asCall) as.call(c(as.name("list"), ret)) else ret
}
formula_margin_check <- function(j, m) {
stopifnot(is.logical(j))
!any(m[!j, j], na.rm = TRUE)
}