0

I am trying to execute linear-mixed models with 3 factors as groups. However running the syntax, I just get 2 terms: time:grup_int2 or time:grup_int3 as you can see in the picture. There should be another term for time:grup_int1. My reference factor group is the level 1, however I don't know how to relevel in the loop, and I intended to post a example with a database provided by the package. The problem is with the database provided by the package works... and I get 3 terms

enter image description here

My code with genes_wide dataframe (plenty of NAs)


#pivot wider to make horizontal data

genes_wide <- genes1 %>% pivot_wider(names_from = gen, values_from = dCt)
genes_wide$time [genes_wide$time == "1"] <- 1 
genes_wide$time [genes_wide$time == "3"] <- 3 
genes_wide$time<-as.numeric(genes_wide$time)

# if they need to be relevel
#genes_wide$grup_int <- relevel(genes_wide$grup_int, "2") 



# names of variables
genes.names <- colnames(genes_wide)[13:67]
no.genes <- length(genes.names)

# create a named list to hold the fitted models
gene.list <- as.list(1:no.genes)
names(gene.list) <- genes.names

# loop over gene names

for(i in genes.names){ 

# create temporary data matrix and model formula
#1st: adjusting per several variables (edad0 = numeric, grup_int = factor, sexo = factor, time = numeric, id = numeric (need to change?)
#2nd: adjusting just per group

tmp <- genes_wide[, c(i, "grup_int","time", "id", "edad0", "sexo")]
fml <- as.formula( paste( i, "~", paste(c("time", "time:grup_int", "edad0", "sexo"), collapse="+")))
  
gene.list[[i]] <- lme(fml, random= ~ time|id  , method="REML", data=tmp, control = lmeControl(opt = "optim"), na.action = na.omit)


}



# List of data.frames/tibble

result <- lapply(gene.list, function(x) tidy(x))

# Solution 2

result <- dplyr::bind_rows(result, .id = "var")
result2   <- result %>%  dplyr::select( var, term, p.value)

**

Example database

** However running the exact syntax with an example database (+ 5 variables to iterate) from the nlme package, the result is different

#adding variables to df problem

BodyWeight$var1 <- rnorm(176, 47 ,3)
BodyWeight$var2 <- rnorm(176, 50 ,10)
BodyWeight$var3 <- rnorm(176, 68, 7)
BodyWeight$var4 <- rnorm(176, 150 ,10)
BodyWeight$var5 <- rnorm(176, 140, 7)

# names of variables
var.names <- colnames(BodyWeight)[5:9]
no.var <- length(var.names)

# create a named list to hold the fitted models
var.list <- as.list(1:no.var)
names(var.list) <- var.names

# loop over gene names

for(i in var.names){ 



  
#2nd: adjusting just per group

tmp <- BodyWeight[, c(i, "Diet","Time", "Rat", "weight")]
fml <- as.formula( paste( i, "~", paste(c("Time:Diet", "weight"), collapse="+")))
  
var.list[[i]] <- lme(fml, random= ~ Time|Rat  , method="REML", data=tmp, control = lmeControl(opt = "optim"), na.action = na.omit)


}



# List of data.frames/tibble

result <- lapply(var.list, function(x) tidy(x))

# Solution 

result <- dplyr::bind_rows(result, .id = "var")
result2   <- result %>%  dplyr::select( var, term, p.value)

As you can see I get 3 terms:

  • Time:Diet1
  • Time:Diet2
  • Time:Diet3

enter image description here

My original database is made up by more variables, but time is numeric, as well as grouping variables are factor (grup_int and Diet respectively). The original database contains NAs I was wondering if it is possible to relevel to get the interaction remaining

UPDATE I don't know how to put in words, what is happening when I get this:

enter image description here

What does it mean? time:grup_int1 is significant but time:grup_int2 not? Is time:grup_int1 the meaning for statistical differences between group 1 ant the others groups (2 and 3) across time ? It doesn't make any sense to be significant for one but not for the others. How can I get an overall p-value comparison between groups?

The df used for the loop

