search

rvlenth / emmeans

Estimated marginal means
https://rvlenth.github.io/emmeans/
340 stars 30 forks source link

Ordinal logistic regression [GEE]: ordLORgee; mode="prob" - all estimates have same value, joint_tests() fails. What does it mean? #481

Closed Generalized closed 2 months ago

Generalized commented 2 months ago

Hello,

I have a 2-arm study with 3 repeated observation timepoints and ordinal response at each timepoint. I attach the data on the bottom.

I analyze it via GEE using the multgee and repolr packages.

Let's fit a model

library(multgee)

fitmod <- ordLORgee(formula = ODIPain ~ Visit_nOrd * Arm + Age_centered + Age_centered : Visit_nOrd,
                    data = data, id = PatientId, repeated = as.numeric(Visit_nOrd), LORstr = "uniform")
coef(summary(fitmod))

                                Estimate san.se  san.z Pr(>|san.z|)
beta10                           -1.2882 0.2588 -4.978       0.0000
beta20                           -0.3159 0.2343 -1.348       0.1776
beta30                            1.1004 0.2364  4.654       0.0000
beta40                            2.1376 0.2743  7.794       0.0000
beta50                            3.6202 0.4557  7.944       0.0000
Visit_nOrdMonth 12               -0.4471 0.2660 -1.681       0.0927
Visit_nOrdMonth 20               -0.5431 0.2575 -2.109       0.0349
ArmB                             -0.0502 0.3247 -0.155       0.8771
Age_centered                      0.0139 0.0150  0.929       0.3527
Visit_nOrdMonth 12:ArmB           0.6400 0.3670  1.744       0.0812
Visit_nOrdMonth 20:ArmB           0.8761 0.3510  2.496       0.0126
Visit_nOrdMonth 12:Age_centered   0.0289 0.0151  1.917       0.0552
Visit_nOrdMonth 20:Age_centered   0.0167 0.0127  1.310       0.1902

It converged. Now let's make the grid and emmeans:

> (lord_grid <- qdrg(formula = ODIPain ~ Visit_nOrd * Arm + Age_centered + Age_centered : Visit_nOrd,
+                    data=data, 
+                    coef = fitmod$coefficients, 
+                    vcov = fitmod$robust.variance, 
+                    df = Inf, 
+                    ordinal = list(dim=6),
+                    link = "logit"))

'emmGrid' object with variables:
    Visit_nOrd = Month 6, Month 12, Month 20
    Arm = A, B
    Age_centered = -2.4336e-15
Transformation: “logit” 
> 
> (lord_em <- emmeans(lord_grid, specs = ~Arm * Visit_nOrd , data=data))
 Arm Visit_nOrd emmean    SE  df asymp.LCL asymp.UCL
 A   Month 6    -1.051 0.248 Inf     -1.54    -0.565
 B   Month 6    -1.101 0.521 Inf     -2.12    -0.079
 A   Month 12   -1.498 0.434 Inf     -2.35    -0.647
 B   Month 12   -0.908 0.530 Inf     -1.95     0.131
 A   Month 20   -1.594 0.431 Inf     -2.44    -0.749
 B   Month 20   -0.768 0.522 Inf     -1.79     0.255

Results are given on the logit (not the response) scale. 
Confidence level used: 0.95

And test some contrasts I was interested in:

> update(contrast(lord_em,
+                 list("Month 6  : A vs. B" = c(-1, 1, 0, 0, 0, 0),
+                      "Month 12 : A vs. B" = c( 0, 0,-1, 1, 0, 0),
+                      "Month 20 : A vs. B" = c( 0, 0, 0, 0,-1, 1))),
+        adjust="none", level = 0.95, infer = c(TRUE, TRUE))
 contrast           estimate    SE  df asymp.LCL asymp.UCL z.ratio p.value
 Month 6  : A vs. B   -0.050 0.325 Inf    -0.687     0.586  -0.155  0.8770
 Month 12 : A vs. B    0.590 0.365 Inf    -0.126     1.305   1.616  0.1060
 Month 20 : A vs. B    0.826 0.350 Inf     0.139     1.513   2.357  0.0180

Results are given on the log odds ratio (not the response) scale. 
Confidence level used: 0.95 

> joint_tests(lord_grid)
 model term     df1 df2 F.ratio Chisq p.value
 Visit_nOrd       2 Inf   0.278  0.56  0.7570
 Arm              1 Inf   2.544  2.54  0.1110
 Visit_nOrd:Arm   2 Inf   3.134  6.27  0.0440

