New error when estimating emmeans with mode mean.class from clmm fit

[qdread]

Hi Russ, about a year ago you helped me out with issue #387 to fix a bug when estimating marginal means with mode = 'mean.class' from a model fit with ordinal::clmm() that has nested fixed effects, and random effects. Everything worked great. In the interim, I updated emmeans to the latest release and then had to revisit the code from a year ago. Now, it returns a new error, specifically Error in object@V[disp, disp, drop = FALSE] : (subscript) logical subscript too long. Note: the ordinal package has not had a new release since last year so I am still using the same version of ordinal that interacted correctly with emmeans when I ran the code a year ago.

Here is a modified reproducible example from that previous issue.

library(emmeans)
library(ordinal)

example_data <- structure(list(Genotype = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 
                                                      1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
                                                      2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 
                                                      1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
                                                      2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 
                                                      3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                      1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
                                                      2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), levels = c("g1", 
                                                                                                                              "g2", "g3"), class = "factor"), Strain = structure(c(3L, 3L, 
                                                                                                                                                                                   3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 
                                                                                                                                                                                   1L, 1L, 1L, 1L, 2L, 2L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 
                                                                                                                                                                                   6L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 
                                                                                                                                                                                   1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 5L, 5L, 
                                                                                                                                                                                   5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 3L, 3L, 3L, 
                                                                                                                                                                                   3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 
                                                                                                                                                                                   1L, 1L, 1L, 2L, 2L, 2L, 2L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 
                                                                                                                                                                                   6L, 6L, 6L), levels = c("a", "b", "c", "d", "e", "f"), class = "factor"), 
                               trial = c("1", "1", "1", "1", "1", "1", "1", "1", "1", "1", 
                                         "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", 
                                         "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", 
                                         "1", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", 
                                         "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", 
                                         "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", 
                                         "2", "2", "2", "2", "2", "2", "2", "2", "2", "3", "3", "3", 
                                         "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", 
                                         "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", 
                                         "3", "3", "3", "3", "3", "3", "3", "3", "3", "3", "3"), score_ordinal = structure(c(2L, 
                                                                                                                             2L, 5L, 6L, 3L, 1L, 4L, 3L, 3L, 2L, 2L, 2L, 2L, 3L, 2L, 2L, 
                                                                                                                             5L, 6L, 6L, 1L, 4L, 4L, 6L, 2L, 2L, 5L, 6L, 5L, 2L, 5L, 3L, 
                                                                                                                             6L, 6L, 5L, 5L, 3L, 3L, 4L, 4L, 1L, 5L, 1L, 3L, 3L, 4L, 5L, 
                                                                                                                             2L, 4L, 3L, 4L, 4L, 6L, 3L, 2L, 2L, 4L, 4L, 3L, 6L, 3L, 3L, 
                                                                                                                             6L, 5L, 3L, 2L, 5L, 4L, 4L, 4L, 6L, 3L, 3L, 3L, 5L, 3L, 4L, 
                                                                                                                             3L, 2L, 5L, 2L, 6L, 4L, 4L, 2L, 6L, 6L, 4L, 2L, 3L, 5L, 6L, 
                                                                                                                             6L, 4L, 5L, 4L, 4L, 4L, 5L, 6L, 6L, 5L, 6L, 4L, 3L, 3L, 2L, 
                                                                                                                             5L, 2L, 3L, 3L, 5L, 5L, 5L, 5L, 4L, 6L, 6L), levels = c("0", 
                                                                                                                                                                                     "1", "2", "3", "4", "5"), class = c("ordered", "factor")) 
                               ), row.names = c(NA, 
                                                                                                                                 -117L), class = "data.frame")

ordinal_fit <- clmm(score_ordinal ~ Genotype + Strain + (1|trial), link = 'logit', data = example_data)

emmeans(ordinal_fit, ~ Genotype, mode = 'mean.class')

The output indicates that the nested structure was correctly detected but the error is new and didn't occur using the release of emmeans that was current in November 2022.

NOTE: A nesting structure was detected in the fitted model:
    Strain %in% Genotype
Error in object@V[disp, disp, drop = FALSE] : 
  (subscript) logical subscript too long

Thanks in advance for your help!

[rvlenth]

Oh my. I see where this is happening, and tried to do a patch, but it screws up something that comes later. I'll have to look at when and why I made the change that affects you. It's very mysterious, and am not sure how long this will take.

[rvlenth]

The really discouraging thing is that the code where the error occurs is in the fix for #287, and that the example code there that this supposedly fixes is once again broken! It creates an error in exactly the same line that you are experiencing here.

[rvlenth]

Well, I feel dumb, but it was a simple fix... Now it works again:

> emmeans(ordinal_fit, ~ Genotype, mode = 'mean.class')
NOTE: A nesting structure was detected in the fitted model:
    Strain %in% Genotype
 Genotype mean.class    SE  df asymp.LCL asymp.UCL
 g2             4.19 0.297 Inf      3.61      4.77
 g1             3.48 0.255 Inf      2.97      3.98
 g3             4.11 0.263 Inf      3.59      4.63

Results are averaged over the levels of: Strain 
Confidence level used: 0.95 

I had coded V[disp, disp] when I meant V[estble, estble]. It is very unclear to me why this ever worked with a nested model. I must have been out of synch with things somehow.

[qdread]

Thanks, Russ!

[rvlenth]

Gonna go ahead and close this. I guess it may reopen when the next gremlin comes along, but hoping that isn't soon.