structure(list(weight = c(240, 250, 255, 260, 262, 258, 266, 
266, 265, 272, 278, 225, 230, 230, 232, 240, 240, 243, 244, 238, 
247, 245, 245, 250, 250, 255, 262, 265, 267, 267, 264, 268, 269, 
260, 255, 255, 265, 265, 268, 270, 272, 274, 273, 275, 255, 260, 
255, 270, 270, 273, 274, 273, 276, 278, 280, 260, 265, 270, 275, 
275, 277, 278, 278, 284, 279, 281, 275, 275, 260, 270, 273, 274, 
276, 271, 282, 281, 284, 245, 255, 260, 268, 270, 265, 265, 267, 
273, 274, 278, 410, 415, 425, 428, 438, 443, 442, 446, 456, 468, 
478, 405, 420, 430, 440, 448, 460, 458, 464, 475, 484, 496, 445, 
445, 450, 452, 455, 455, 451, 450, 462, 466, 472, 555, 560, 565, 
580, 590, 597, 595, 595, 612, 618, 628, 470, 465, 475, 485, 487, 
493, 493, 504, 507, 518, 525, 535, 525, 530, 533, 535, 540, 525, 
530, 543, 544, 559, 520, 525, 530, 540, 543, 546, 538, 544, 553, 
555, 548, 510, 510, 520, 515, 530, 538, 535, 542, 550, 553, 569
), Time = c(1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64), Rat = structure(c(4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 6L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 
7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 5L, 5L, 5L, 5L, 
5L, 5L, 5L, 5L, 5L, 5L, 5L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 
10L, 10L, 10L, 10L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 
11L, 11L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 12L, 12L, 
12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 13L, 13L, 13L, 13L, 
13L, 13L, 13L, 13L, 13L, 13L, 13L, 15L, 15L, 15L, 15L, 15L, 15L, 
15L, 15L, 15L, 15L, 15L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 
14L, 14L, 14L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 
16L), levels = c("2", "3", "4", "1", "8", "5", "6", "7", "11", 
"9", "10", "12", "13", "15", "14", "16"), class = c("ordered", 
"factor")), Diet = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 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, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 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, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), levels = c("1", "2", "3"), class = "factor"), 
    var1 = c(46.804, 45.792, 43.572, 50.587, 52.268, 45.402, 
    45.996, 48.566, 47.933, 48.336, 43.964, 52.441, 45.077, 45.699, 
    51.348, 47.888, 50.245, 48.111, 42.373, 44.231, 46.461, 45.639, 
    44.993, 50.391, 47.997, 48.368, 44.939, 47.805, 49.195, 51.419, 
    50.242, 45.856, 48.115, 48.053, 44.318, 50.493, 40.48, 47.418, 
    49.441, 49.414, 49.776, 49.043, 41.456, 48.52, 46.079, 43.058, 
    43.649, 43.967, 44.799, 48.081, 50.728, 51.207, 50.01, 47.346, 
    43.816, 52.575, 47.962, 49.061, 45.438, 43.822, 43.022, 48.864, 
    41.669, 45.58, 42.206, 42.788, 43.988, 46.636, 45.651, 42.502, 
    42.93, 49.196, 48.314, 46.412, 49.436, 47.193, 48.806, 45.142, 
    48.4, 48.419, 47.937, 42.605, 46.027, 45.325, 43.811, 45.922, 
    47.262, 46.566, 47.016, 46.111, 50.573, 47.436, 47.859, 47.839, 
    48.239, 46.762, 46.581, 55.9, 42.816, 49.636, 47.514, 39.489, 
    46.513, 49.323, 51.504, 45.558, 42.402, 45.153, 48.632, 47.908, 
    46.033, 44.963, 49.312, 45.733, 51.149, 45.755, 47.247, 46.674, 
    49.121, 43.411, 43.735, 44.777, 44.942, 48.207, 50.962, 47.1, 
    44.167, 50.007, 48.531, 46.749, 47.69, 47.396, 45.074, 49.203, 
    42.869, 47.713, 46.186, 52.191, 46.809, 47.578, 45.365, 53.228, 
    45.896, 45.428, 42.702, 47.328, 46.108, 49.857, 45.128, 43.772, 
    46.479, 51.573, 51.255, 45.33, 45.651, 48.521, 49.006, 51.906, 
    48.924, 53.693, 48.555, 48.367, 49.659, 47.839, 47.829, 48.212, 
    45.867, 51.965, 48.275, 49.091, 44.54, 44.978, 50.431, 48.915, 
    51.569, 46.008), var2 = c(53.273, 35.7, 44.406, 62.644, 42.971, 
    49.393, 52.585, 30.572, 52.156, 58.566, 44.631, 50.967, 48.583, 
    60.878, 50.712, 48.741, 47.584, 54.774, 62.288, 49.952, 33.225, 
    39.251, 53.158, 46.74, 46.796, 40.343, 58.491, 48.12, 41.854, 
    55.315, 44.553, 55.674, 47.404, 56.085, 53.443, 56.382, 74.641, 
    44.822, 73.191, 46.409, 66.453, 60.929, 51.394, 56.812, 53.244, 
    48.993, 54.763, 35.624, 63.346, 38.389, 48.288, 61.422, 45.707, 
    50.213, 42.629, 64.381, 39.218, 58.923, 61.569, 39.279, 54.748, 
    63.149, 54.639, 56.494, 66.907, 36.303, 52.297, 46.288, 62.852, 
    53.063, 43.096, 70.983, 49.175, 43.663, 44.456, 55.17, 51.383, 
    56.945, 51.471, 49.053, 63.062, 46.386, 48.094, 43.282, 35.509, 
    48.079, 48.223, 34.605, 42.669, 40.146, 53.089, 40.911, 39.396, 
    39.507, 41.481, 41.729, 28.367, 62.781, 59.154, 50.788, 84.282, 
    58.214, 53.404, 46.312, 36.983, 61.712, 53.109, 37.285, 50.834, 
    32.306, 47.969, 36.232, 46.744, 37.876, 50.369, 59.684, 40.427, 
    44.2, 46.955, 56.388, 46.667, 46.373, 49.932, 60.3, 27.732, 
    56.772, 54.569, 45.164, 37.747, 57.827, 49.998, 23.109, 53.339, 
    51.024, 58.875, 55.315, 40.637, 37.781, 47.209, 57.51, 46.956, 
    66.715, 51.034, 52.54, 60.342, 34.439, 58.08, 64.298, 56.374, 
    64.829, 49.392, 48.939, 38.746, 40.048, 46.856, 61.697, 51.253, 
    54.024, 61.166, 45.834, 43.578, 39.212, 67.632, 35.208, 23.676, 
    50.577, 46.832, 42.866, 57.808, 53.292, 42.474, 62.779, 34.747, 
    55.367, 27.928, 57.481), var3 = c(67.483, 67.977, 70.749, 
    76.851, 77.351, 59.821, 67.203, 72.635, 67.767, 67.627, 64.729, 
    70.388, 71.889, 61.95, 54.182, 71.049, 56.33, 68.982, 70.615, 
    68.416, 60.661, 73.205, 68.086, 72.611, 62.082, 74.026, 69.073, 
    67.488, 79.656, 79.612, 71.016, 61.284, 76.297, 62.713, 79.335, 
    68.155, 67.184, 73.889, 61.59, 68.273, 57.056, 66.95, 60.504, 
    77.121, 73.056, 60.739, 75.891, 74.363, 82.829, 86.902, 59.38, 
    77.63, 67.613, 72.666, 68.903, 71.631, 68.86, 57.653, 65.119, 
    60.825, 63.301, 64.945, 65.89, 80.092, 79.951, 65.475, 66.491, 
    72.355, 71.909, 70.339, 83.508, 59.686, 62.102, 67.291, 52.908, 
    70.508, 63.2, 66.62, 75.138, 72.311, 66.394, 64.812, 74.063, 
    51.087, 55.095, 56.906, 55.88, 71.916, 81.805, 67.339, 66.459, 
    70.005, 82.891, 65.456, 86.611, 65.401, 72.647, 59.198, 73.558, 
    66.536, 69.83, 77.101, 65.494, 62.129, 71.86, 85.406, 70.619, 
    79.225, 71.806, 62.86, 79.615, 65.049, 63.472, 63.46, 80.215, 
    74.033, 58.949, 69.519, 66.39, 68.751, 69.356, 75.201, 54.815, 
    50.592, 76.019, 62.177, 66.677, 69.495, 82.428, 70.743, 59.458, 
    75.326, 69.615, 67.157, 69.793, 66.205, 75.328, 76.311, 64.36, 
    62.207, 74.16, 51.799, 66.378, 67.964, 66.937, 59.779, 59.735, 
    64.992, 66.817, 60.595, 65.388, 50.404, 61.742, 65.789, 77.047, 
    56.86, 69.756, 66.338, 57.149, 63.407, 73.601, 66.042, 67.421, 
    51.208, 47.265, 69.465, 63.471, 75.845, 74.938, 69.009, 81.807, 
    81.195, 64.288, 77.085, 64.698, 65.747), var4 = c(143.048, 
    139.381, 152.849, 163.726, 155.712, 146.957, 140.231, 148.321, 
    141.917, 156.873, 165.845, 158.824, 157.305, 144.521, 145.253, 
    134.892, 153.597, 124.482, 153.096, 134.954, 171.238, 147.43, 
    150.983, 129.945, 150.652, 142.259, 156.631, 152.215, 152.196, 
    153.54, 139.712, 168.495, 152.043, 148.538, 164.067, 172.445, 
    144.175, 144.015, 149.093, 135.921, 156.787, 148.292, 144.958, 
    137.332, 148.185, 151.678, 160.774, 168.479, 175.839, 148.211, 
    163.747, 145.416, 131.502, 154.997, 152.35, 125.955, 167.116, 
    141.734, 135.958, 162.752, 139.482, 149.315, 144.425, 144.654, 
    160.05, 153.767, 148.532, 143.102, 135.553, 143.305, 153.021, 
    127.982, 165.047, 144.142, 136.568, 145.793, 143.913, 158.002, 
    154.043, 167.63, 165.472, 143.681, 166.417, 152.012, 145.843, 
    144.448, 138.782, 141.475, 150.717, 140.765, 164.657, 133.655, 
    151.908, 160.084, 138.25, 155.93, 142.756, 170.935, 128.346, 
    149.986, 154.099, 141.688, 150.616, 149.999, 139.195, 160.921, 
    148.261, 144.611, 154.183, 157.2, 130.512, 155.004, 125.064, 
    124.813, 150.074, 145.615, 152.238, 144.241, 149.079, 159.354, 
    153.398, 146.609, 155.032, 145.522, 143.296, 139.502, 148.944, 
    143.446, 155.094, 148.208, 154.417, 153.249, 154.544, 140.746, 
    156.798, 146.069, 134.905, 145.318, 137.198, 142.087, 137.117, 
    151.782, 147.245, 148.583, 152.535, 151.265, 152.272, 162.014, 
    150.055, 154.343, 139.717, 148.752, 140.933, 146.761, 145.311, 
    146.361, 155.264, 154.984, 168.623, 129.31, 147.492, 139.828, 
    161.22, 145.89, 139.146, 138.104, 148.543, 155.836, 143.594, 
    159.32, 163.648, 152.058, 160.412, 133, 155.305, 144.094), 
    var5 = c(135.794, 136.92, 131.548, 132.091, 143.587, 130.393, 
    133.678, 143.005, 136.036, 128.161, 134.036, 144.237, 125.767, 
    133.32, 137.047, 140.3, 137.92, 144.883, 140.901, 140.632, 
    136.315, 138.138, 140.855, 124.721, 143.153, 144.514, 142.092, 
    144.765, 139.205, 134.905, 141.815, 122.992, 138.249, 125.918, 
    125.236, 136.503, 137.462, 130.757, 135.141, 135.413, 134.991, 
    137.143, 153.329, 148.779, 138.249, 133.978, 132.454, 139.955, 
    122.84, 135.984, 140.941, 134.203, 146.091, 146.766, 136.259, 
    151.474, 139.683, 132.913, 142.39, 134.625, 124.717, 145.74, 
    135.647, 135.119, 144.182, 133.961, 130.387, 141.314, 129.26, 
    127.154, 154.74, 142.327, 142.23, 139.433, 148.245, 132.123, 
    127.374, 134.518, 133.422, 143.66, 129.306, 138.31, 147.787, 
    137.415, 138.686, 143.896, 140.973, 128.592, 139.53, 146.29, 
    145, 129.926, 132.338, 155.193, 132.286, 140.86, 135.343, 
    137.037, 140.147, 131.872, 137.858, 142.137, 131.587, 121.455, 
    142.264, 138.533, 138.769, 136.296, 141.055, 133.475, 148.136, 
    136.175, 143.947, 152.026, 137.649, 142.466, 144.819, 153.745, 
    138.95, 138.146, 133.404, 147.928, 145.928, 146.421, 145.156, 
    134.289, 127.955, 148.223, 129.987, 139.295, 141.035, 150.867, 
    136.412, 139.705, 142.013, 129.377, 143.521, 145.595, 142.936, 
    137.703, 133.344, 144.46, 143.001, 137.017, 149.985, 144.242, 
    158.029, 151.835, 133.262, 123.893, 136.338, 143.951, 136.71, 
    136.007, 147.102, 134.472, 128.058, 145.521, 142.379, 135.869, 
    144.135, 141.297, 149.133, 134.008, 138.566, 117.773, 151.701, 
    140.288, 144.17, 140.886, 151.585, 136.438, 150.399, 140.535, 
    142.601, 136.449)), row.names = c(NA, 176L), class = c("nfnGroupedData", 
"nfGroupedData", "groupedData", "data.frame"), outer = ~Diet, formula = weight ~ 
    Time | Rat, labels = list(x = "Time", y = "Body weight"), units = list(
    x = "(days)", y = "(g)"), FUN = structure(function (x) 
max(x, na.rm = TRUE), source = "function(x) max(x, na.rm = TRUE)"), order.groups = TRUE)