OK. Now let's switch to other modes:

mode = "cum.prob" emmeans worked, yet joint_tests failed. What could cause that?

> 
> (lord_grid <- qdrg(formula = ODIPain ~ Visit_nOrd * Arm + Age_centered + Age_centered : Visit_nOrd,
+                    data=data, 
+                    coef = fitmod$coefficients, 
+                    vcov = fitmod$robust.variance, 
+                    df = Inf, 
+                    ordinal = list(dim=6, mode="cum.prob")))
'emmGrid' object with variables:
    Visit_nOrd = Month 6, Month 12, Month 20
    Arm = A, B
    Age_centered = -2.4336e-15
    cut = multivariate response levels: 1|2, 2|3, 3|4, 4|5, 5|6
> 
> (lord_em <- emmeans(lord_grid, specs = ~Arm * Visit_nOrd , data=data))
 Arm Visit_nOrd emmean    SE  df asymp.LCL asymp.UCL
 A   Month 6     1.051 0.248 Inf     0.565      1.54
 B   Month 6     1.101 0.521 Inf     0.079      2.12
 A   Month 12    1.498 0.434 Inf     0.647      2.35
 B   Month 12    0.908 0.530 Inf    -0.131      1.95
 A   Month 20    1.594 0.431 Inf     0.749      2.44
 B   Month 20    0.768 0.522 Inf    -0.255      1.79

Results are averaged over the levels of: cut 
Confidence level used: 0.95 
> 
> update(contrast(lord_em,
+                 list("Month 6  : A vs. B" = c(-1, 1, 0, 0, 0, 0),
+                      "Month 12 : A vs. B" = c( 0, 0,-1, 1, 0, 0),
+                      "Month 20 : A vs. B" = c( 0, 0, 0, 0,-1, 1))),
+        adjust="none", level = 0.95, infer = c(TRUE, TRUE))
 contrast           estimate    SE  df asymp.LCL asymp.UCL z.ratio p.value
 Month 6  : A vs. B    0.050 0.325 Inf    -0.586     0.687   0.155  0.8770
 Month 12 : A vs. B   -0.590 0.365 Inf    -1.305     0.126  -1.616  0.1060
 Month 20 : A vs. B   -0.826 0.350 Inf    -1.513    -0.139  -2.357  0.0180

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

> joint_tests(lord_grid)
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 2.87943e-19
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 1.65641e-18
Error in rbind(deparse.level, ...) : 
  numbers of columns of arguments do not match

This time joint_tests() failed.

mode = "prob"

> (lord_grid <- qdrg(formula = ODIPain ~ Visit_nOrd * Arm + Age_centered + Age_centered : Visit_nOrd,
+                    data=data, 
+                    coef = fitmod$coefficients, 
+                    vcov = fitmod$robust.variance, 
+                    df = Inf, 
+                    ordinal = list(dim=6, mode="prob")))
'emmGrid' object with variables:
    Visit_nOrd = Month 6, Month 12, Month 20
    Arm = A, B
    Age_centered = -2.4336e-15
    ODIPain = multivariate response levels: 1, 2, 3, 4, 5, 6
> 
> (lord_em <- emmeans(lord_grid, specs = ~Arm * Visit_nOrd , data=data))
 Arm Visit_nOrd  prob SE  df asymp.LCL asymp.UCL
 A   Month 6    0.167  0 Inf     0.167     0.167
 B   Month 6    0.167  0 Inf     0.167     0.167
 A   Month 12   0.167  0 Inf     0.167     0.167
 B   Month 12   0.167  0 Inf     0.167     0.167
 A   Month 20   0.167  0 Inf     0.167     0.167
 B   Month 20   0.167  0 Inf     0.167     0.167

Results are averaged over the levels of: ODIPain 
Confidence level used: 0.95 
> 
> update(contrast(lord_em,
+                 list("Month 6  : A vs. B" = c(-1, 1, 0, 0, 0, 0),
+                      "Month 12 : A vs. B" = c( 0, 0,-1, 1, 0, 0),
+                      "Month 20 : A vs. B" = c( 0, 0, 0, 0,-1, 1))),
+        adjust="none", level = 0.95, infer = c(TRUE, TRUE))
 contrast           estimate SE  df asymp.LCL asymp.UCL z.ratio p.value
 Month 6  : A vs. B        0  0 Inf         0         0     NaN     NaN
 Month 12 : A vs. B        0  0 Inf         0         0     NaN     NaN
 Month 20 : A vs. B        0  0 Inf         0         0     NaN     NaN

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

