## These are the functions for fitting the multilevel index using lme4,
## calculating the percentage of the variance at each level and
## also the holdback scores
.mid <- function(data, vars, levels, expected = FALSE, nsims = 100) {
if (is.character((levels))) {
ifelse (all(levels %in% names(data)), levels <- match(levels, names(data)),
stop("Higher level grouping variable not found"))
}
if (anyNA(data[,levels])) stop("Data contain NAs")
if (any(sapply(data[, levels], is.numeric)))
warning("The levels contain numeric data")
id <- .idx(data, vars, expected, nsims)
ols <- attr(id, "ols")
vv <- attr(id, "vars")
df <- attr(id, "data")
lvls <- names(data)[levels]
data <- data.frame(y = data[, vars[1]], x = data[, vars[2]], data[, levels])
names(data) <- c("y", "x", lvls)
data$y <- data$y / sum(data$y)
data$x <- data$x / sum(data$x)
f <- "y ~ 0"
for(k in 3: ncol(data)) {
f <- paste(f, "+", paste0("(1|", names(data)[k],")"))
}
mlm <- lme4::lmer(f, data=data, offset=data$x)
attributes(id) <- list(ols = ols, vars = vv, levels = lvls,
mlm = mlm, variance = .varshare(mlm),
holdback = holdback(mlm),
data = df)
class(id) <- "index"
return(id)
}
.mlvar <- function(mlm) {
vv <- lme4::VarCorr(mlm)
res <- attr(vv, "sc")
for(i in 1: length(vv)) {
res <- c(res,attr(vv[[i]], "stddev"))
}
res <- res^2
names(res) <- c("Base",names(vv))
return(res)
}
.varshare <- function(mlm) {
mlvar <- .mlvar(mlm)
mlvar <- mlvar/sum(mlvar)*100
return(round(mlvar,2))
}
#' Holdback scores
#'
#' Calculates the holdback scores for a multilevel index of dissimilarity
#'
#' For the index of dissimilarity (ID), the residuals are the differences
#' between the share of the Y population and the share of the X population per
#' neighbourhood. For the multilevel index, the residuals are estimated at and
#' partitioned between each level of the model. The holdback scores take each
#' level in the model in turn and set the residuals (the effects) at that
#' level to zero, then recalculating the ID on that basis and recording
#' the percentage change in the original value that occurs. The holdback
#' scores are calculated automatically as part of the function \code{\link{id}}
#' and can be viewed through \code{print(index)}, where \code{index} is the
#' object returned by the function, or as \code{attr(index, "holdback")}.
#'
#' @param mlm an object of class \code{lmerMod} generated by the
#' \code{lme4} package
#' @return a numeric vector containing the holdback scores
#' @seealso \code{\link{id}} \code{\link{print.index}}
holdback <- function(mlm) {
rvals <- .rvals(mlm)
k <- ncol(rvals)
hb <- rep(NA, k)
rwsm <- rowSums(rvals)
ID <- sum(abs(rwsm))
for(i in 1:k) {
ifelse(k > 2, rwsm <- rowSums(rvals[,-i]), rwsm <- rvals[,-i])
hb[i] <- sum(abs(rwsm))/ID - 1
}
names(hb) <- colnames(rvals)
return(round(hb*100,2))
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.