1 Answers1

0

The problem (i.e., difference between what you get in your original model and what you get in the other example you posted) is not because of the levelling of the factors, but because of a difference in model specification. In the first model results you posted you include Time + Time:Diet - the "main effect" of Time and the product regressors with Diet. However in the second set of models - the ones with the dummy data you provided - you do not include the main effect in the model. This is what produces the difference. I show the difference below. First, load the data.

BodyWeight <- structure(list(weight = c(240, 250, 255, 260, 262, 258, 266, 
                               266, 265, 272, 278, 225, 230, 230, 232, 240, 240, 243, 244, 238, 
                               247, 245, 245, 250, 250, 255, 262, 265, 267, 267, 264, 268, 269, 
                               260, 255, 255, 265, 265, 268, 270, 272, 274, 273, 275, 255, 260, 
                               255, 270, 270, 273, 274, 273, 276, 278, 280, 260, 265, 270, 275, 
                               275, 277, 278, 278, 284, 279, 281, 275, 275, 260, 270, 273, 274, 
                               276, 271, 282, 281, 284, 245, 255, 260, 268, 270, 265, 265, 267, 
                               273, 274, 278, 410, 415, 425, 428, 438, 443, 442, 446, 456, 468, 
                               478, 405, 420, 430, 440, 448, 460, 458, 464, 475, 484, 496, 445, 
                               445, 450, 452, 455, 455, 451, 450, 462, 466, 472, 555, 560, 565, 
                               580, 590, 597, 595, 595, 612, 618, 628, 470, 465, 475, 485, 487, 
                               493, 493, 504, 507, 518, 525, 535, 525, 530, 533, 535, 540, 525, 
                               530, 543, 544, 559, 520, 525, 530, 540, 543, 546, 538, 544, 553, 
                               555, 548, 510, 510, 520, 515, 530, 538, 535, 542, 550, 553, 569
), Time = c(1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
            22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
            50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
            22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
            50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
            22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
            50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
            22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
            50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 
            22, 29, 36, 43, 44, 50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 
            50, 57, 64, 1, 8, 15, 22, 29, 36, 43, 44, 50, 57, 64), Rat = structure(c(4L, 
                                                                                     4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                     1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
                                                                                     3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 6L, 6L, 6L, 6L, 6L, 
                                                                                     6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 
                                                                                     7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 5L, 5L, 5L, 5L, 
                                                                                     5L, 5L, 5L, 5L, 5L, 5L, 5L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 
                                                                                     10L, 10L, 10L, 10L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 
                                                                                     11L, 11L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 12L, 12L, 
                                                                                     12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 13L, 13L, 13L, 13L, 
                                                                                     13L, 13L, 13L, 13L, 13L, 13L, 13L, 15L, 15L, 15L, 15L, 15L, 15L, 
                                                                                     15L, 15L, 15L, 15L, 15L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 
                                                                                     14L, 14L, 14L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 
                                                                                     16L), levels = c("2", "3", "4", "1", "8", "5", "6", "7", "11", 
                                                                                                      "9", "10", "12", "13", "15", "14", "16"), class = c("ordered", 
                                                                                                                                                          "factor")), Diet = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                         1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                         1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                         1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
                                                                                                                                                                                         1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 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, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
                                                                                                                                                                                         2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 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, 3L, 3L, 3L, 3L, 3L, 
                                                                                                                                                                                         3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
                                                                                                                                                                                         3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), levels = c("1", "2", "3"), class = "factor"), 