> joint_tests(lord_grid)
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 2.87897e-17
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 1.89748e-17
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 1.73311e-17
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 3.28385e-19
Error in rbind(deparse.level, ...) : 
  numbers of columns of arguments do not match

What does it mean to have all estimated equal to same value 0.167 here? joint_tests() failed as well. The behaviour of both tested packages is the same here.

The data:

> dput(data)
structure(list(PatientId = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 
3L, 3L, 3L, 4L, 4L, 4L, 5L, 5L, 5L, 6L, 6L, 6L, 7L, 7L, 7L, 8L, 
8L, 8L, 9L, 9L, 9L, 10L, 10L, 10L, 11L, 11L, 11L, 12L, 12L, 12L, 
13L, 13L, 13L, 14L, 14L, 14L, 15L, 15L, 15L, 16L, 16L, 16L, 17L, 
17L, 17L, 18L, 18L, 18L, 19L, 19L, 19L, 20L, 20L, 20L, 21L, 21L, 
21L, 22L, 22L, 22L, 23L, 23L, 23L, 24L, 24L, 24L, 25L, 25L, 25L, 
26L, 26L, 26L, 27L, 27L, 27L, 28L, 28L, 28L, 29L, 29L, 29L, 30L, 
30L, 30L, 31L, 31L, 31L, 32L, 32L, 32L, 33L, 33L, 33L, 34L, 34L, 
34L, 35L, 35L, 35L, 36L, 36L, 36L, 37L, 37L, 37L, 38L, 38L, 38L, 
39L, 39L, 39L, 40L, 40L, 40L, 41L, 41L, 41L, 42L, 42L, 42L, 43L, 
43L, 43L, 44L, 44L, 44L, 45L, 45L, 45L, 46L, 46L, 46L, 47L, 47L, 
47L, 48L, 48L, 48L, 49L, 49L, 49L, 50L, 50L, 50L, 51L, 51L, 51L, 
52L, 52L, 52L, 53L, 53L, 53L, 54L, 54L, 54L, 55L, 55L, 55L, 56L, 
56L, 56L, 57L, 57L, 57L, 58L, 58L, 58L, 59L, 59L, 59L, 60L, 60L, 
60L, 61L, 61L, 61L, 62L, 62L, 62L, 63L, 63L, 63L, 64L, 64L, 64L, 
65L, 65L, 65L, 66L, 66L, 66L, 67L, 67L, 67L, 68L, 68L, 68L, 69L, 
69L, 69L, 70L, 70L, 70L, 71L, 71L, 71L, 72L, 72L, 72L, 73L, 73L, 
73L, 74L, 74L, 74L, 75L, 75L, 75L, 76L, 76L, 76L, 77L, 77L, 77L, 
78L, 78L, 78L, 79L, 79L, 79L, 80L, 80L, 80L, 81L, 81L, 81L, 82L, 
82L, 82L, 83L, 83L, 83L, 84L, 84L, 84L, 85L, 85L, 85L, 86L, 86L, 
86L, 87L, 87L, 87L, 88L, 88L, 88L, 89L, 89L, 89L, 90L, 90L, 90L, 
91L, 91L, 91L, 92L, 92L, 92L, 93L, 93L, 93L, 94L, 94L, 94L, 95L, 
95L, 95L, 96L, 96L, 96L, 97L, 97L, 97L, 98L, 98L, 98L, 99L, 99L, 
99L, 100L, 100L, 100L, 101L, 101L, 101L, 102L, 102L, 102L, 103L, 
103L, 103L, 104L, 104L, 104L, 105L, 105L, 105L, 106L, 106L, 106L, 
107L, 107L, 107L, 108L, 108L, 108L, 109L, 109L, 109L, 110L, 110L, 
110L, 111L, 111L, 111L), levels = c("4x9x110004", "4x9x110007", 
"4x9x110008", "4x9x110010", "4x9x110012", "4x9x110014", "4x9x110017", 
"4x9x110018", "4x9x110019", "4x9x110021", "4x9x110023", "4x9x110024", 
"4x9x110029", "4x9x110030", "4x9x110031", "4x9x110033", "4x9x110034", 
"4x9x110036", "4x9x110037", "4x9x110038", "4x9x110041", "4x9x110042", 
"4x9x110043", "4x9x110044", "4x9x110045", "4x9x110047", "4x9x110048", 
"4x9x110049", "4x9x110050", "4x9x110051", "4x9x110052", "4x9x110054", 
"4x9x110055", "4x9x110056", "4x9x110057", "4x9x110058", "4x9x110059", 
"4x9x110061", "4x9x110062", "4x9x110063", "4x9x110064", "4x9x110066", 
"4x9x110068", "4x9x110069", "4x9x110070", "4x9x110071", "4x9x110072", 
"4x9x110073", "4x9x110074", "4x9x110075", "4x9x110076", "4x9x110077", 
"4x9x120001", "4x9x120003", "4x9x120005", "4x9x120009", "4x9x120014", 
"4x9x120015", "4x9x120017", "4x9x120020", "4x9x120021", "4x9x120022", 
"4x9x120026", "4x9x120027", "4x9x120028", "4x9x120029", "4x9x120030", 
"4x9x120031", "4x9x120032", "4x9x120033", "4x9x120034", "4x9x120035", 
"4x9x120036", "4x9x120037", "4x9x120038", "4x9x120039", "4x9x120040", 
"4x9x120041", "4x9x120042", "4x9x120043", "4x9x120045", "4x9x120047", 
"4x9x120048", "4x9x120049", "4x9x120050", "4x9x120051", "4x9x120052", 
"4x9x130001", "4x9x130004", "4x9x130005", "4x9x130006", "4x9x130009", 
"4x9x130012", "4x9x130013", "4x9x130014", "4x9x130015", "4x9x130016", 
"4x9x130017", "4x9x130018", "4x9x130019", "4x9x130020", "4x9x130022", 
"4x9x130023", "4x9x130024", "4x9x130025", "4x9x130026", "4x9x130027", 
"4x9x130028", "4x9x130030", "4x9x130031", "4x9x130032"), class = "factor"), 
    ODIPain = structure(c(3L, 2L, 2L, 1L, 1L, 3L, 2L, 4L, 1L, 
    3L, 1L, 1L, 4L, 3L, 4L, 3L, 2L, 3L, 4L, 4L, 4L, 2L, 3L, 2L, 
    3L, 3L, 4L, 5L, 3L, 4L, 2L, 1L, 1L, 5L, 4L, 4L, 5L, 3L, 2L, 
    6L, 2L, 3L, 2L, 3L, 3L, 3L, 3L, 4L, 1L, 5L, 4L, 2L, 2L, 2L, 
    3L, 6L, 4L, 3L, 3L, 2L, 1L, 1L, 1L, 3L, 4L, 5L, 3L, 3L, 4L, 
    3L, 5L, 6L, 2L, 5L, 4L, 5L, 5L, 5L, 5L, 2L, 3L, 2L, 3L, 4L, 
    5L, 3L, 3L, 3L, 3L, 3L, 5L, 1L, 3L, 1L, 4L, 4L, 3L, 3L, 4L, 
    3L, 3L, 1L, 3L, 3L, 3L, 1L, 1L, 2L, 5L, 4L, 5L, 3L, 4L, 4L, 
    5L, 4L, 4L, 1L, 2L, 3L, 1L, 2L, 1L, 2L, 3L, 5L, 1L, 1L, 1L, 
    2L, 1L, 1L, 3L, 3L, 3L, 2L, 2L, 3L, 3L, 2L, 4L, 1L, 3L, 1L, 
    2L, 2L, 1L, 1L, 1L, 1L, 2L, 3L, 1L, 4L, 1L, 1L, 6L, 6L, 6L, 
    4L, 4L, 3L, 3L, 5L, 3L, 3L, 2L, 1L, 3L, 1L, 1L, 3L, 1L, 1L, 
    1L, 4L, 2L, 4L, 5L, 6L, 5L, 4L, 3L, 4L, 5L, 5L, 3L, 3L, 1L, 
    3L, 4L, 3L, 1L, 1L, 3L, 3L, 2L, 2L, 1L, 1L, 2L, 4L, 6L, 4L, 
    3L, 3L, 2L, 3L, 3L, 4L, 2L, 3L, 3L, 3L, 4L, 4L, 3L, 3L, 5L, 
    3L, 3L, 4L, 3L, 1L, 3L, 3L, 5L, 4L, 2L, 1L, 1L, 2L, 2L, 2L, 
    1L, 2L, 4L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 3L, 1L, 3L, 3L, 3L, 
    4L, 5L, 3L, 4L, 3L, 3L, 3L, 3L, 1L, 2L, 1L, 3L, 2L, 1L, 2L, 
    2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 3L, 1L, 3L, 4L, 1L, 1L, 1L, 
    3L, 2L, 2L, 1L, 1L, 1L, 3L, 3L, 2L, 1L, 1L, 1L, 3L, 4L, 3L, 
    4L, 6L, 6L, 4L, 4L, 3L, 3L, 5L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 
    2L, 2L, 3L, 3L, 2L, 3L, 3L, 5L, 4L, 3L, 5L, 3L, 3L, 3L, 3L, 
    5L, 5L, 4L, 2L, 3L, 2L, 2L, 1L, 2L), levels = c("0", "1", 
    "2", "3", "4", "5"), class = c("ordered", "factor")), Arm = structure(c(2L, 
    2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 
    1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 
    1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L), levels = c("A", "B"), class = "factor"), Age_centered = c(-6.15315315315316, 
    -6.15315315315316, -6.15315315315316, 12.8468468468468, 12.8468468468468, 
    12.8468468468468, -9.15315315315316, -9.15315315315316, -9.15315315315316, 
    14.8468468468468, 14.8468468468468, 14.8468468468468, 12.8468468468468, 
    12.8468468468468, 12.8468468468468, 2.84684684684684, 2.84684684684684, 
    2.84684684684684, -10.1531531531532, -10.1531531531532, -10.1531531531532, 
    -18.1531531531532, -18.1531531531532, -18.1531531531532, 
    -1.15315315315316, -1.15315315315316, -1.15315315315316, 
    8.84684684684684, 8.84684684684684, 8.84684684684684, -17.1531531531532, 
    -17.1531531531532, -17.1531531531532, 13.8468468468468, 13.8468468468468, 
    13.8468468468468, 9.84684684684684, 9.84684684684684, 9.84684684684684, 
    17.8468468468468, 17.8468468468468, 17.8468468468468, -19.1531531531532, 
    -19.1531531531532, -19.1531531531532, -7.15315315315316, 
    -7.15315315315316, -7.15315315315316, -10.1531531531532, 
    -10.1531531531532, -10.1531531531532, -19.1531531531532, 
    -19.1531531531532, -19.1531531531532, -7.15315315315316, 
    -7.15315315315316, -7.15315315315316, 0.846846846846844, 
    0.846846846846844, 0.846846846846844, -17.1531531531532, 
    -17.1531531531532, -17.1531531531532, 5.84684684684684, 5.84684684684684, 
    5.84684684684684, -25.1531531531532, -25.1531531531532, -25.1531531531532, 
    -1.15315315315316, -1.15315315315316, -1.15315315315316, 
    -15.1531531531532, -15.1531531531532, -15.1531531531532, 
    4.84684684684684, 4.84684684684684, 4.84684684684684, 1.84684684684684, 
    1.84684684684684, 1.84684684684684, 12.8468468468468, 12.8468468468468, 
    12.8468468468468, -11.1531531531532, -11.1531531531532, -11.1531531531532, 
    5.84684684684684, 5.84684684684684, 5.84684684684684, -6.15315315315316, 
    -6.15315315315316, -6.15315315315316, -0.153153153153156, 
    -0.153153153153156, -0.153153153153156, 20.8468468468468, 
    20.8468468468468, 20.8468468468468, 5.84684684684684, 5.84684684684684, 
    5.84684684684684, -0.153153153153156, -0.153153153153156, 
    -0.153153153153156, 12.8468468468468, 12.8468468468468, 12.8468468468468, 
    -19.1531531531532, -19.1531531531532, -19.1531531531532, 
    -11.1531531531532, -11.1531531531532, -11.1531531531532, 
    1.84684684684684, 1.84684684684684, 1.84684684684684, 0.846846846846844, 
    0.846846846846844, 0.846846846846844, -21.1531531531532, 
    -21.1531531531532, -21.1531531531532, 9.84684684684684, 9.84684684684684, 
    9.84684684684684, 15.8468468468468, 15.8468468468468, 15.8468468468468, 
    14.8468468468468, 14.8468468468468, 14.8468468468468, -12.1531531531532, 
    -12.1531531531532, -12.1531531531532, -11.1531531531532, 
    -11.1531531531532, -11.1531531531532, -9.15315315315316, 
    -9.15315315315316, -9.15315315315316, 5.84684684684684, 5.84684684684684, 
    5.84684684684684, -4.15315315315316, -4.15315315315316, -4.15315315315316, 
    12.8468468468468, 12.8468468468468, 12.8468468468468, 1.84684684684684, 
    1.84684684684684, 1.84684684684684, -7.15315315315316, -7.15315315315316, 
    -7.15315315315316, -3.15315315315316, -3.15315315315316, 
    -3.15315315315316, 7.84684684684684, 7.84684684684684, 7.84684684684684, 
    0.846846846846844, 0.846846846846844, 0.846846846846844, 
    -4.15315315315316, -4.15315315315316, -4.15315315315316, 
    5.84684684684684, 5.84684684684684, 5.84684684684684, -0.153153153153156, 
    -0.153153153153156, -0.153153153153156, 1.84684684684684, 
    1.84684684684684, 1.84684684684684, -7.15315315315316, -7.15315315315316, 
    -7.15315315315316, 1.84684684684684, 1.84684684684684, 1.84684684684684, 
    -9.15315315315316, -9.15315315315316, -9.15315315315316, 
    6.84684684684684, 6.84684684684684, 6.84684684684684, 9.84684684684684, 
    9.84684684684684, 9.84684684684684, 17.8468468468468, 17.8468468468468, 
    17.8468468468468, 5.84684684684684, 5.84684684684684, 5.84684684684684, 
    9.84684684684684, 9.84684684684684, 9.84684684684684, -10.1531531531532, 
    -10.1531531531532, -10.1531531531532, -5.15315315315316, 
    -5.15315315315316, -5.15315315315316, 18.8468468468468, 18.8468468468468, 
    18.8468468468468, 21.8468468468468, 21.8468468468468, 21.8468468468468, 
    -0.153153153153156, -0.153153153153156, -0.153153153153156, 
    2.84684684684684, 2.84684684684684, 2.84684684684684, -8.15315315315316, 
    -8.15315315315316, -8.15315315315316, -5.15315315315316, 
    -5.15315315315316, -5.15315315315316, 5.84684684684684, 5.84684684684684, 
    5.84684684684684, 2.84684684684684, 2.84684684684684, 2.84684684684684, 
    -15.1531531531532, -15.1531531531532, -15.1531531531532, 
    2.84684684684684, 2.84684684684684, 2.84684684684684, 25.8468468468468, 
    25.8468468468468, 25.8468468468468, -11.1531531531532, -11.1531531531532, 
    -11.1531531531532, 27.8468468468468, 27.8468468468468, 27.8468468468468, 
    2.84684684684684, 2.84684684684684, 2.84684684684684, 20.8468468468468, 
    20.8468468468468, 20.8468468468468, -0.153153153153156, -0.153153153153156, 
    -0.153153153153156, -2.15315315315316, -2.15315315315316, 
    -2.15315315315316, 12.8468468468468, 12.8468468468468, 12.8468468468468, 
    -0.153153153153156, -0.153153153153156, -0.153153153153156, 
    0.846846846846844, 0.846846846846844, 0.846846846846844, 
    11.8468468468468, 11.8468468468468, 11.8468468468468, -8.15315315315316, 
    -8.15315315315316, -8.15315315315316, 3.84684684684684, 3.84684684684684, 
    3.84684684684684, 22.8468468468468, 22.8468468468468, 22.8468468468468, 
    5.84684684684684, 5.84684684684684, 5.84684684684684, 12.8468468468468, 
    12.8468468468468, 12.8468468468468, 4.84684684684684, 4.84684684684684, 
    4.84684684684684, 11.8468468468468, 11.8468468468468, 11.8468468468468, 
    -5.15315315315316, -5.15315315315316, -5.15315315315316, 
    -17.1531531531532, -17.1531531531532, -17.1531531531532, 
    -7.15315315315316, -7.15315315315316, -7.15315315315316, 
    -16.1531531531532, -16.1531531531532, -16.1531531531532, 
    0.846846846846844, 0.846846846846844, 0.846846846846844, 
    -13.1531531531532, -13.1531531531532, -13.1531531531532, 
    -13.1531531531532, -13.1531531531532, -13.1531531531532, 
    -19.1531531531532, -19.1531531531532, -19.1531531531532, 
    -15.1531531531532, -15.1531531531532, -15.1531531531532, 
    -11.1531531531532, -11.1531531531532, -11.1531531531532, 
    -4.15315315315316, -4.15315315315316, -4.15315315315316, 
    -0.153153153153156, -0.153153153153156, -0.153153153153156, 
    -4.15315315315316, -4.15315315315316, -4.15315315315316, 
    -7.15315315315316, -7.15315315315316, -7.15315315315316), 
    Visit_nOrd = structure(c(1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 
    3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L), levels = c("Month 6", 
    "Month 12", "Month 20"), class = "factor"), ODIPain_n = c(3, 
    2, 2, 1, 1, 3, 2, 4, 1, 3, 1, 1, 4, 3, 4, 3, 2, 3, 4, 4, 
    4, 2, 3, 2, 3, 3, 4, 5, 3, 4, 2, 1, 1, 5, 4, 4, 5, 3, 2, 
    6, 2, 3, 2, 3, 3, 3, 3, 4, 1, 5, 4, 2, 2, 2, 3, 6, 4, 3, 
    3, 2, 1, 1, 1, 3, 4, 5, 3, 3, 4, 3, 5, 6, 2, 5, 4, 5, 5, 
    5, 5, 2, 3, 2, 3, 4, 5, 3, 3, 3, 3, 3, 5, 1, 3, 1, 4, 4, 
    3, 3, 4, 3, 3, 1, 3, 3, 3, 1, 1, 2, 5, 4, 5, 3, 4, 4, 5, 
    4, 4, 1, 2, 3, 1, 2, 1, 2, 3, 5, 1, 1, 1, 2, 1, 1, 3, 3, 
    3, 2, 2, 3, 3, 2, 4, 1, 3, 1, 2, 2, 1, 1, 1, 1, 2, 3, 1, 
    4, 1, 1, 6, 6, 6, 4, 4, 3, 3, 5, 3, 3, 2, 1, 3, 1, 1, 3, 
    1, 1, 1, 4, 2, 4, 5, 6, 5, 4, 3, 4, 5, 5, 3, 3, 1, 3, 4, 
    3, 1, 1, 3, 3, 2, 2, 1, 1, 2, 4, 6, 4, 3, 3, 2, 3, 3, 4, 
    2, 3, 3, 3, 4, 4, 3, 3, 5, 3, 3, 4, 3, 1, 3, 3, 5, 4, 2, 
    1, 1, 2, 2, 2, 1, 2, 4, 2, 1, 2, 2, 2, 2, 1, 3, 1, 3, 3, 
    3, 4, 5, 3, 4, 3, 3, 3, 3, 1, 2, 1, 3, 2, 1, 2, 2, 2, 1, 
    2, 2, 2, 2, 2, 3, 1, 3, 4, 1, 1, 1, 3, 2, 2, 1, 1, 1, 3, 
    3, 2, 1, 1, 1, 3, 4, 3, 4, 6, 6, 4, 4, 3, 3, 5, 3, 3, 3, 
    3, 1, 1, 1, 2, 2, 3, 3, 2, 3, 3, 5, 4, 3, 5, 3, 3, 3, 3, 
    5, 5, 4, 2, 3, 2, 2, 1, 2)), row.names = c(56L, 167L, 278L, 
57L, 168L, 279L, 1L, 112L, 223L, 58L, 169L, 280L, 2L, 113L, 224L, 
59L, 170L, 281L, 3L, 114L, 225L, 60L, 171L, 282L, 61L, 172L, 
283L, 62L, 173L, 284L, 4L, 115L, 226L, 5L, 116L, 227L, 6L, 117L, 
228L, 63L, 174L, 285L, 64L, 175L, 286L, 7L, 118L, 229L, 8L, 119L, 
230L, 65L, 176L, 287L, 9L, 120L, 231L, 66L, 177L, 288L, 67L, 
178L, 289L, 10L, 121L, 232L, 68L, 179L, 290L, 11L, 122L, 233L, 
12L, 123L, 234L, 69L, 180L, 291L, 13L, 124L, 235L, 14L, 125L, 
236L, 70L, 181L, 292L, 71L, 182L, 293L, 72L, 183L, 294L, 15L, 
126L, 237L, 16L, 127L, 238L, 17L, 128L, 239L, 73L, 184L, 295L, 
18L, 129L, 240L, 19L, 130L, 241L, 74L, 185L, 296L, 20L, 131L, 
242L, 75L, 186L, 297L, 76L, 187L, 298L, 77L, 188L, 299L, 21L, 
132L, 243L, 22L, 133L, 244L, 78L, 189L, 300L, 23L, 134L, 245L, 
24L, 135L, 246L, 79L, 190L, 301L, 25L, 136L, 247L, 80L, 191L, 
302L, 26L, 137L, 248L, 81L, 192L, 303L, 27L, 138L, 249L, 82L, 
193L, 304L, 28L, 139L, 250L, 83L, 194L, 305L, 29L, 140L, 251L, 
84L, 195L, 306L, 85L, 196L, 307L, 86L, 197L, 308L, 87L, 198L, 
309L, 30L, 141L, 252L, 88L, 199L, 310L, 89L, 200L, 311L, 31L, 
142L, 253L, 32L, 143L, 254L, 33L, 144L, 255L, 90L, 201L, 312L, 
91L, 202L, 313L, 92L, 203L, 314L, 34L, 145L, 256L, 35L, 146L, 
257L, 36L, 147L, 258L, 37L, 148L, 259L, 38L, 149L, 260L, 93L, 
204L, 315L, 94L, 205L, 316L, 95L, 206L, 317L, 39L, 150L, 261L, 
96L, 207L, 318L, 97L, 208L, 319L, 40L, 151L, 262L, 98L, 209L, 
320L, 41L, 152L, 263L, 42L, 153L, 264L, 99L, 210L, 321L, 43L, 
154L, 265L, 100L, 211L, 322L, 44L, 155L, 266L, 45L, 156L, 267L, 
101L, 212L, 323L, 46L, 157L, 268L, 102L, 213L, 324L, 47L, 158L, 
269L, 48L, 159L, 270L, 103L, 214L, 325L, 49L, 160L, 271L, 104L, 
215L, 326L, 50L, 161L, 272L, 105L, 216L, 327L, 51L, 162L, 273L, 
52L, 163L, 274L, 106L, 217L, 328L, 107L, 218L, 329L, 53L, 164L, 
275L, 54L, 165L, 276L, 108L, 219L, 330L, 55L, 166L, 277L, 109L, 
220L, 331L, 110L, 221L, 332L, 111L, 222L, 333L), class = "data.frame")
rvlenth commented 2 months ago

