how to use dotplot after clmm

Question

I have tried to plot random effects after clmm, but I got error message: "Error in sort.list(y):'x' must be atomic for `sort.list' Have you called 'sort' an a list?" The following codes are representative of my actual codes.

library(ordinal)
library(lattice)

###I am using the wine data in the ordinal package
d<-data.frame(wine)
result<-clmm(rating~1+temp+contact+(1+temp|judge), data=d)
###In my actual codes, I put "as.ordered(rating)" instead of "rating".

Then, I am trying to plot random effects of temp and judge:

dotplot(ranef(result, condVar=TRUE))

Then, the error message appears: "Error in sort.list(y):'x' must be atomic for `sort.list' Have you called 'sort' on a list?"

My best guess is that rating is implemented as ordered under clmm, which seems for me to make sense given the error message. Yet, the guess is highly speculative, and I have no idea how to handle this situation. Specifically, what I would like to plot is the random effects of temp and judge (intercept) with their CIs. Please refer to the following plot that I generated using the codes below

result2<-lmer(as.numeric(rating)~1+temp+contact+(1+temp|judge), data=d)
dotplot(ranef(result2, condVar=TRUE))

If you could give any comment or suggestion on how to circumvent the seemingly-conflicting situation in which using clmm prevents me from using dotplot using the code above, it will be really appreciated.

Answer 1

dotplot does not know what to do with the result from ranef(clmm(...)) as it has no associated method for this result.

The 'simple' answer is that you can set the class of the result (to the same class as provided by ranef(lm34::lmer(...))) to circumvent this:

r1 <- ranef(result, condVar=TRUE)
class(r1)="ranef.mer"
dotplot(r1)

This will get you a dotplot, but without the conditional standard error bars.

To also have the conditional error bars plotted need a little more effort. The problem is that lme4::dotplot.ranef.mer looks to find the errors in an attribute named postVar, whereas clmm provides them in a differently formatted attribute called condVar. Fortunately, it is relatively simple to edit dotplot.ranef.mer to deal also with clmm objects. we can do this by adding one line to the function:

if (!is.null(cv <- attr(xt, "condVar")))  se <- as.vector(unlist(cv))

Then we can get the plot using

result<-clmm(rating~1+temp+contact+(1+temp|judge), data=d)
r1 <- ranef(result, condVar=TRUE)
dp(r1, scales = list(x = list(relation = 'free')))

Here's the whole of the function to draw this plot. This is a verbatim copy of dotplot.ranef.mer with the one added line as mentioned above

dp <- function (x, data, main = TRUE, ...) 
{
  prepanel.ci <- function(x, y, se, subscripts, ...) {
    if (is.null(se)) 
      return(list())
    x <- as.numeric(x)
    hw <- 1.96 * as.numeric(se[subscripts])
    list(xlim = range(x - hw, x + hw, finite = TRUE))
  }
  panel.ci <- function(x, y, se, subscripts, pch = 16, horizontal = TRUE, 
                       col = dot.symbol$col, lty = dot.line$lty, lwd = dot.line$lwd, 
                       col.line = dot.line$col, levels.fos = unique(y), groups = NULL, 
                       ...) {
    x <- as.numeric(x)
    y <- as.numeric(y)
    dot.line <- trellis.par.get("dot.line")
    dot.symbol <- trellis.par.get("dot.symbol")
    sup.symbol <- trellis.par.get("superpose.symbol")
    panel.abline(h = levels.fos, col = col.line, lty = lty, 
                 lwd = lwd)
    panel.abline(v = 0, col = col.line, lty = lty, lwd = lwd)
    if (!is.null(se)) {
      se <- as.numeric(se[subscripts])
      panel.segments(x - 1.96 * se, y, x + 1.96 * se, y, 
                     col = "black")
    }
    panel.xyplot(x, y, pch = pch, ...)
  }
  f <- function(nx, ...) {
    xt <- x[[nx]]
    ss <- stack(xt)
    mtit <- if (main) 
      nx
    ss$ind <- factor(as.character(ss$ind), levels = colnames(xt))
    ss$.nn <- rep.int(reorder(factor(rownames(xt)), xt[[1]], 
                              FUN = mean, sort = sort), ncol(xt))
    se <- NULL
    if (!is.null(pv <- attr(xt, "postVar"))) 
      se <- unlist(lapply(1:(dim(pv)[1]), function(i) sqrt(pv[i, i, ])))
#############################################################
# Next line is the inseerted line to deal with clmm objects
#############################################################
    if (!is.null(cv <- attr(xt, "condVar"))) se <- as.vector(unlist(cv))
    dotplot(.nn ~ values | ind, ss, se = se, prepanel = prepanel.ci, 
            panel = panel.ci, xlab = NULL, main = mtit, ...)
  }
  setNames(lapply(names(x), f, ...), names(x))
}

Answer 2

I started out trying to hack the guts of lme4:::dotplot.ranef.mer, but it turned out to be easier to do this with reshape2::melt and ggplot2. Sorry if you have a really strong preference for lattice ...

Note by the way that, technically, these aren't "confidence intervals" (since the predicted conditional modes are not "estimates" in the technical, frequentist sense) ...

Rearrange point estimates and SEs of conditional modes into a useful shape:

library(reshape2)
melt.ranef.clmm <- function(re,cv) {
  reList <- lapply(re,
        function(x) {
           x$id <- reorder(factor(rownames(x)),x[,1])
           return(melt(x,id.var="id"))
        })
  cvList <- lapply(cv,melt,id.var=NULL,value.name="se")
  mm <- Map(cbind,reList,cvList)
  return(mm)
}

Apply this to a model (m1 is a fitted model):

ss <- melt.ranef.clmm(ranef(m1),condVar(m1))

Plot:

library(ggplot2); theme_set(theme_bw())
ggplot(ss[[1]],aes(value,id))+
  geom_errorbarh(aes(xmin=value-1.96*se,xmax=value+1.96*se),
                 height=0)+
  ggtitle(names(ss)[1])+
  geom_point(colour="blue")+
  facet_wrap(~variable,scale="free_x")