var1 = c(46.804, 45.792, 43.572, 50.587, 52.268, 45.402, 
         45.996, 48.566, 47.933, 48.336, 43.964, 52.441, 45.077, 45.699, 
         51.348, 47.888, 50.245, 48.111, 42.373, 44.231, 46.461, 45.639, 
         44.993, 50.391, 47.997, 48.368, 44.939, 47.805, 49.195, 51.419, 
         50.242, 45.856, 48.115, 48.053, 44.318, 50.493, 40.48, 47.418, 
         49.441, 49.414, 49.776, 49.043, 41.456, 48.52, 46.079, 43.058, 
         43.649, 43.967, 44.799, 48.081, 50.728, 51.207, 50.01, 47.346, 
         43.816, 52.575, 47.962, 49.061, 45.438, 43.822, 43.022, 48.864, 
         41.669, 45.58, 42.206, 42.788, 43.988, 46.636, 45.651, 42.502, 
         42.93, 49.196, 48.314, 46.412, 49.436, 47.193, 48.806, 45.142, 
         48.4, 48.419, 47.937, 42.605, 46.027, 45.325, 43.811, 45.922, 
         47.262, 46.566, 47.016, 46.111, 50.573, 47.436, 47.859, 47.839, 
         48.239, 46.762, 46.581, 55.9, 42.816, 49.636, 47.514, 39.489, 
         46.513, 49.323, 51.504, 45.558, 42.402, 45.153, 48.632, 47.908, 
         46.033, 44.963, 49.312, 45.733, 51.149, 45.755, 47.247, 46.674, 
         49.121, 43.411, 43.735, 44.777, 44.942, 48.207, 50.962, 47.1, 
         44.167, 50.007, 48.531, 46.749, 47.69, 47.396, 45.074, 49.203, 
         42.869, 47.713, 46.186, 52.191, 46.809, 47.578, 45.365, 53.228, 
         45.896, 45.428, 42.702, 47.328, 46.108, 49.857, 45.128, 43.772, 
         46.479, 51.573, 51.255, 45.33, 45.651, 48.521, 49.006, 51.906, 
         48.924, 53.693, 48.555, 48.367, 49.659, 47.839, 47.829, 48.212, 
         45.867, 51.965, 48.275, 49.091, 44.54, 44.978, 50.431, 48.915, 
         51.569, 46.008), var2 = c(53.273, 35.7, 44.406, 62.644, 42.971, 
                                   49.393, 52.585, 30.572, 52.156, 58.566, 44.631, 50.967, 48.583, 
                                   60.878, 50.712, 48.741, 47.584, 54.774, 62.288, 49.952, 33.225, 
                                   39.251, 53.158, 46.74, 46.796, 40.343, 58.491, 48.12, 41.854, 
                                   55.315, 44.553, 55.674, 47.404, 56.085, 53.443, 56.382, 74.641, 
                                   44.822, 73.191, 46.409, 66.453, 60.929, 51.394, 56.812, 53.244, 
                                   48.993, 54.763, 35.624, 63.346, 38.389, 48.288, 61.422, 45.707, 
                                   50.213, 42.629, 64.381, 39.218, 58.923, 61.569, 39.279, 54.748, 
                                   63.149, 54.639, 56.494, 66.907, 36.303, 52.297, 46.288, 62.852, 
                                   53.063, 43.096, 70.983, 49.175, 43.663, 44.456, 55.17, 51.383, 
                                   56.945, 51.471, 49.053, 63.062, 46.386, 48.094, 43.282, 35.509, 
                                   48.079, 48.223, 34.605, 42.669, 40.146, 53.089, 40.911, 39.396, 
                                   39.507, 41.481, 41.729, 28.367, 62.781, 59.154, 50.788, 84.282, 
                                   58.214, 53.404, 46.312, 36.983, 61.712, 53.109, 37.285, 50.834, 
                                   32.306, 47.969, 36.232, 46.744, 37.876, 50.369, 59.684, 40.427, 
                                   44.2, 46.955, 56.388, 46.667, 46.373, 49.932, 60.3, 27.732, 
                                   56.772, 54.569, 45.164, 37.747, 57.827, 49.998, 23.109, 53.339, 
                                   51.024, 58.875, 55.315, 40.637, 37.781, 47.209, 57.51, 46.956, 
                                   66.715, 51.034, 52.54, 60.342, 34.439, 58.08, 64.298, 56.374, 
                                   64.829, 49.392, 48.939, 38.746, 40.048, 46.856, 61.697, 51.253, 
                                   54.024, 61.166, 45.834, 43.578, 39.212, 67.632, 35.208, 23.676, 
                                   50.577, 46.832, 42.866, 57.808, 53.292, 42.474, 62.779, 34.747, 
                                   55.367, 27.928, 57.481), var3 = c(67.483, 67.977, 70.749, 
                                                                     76.851, 77.351, 59.821, 67.203, 72.635, 67.767, 67.627, 64.729, 
                                                                     70.388, 71.889, 61.95, 54.182, 71.049, 56.33, 68.982, 70.615, 
                                                                     68.416, 60.661, 73.205, 68.086, 72.611, 62.082, 74.026, 69.073, 
                                                                     67.488, 79.656, 79.612, 71.016, 61.284, 76.297, 62.713, 79.335, 
                                                                     68.155, 67.184, 73.889, 61.59, 68.273, 57.056, 66.95, 60.504, 
                                                                     77.121, 73.056, 60.739, 75.891, 74.363, 82.829, 86.902, 59.38, 
                                                                     77.63, 67.613, 72.666, 68.903, 71.631, 68.86, 57.653, 65.119, 
                                                                     60.825, 63.301, 64.945, 65.89, 80.092, 79.951, 65.475, 66.491, 
                                                                     72.355, 71.909, 70.339, 83.508, 59.686, 62.102, 67.291, 52.908, 
                                                                     70.508, 63.2, 66.62, 75.138, 72.311, 66.394, 64.812, 74.063, 
                                                                     51.087, 55.095, 56.906, 55.88, 71.916, 81.805, 67.339, 66.459, 
                                                                     70.005, 82.891, 65.456, 86.611, 65.401, 72.647, 59.198, 73.558, 
                                                                     66.536, 69.83, 77.101, 65.494, 62.129, 71.86, 85.406, 70.619, 
                                                                     79.225, 71.806, 62.86, 79.615, 65.049, 63.472, 63.46, 80.215, 
                                                                     74.033, 58.949, 69.519, 66.39, 68.751, 69.356, 75.201, 54.815, 
                                                                     50.592, 76.019, 62.177, 66.677, 69.495, 82.428, 70.743, 59.458, 
                                                                     75.326, 69.615, 67.157, 69.793, 66.205, 75.328, 76.311, 64.36, 
                                                                     62.207, 74.16, 51.799, 66.378, 67.964, 66.937, 59.779, 59.735, 
                                                                     64.992, 66.817, 60.595, 65.388, 50.404, 61.742, 65.789, 77.047, 
                                                                     56.86, 69.756, 66.338, 57.149, 63.407, 73.601, 66.042, 67.421, 
                                                                     51.208, 47.265, 69.465, 63.471, 75.845, 74.938, 69.009, 81.807, 
                                                                     81.195, 64.288, 77.085, 64.698, 65.747), var4 = c(143.048, 
                                                                                                                       139.381, 152.849, 163.726, 155.712, 146.957, 140.231, 148.321, 
                                                                                                                       141.917, 156.873, 165.845, 158.824, 157.305, 144.521, 145.253, 
                                                                                                                       134.892, 153.597, 124.482, 153.096, 134.954, 171.238, 147.43, 
                                                                                                                       150.983, 129.945, 150.652, 142.259, 156.631, 152.215, 152.196, 
                                                                                                                       153.54, 139.712, 168.495, 152.043, 148.538, 164.067, 172.445, 
                                                                                                                       144.175, 144.015, 149.093, 135.921, 156.787, 148.292, 144.958, 
                                                                                                                       137.332, 148.185, 151.678, 160.774, 168.479, 175.839, 148.211, 
                                                                                                                       163.747, 145.416, 131.502, 154.997, 152.35, 125.955, 167.116, 
                                                                                                                       141.734, 135.958, 162.752, 139.482, 149.315, 144.425, 144.654, 
                                                                                                                       160.05, 153.767, 148.532, 143.102, 135.553, 143.305, 153.021, 
                                                                                                                       127.982, 165.047, 144.142, 136.568, 145.793, 143.913, 158.002, 
                                                                                                                       154.043, 167.63, 165.472, 143.681, 166.417, 152.012, 145.843, 
                                                                                                                       144.448, 138.782, 141.475, 150.717, 140.765, 164.657, 133.655, 
                                                                                                                       151.908, 160.084, 138.25, 155.93, 142.756, 170.935, 128.346, 
                                                                                                                       149.986, 154.099, 141.688, 150.616, 149.999, 139.195, 160.921, 
                                                                                                                       148.261, 144.611, 154.183, 157.2, 130.512, 155.004, 125.064, 
                                                                                                                       124.813, 150.074, 145.615, 152.238, 144.241, 149.079, 159.354, 
                                                                                                                       153.398, 146.609, 155.032, 145.522, 143.296, 139.502, 148.944, 
                                                                                                                       143.446, 155.094, 148.208, 154.417, 153.249, 154.544, 140.746, 
                                                                                                                       156.798, 146.069, 134.905, 145.318, 137.198, 142.087, 137.117, 
                                                                                                                       151.782, 147.245, 148.583, 152.535, 151.265, 152.272, 162.014, 
                                                                                                                       150.055, 154.343, 139.717, 148.752, 140.933, 146.761, 145.311, 
                                                                                                                       146.361, 155.264, 154.984, 168.623, 129.31, 147.492, 139.828, 
                                                                                                                       161.22, 145.89, 139.146, 138.104, 148.543, 155.836, 143.594, 
                                                                                                                       159.32, 163.648, 152.058, 160.412, 133, 155.305, 144.094), 