It means that ODIPain is the response variable, which has 6 levels. Since for each combination of other factors, the estimates in prob mode must sum to 1, then the average is always 1/6.

I suggest not averaging over ODIPain. Similarly, in other modes where cut is a generated factor, it is fairly meaningless to average over cut. If you want to filter-out the ordinal categories in your estimates, you should use mode = "latent" or mode = "mean.class".

Generalized commented 2 months ago

Ah! Got it, it was simple yet I couldn't get it until you showed me that, thank you.

Regarding the joint_tests()

> update(contrast(lord_em,
+                 list("Month 6  : A vs. B" = c(-1, 1, 0, 0, 0, 0),
+                      "Month 12 : A vs. B" = c( 0, 0,-1, 1, 0, 0),
+                      "Month 20 : A vs. B" = c( 0, 0, 0, 0,-1, 1))),
+        adjust="none", level = 0.95, infer = c(TRUE, TRUE))
 contrast           estimate    SE  df asymp.LCL asymp.UCL z.ratio p.value
 Month 6  : A vs. B    0.050 0.325 Inf    -0.586     0.687   0.155  0.8770
 Month 12 : A vs. B   -0.590 0.365 Inf    -1.305     0.126  -1.616  0.1060
 Month 20 : A vs. B   -0.826 0.350 Inf    -1.513    -0.139  -2.357  0.0180

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

> joint_tests(lord_grid)
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 2.87943e-19
Error in solve.default(zcov, z) : 
  system is computationally singular: reciprocal condition number = 1.65641e-18
Error in rbind(deparse.level, ...) : 
  numbers of columns of arguments do not match

I'm wondering what caused this, since the contrasts estimates showed reasonable values...

rvlenth commented 2 months ago

Same answer. You averaged over cut. It creates linear dependence.

Generalized commented 2 months ago

Right, system is computationally singular: reciprocal condition number = 1.65641e-18. Thank you very much.