var5 = c(135.794, 136.92, 131.548, 132.091, 143.587, 130.393, 
         133.678, 143.005, 136.036, 128.161, 134.036, 144.237, 125.767, 
         133.32, 137.047, 140.3, 137.92, 144.883, 140.901, 140.632, 
         136.315, 138.138, 140.855, 124.721, 143.153, 144.514, 142.092, 
         144.765, 139.205, 134.905, 141.815, 122.992, 138.249, 125.918, 
         125.236, 136.503, 137.462, 130.757, 135.141, 135.413, 134.991, 
         137.143, 153.329, 148.779, 138.249, 133.978, 132.454, 139.955, 
         122.84, 135.984, 140.941, 134.203, 146.091, 146.766, 136.259, 
         151.474, 139.683, 132.913, 142.39, 134.625, 124.717, 145.74, 
         135.647, 135.119, 144.182, 133.961, 130.387, 141.314, 129.26, 
         127.154, 154.74, 142.327, 142.23, 139.433, 148.245, 132.123, 
         127.374, 134.518, 133.422, 143.66, 129.306, 138.31, 147.787, 
         137.415, 138.686, 143.896, 140.973, 128.592, 139.53, 146.29, 
         145, 129.926, 132.338, 155.193, 132.286, 140.86, 135.343, 
         137.037, 140.147, 131.872, 137.858, 142.137, 131.587, 121.455, 
         142.264, 138.533, 138.769, 136.296, 141.055, 133.475, 148.136, 
         136.175, 143.947, 152.026, 137.649, 142.466, 144.819, 153.745, 
         138.95, 138.146, 133.404, 147.928, 145.928, 146.421, 145.156, 
         134.289, 127.955, 148.223, 129.987, 139.295, 141.035, 150.867, 
         136.412, 139.705, 142.013, 129.377, 143.521, 145.595, 142.936, 
         137.703, 133.344, 144.46, 143.001, 137.017, 149.985, 144.242, 
         158.029, 151.835, 133.262, 123.893, 136.338, 143.951, 136.71, 
         136.007, 147.102, 134.472, 128.058, 145.521, 142.379, 135.869, 
         144.135, 141.297, 149.133, 134.008, 138.566, 117.773, 151.701, 
         140.288, 144.17, 140.886, 151.585, 136.438, 150.399, 140.535, 
         142.601, 136.449)), row.names = c(NA, 176L), class = c("nfnGroupedData", 
                                                                "nfGroupedData", "groupedData", "data.frame"), outer = ~Diet, formula = weight ~ 
  Time | Rat, labels = list(x = "Time", y = "Body weight"), units = list(
    x = "(days)", y = "(g)"), FUN = structure(function (x) 
      max(x, na.rm = TRUE), source = "function(x) max(x, na.rm = TRUE)"), order.groups = TRUE)

Here is the model without the main effect of Time. You can see, you get the desired result - coefficients for each value of Diet.

library(name)

mod0 <- lme(var1 ~ Time:Diet + weight, random= ~ Time|Rat  , method="REML", data=BodyWeight, control = lmeControl(opt = "optim"), na.action = na.omit)
summary(mod0)
#> Linear mixed-effects model fit by REML
#>   Data: BodyWeight 
#>        AIC      BIC    logLik
#>   905.3992 933.6741 -443.6996
#> 
#> Random effects:
#>  Formula: ~Time | Rat
#>  Structure: General positive-definite, Log-Cholesky parametrization
#>             StdDev     Corr  
#> (Intercept) 0.73280912 (Intr)
#> Time        0.02632007 -0.912
#> Residual    2.73566273       
#> 
#> Fixed effects:  var1 ~ Time:Diet + weight 
#>                Value Std.Error  DF  t-value p-value
#> (Intercept) 46.94072 1.3142431 156 35.71692  0.0000
#> weight       0.00016 0.0033845 156  0.04685  0.9627
#> Time:Diet1  -0.00470 0.0160669 156 -0.29282  0.7701
#> Time:Diet2   0.00369 0.0194774 156  0.18956  0.8499
#> Time:Diet3   0.02582 0.0213189 156  1.21119  0.2277
#>  Correlation: 
#>            (Intr) weight Tm:Dt1 Tm:Dt2
#> weight     -0.940                     
#> Time:Diet1 -0.698  0.497              
#> Time:Diet2  0.311 -0.533  0.109       
#> Time:Diet3  0.422 -0.635  0.026  0.620
#> 
#> Standardized Within-Group Residuals:
#>         Min          Q1         Med          Q3         Max 
#> -2.75835293 -0.70479591  0.02987446  0.58694850  3.08058470 
#> 
#> Number of Observations: 176
#> Number of Groups: 16

Now, let's look at the model specified with the Time main effect as well. Here, you can see that you get the Time coefficient and one for all but the reference level of Diet interacted with Time. The omission comes because Time = Time:Diet1 + Time:Diet2 + Time:Diet3, so these four terms are perfectly collinear.

mod1 <- lme(var1 ~ Time + Time:Diet + weight, random= ~ Time|Rat  , method="REML", data=BodyWeight, control = lmeControl(opt = "optim"), na.action = na.omit)
summary(mod1)
#> Linear mixed-effects model fit by REML
#>   Data: BodyWeight 
#>        AIC      BIC    logLik
#>   905.3992 933.6741 -443.6996
#> 
#> Random effects:
#>  Formula: ~Time | Rat
#>  Structure: General positive-definite, Log-Cholesky parametrization
#>             StdDev     Corr  
#> (Intercept) 0.73280912 (Intr)
#> Time        0.02632007 -0.912
#> Residual    2.73566273       
#> 
#> Fixed effects:  var1 ~ Time + Time:Diet + weight 
#>                Value Std.Error  DF  t-value p-value
#> (Intercept) 46.94072 1.3142431 156 35.71692  0.0000
#> Time        -0.00470 0.0160669 156 -0.29282  0.7701
#> weight       0.00016 0.0033845 156  0.04685  0.9627
#> Time:Diet2   0.00840 0.0238626 156  0.35188  0.7254
#> Time:Diet3   0.03053 0.0263585 156  1.15811  0.2486
#>  Correlation: 
#>            (Intr) Time   weight Tm:Dt2
#> Time       -0.698                     
#> weight     -0.940  0.497              
#> Time:Diet2  0.724 -0.585 -0.770       
#> Time:Diet3  0.767 -0.588 -0.816  0.752
#> 
#> Standardized Within-Group Residuals:
#>         Min          Q1         Med          Q3         Max 
#> -2.75835293 -0.70479591  0.02987446  0.58694850  3.08058470 
#> 
#> Number of Observations: 176
#> Number of Groups: 16

You could recover all the same information from the second model. The effect of Time for Diet1 is just the coefficient on Time. The effect of Time for Diet2 is the effect of Time plus the effect of Time:Diet2 and its standard error is the square root of the variance of the time effect, plus the variance of the interaction coefficient plus two times the covariance of the two coefficients (this comes from the equation for the calculation of the variance of a linear combination of random variables).

b <- fixef(mod1)
V <- vcov(mod1)

## Time:Diet1
d1 <- c(b[2], sqrt(V[2,2]))
d1 <- c(d1, d1[1]/d1[2])
d1 <- c(d1, 2*pt(abs(d1[3]), 156, lower.tail=FALSE))

## Time:Diet2
d2 <- c(b[2] + b[4], sqrt(V[2,2] + V[4,4] + 2*V[2,4]))
d2 <- c(d2, d2[1]/d2[2])
d2 <- c(d2, 2*pt(abs(d2[3]), 156, lower.tail=FALSE))

## Time:Diet2
d3 <- c(b[2] + b[5], sqrt(V[2,2] + V[5,5] + 2*V[2,5]))
d3 <- c(d3, d3[1]/d3[2])
d3 <- c(d3, 2*pt(abs(d3[3]), 156, lower.tail=FALSE))

d <- rbind(d1, d2, d3)
colnames(d) <- c("Estimate", "SE", "t-stat", "p-value")
rownames(d) <- c("Time:Diet1", "Time:Diet2", "Time:Diet3")
d
#>                Estimate         SE     t-stat   p-value
#> Time:Diet1 -0.004704694 0.01606695 -0.2928182 0.7700503
#> Time:Diet2  0.003692093 0.01947738  0.1895580 0.8499019
#> Time:Diet3  0.025821390 0.02131895  1.2111944 0.2276525

Compare this result with the one from the first model above you can see that they are the same:

#>             Estimate        SE  DF   t-stat   p-value
#> Time:Diet1  -0.00470 0.0160669 156 -0.29282  0.7701
#> Time:Diet2   0.00369 0.0194774 156  0.18956  0.8499
#> Time:Diet3   0.02582 0.0213189 156  1.21119  0.2277

Created on 2023-02-01 by the reprex package (v2.0.1)

So, the answer to how you get a coefficient for the interaction of each level of Diet and Time is to leave out the main effect of Time. That said, I would encourage you to include the main effect of Diet in the model as well - this would allow the intercepts to vary according to Diet.

DaveArmstrong
  • 18,377
  • 2
  • 13
  • 25
  • Yes, the two model specifications: `Time + Time:Diet` and `Time:Diet` are in essence the same. The difference will actually be in the other coefficients: `time:grup_int2` in the first specification will be the difference between the effect of timefor `grup_int1` and `grup_int2`, such that the effect of `time` for `grup_int2` is the coefficient for `time` plus the coefficient for `time:grup_int2`. In the second specification, the coefficient for `time:grup_int2` is the effect of `time` for `grup_int2`. – DaveArmstrong Feb 01 '23 at 13:21
  • My comment above notwithstanding, while `time:grup_int` and `time + time:grup_int` are the same, `time:grup_int` and `grup_int + time:grup_int` are not the same. – DaveArmstrong Feb 01 '23 at 13:23
  • My previous comment wasn't intended to be posted. The term time:grup_int2 measures the interaction time:grup_int2 with respect the other groups, I guess. I'm not pretty sure changing it by grup_int + time:grup_int, what difference I am making in the theoretical aspect of the model – Javier Hernando Feb 01 '23 at 13:45
  • I was told that maybe my model will be better adjusted with time + time:grup_int because being a clinical trial, at time 0 the intervention cannot possibly cause any effect. And the time term alone could account for that. Does it make any sense? – Javier Hernando Feb 01 '23 at 13:54
  • It should be that the two specifications `time + grup_int + time:grup_int` and `grup_int + time:grup_int` provide the same model fit. The effect of the intervention at time=0 is captured by the inclusion of the main effect of `grup_int` not by the main effect of `time`. – DaveArmstrong Feb 01 '23 at 14:27
  • I have been thinking about this a lot. Shouldn't it be the most pure way to do it: gen ~ time:group + age + sex, data = tmp, random = ~time | id, method = "REML" I am just putting how groups evolve across time. And I get the following terms time:group1 time:group2 time:group3 Isn't it time:group1 how the group1 behaves across time in respect of the other 2 groups? This is may interest, if there are significant differences. Is this approach too complex, and there is no need for mixed models and just need basic ANOVA? – Javier Hernando Mar 17 '23 at 12:55