rr-marginaleffects

Plotting Population-level predictions with marginaleffects

I was trying to replicate this blog post from Vincent. to represent the "group" effect at population level of my model.

Edit: I have now realized that I need to treat time as continuous here.

library(patchwork)
library(marginaleffects)
library(lme4)
library(ggplot2)

m1 <- lmer (y~ time_cont*group + cov + (1 | ID), data = dat)
summary(m1)



pred <- predictions(  m1,
  newdata = datagrid(ID = NA,
                     group = c("CONTROL", "INT"),
                     time_cont = 0:4),
  include_random = FALSE)




ggplot(pred, aes(x = time_cont, y = predicted, ymin = conf.low, ymax = conf.high)) +
  geom_ribbon(alpha = .1, fill = "red") +
  geom_line() +
  facet_wrap(~ group) +
  labs(title = "Population-level trajectories")

Which produces this graph: enter image description here

The only doubt remaining is if it is possible to produce a single graph with the lines and 95% CI of both groups?

Thanks in advance.

Data below

dat <- structure(list(group = c("CONTROL", "CONTROL", "CONTROL", "INT", 
"INT", "INT", "CONTROL", "CONTROL", "CONTROL", "INT", 
"INT", "INT", "CONTROL", "CONTROL", "CONTROL", "INT", 
"INT", "INT", "INT", "INT", "INT", "INT", "INT", 
"INT", "CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", 
"CONTROL", "INT", "INT", "INT", "INT", "INT", 
"INT", "INT", "INT", "INT", "INT", "INT", "INT", 
"CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", 
"CONTROL", "CONTROL", "CONTROL", "INT", "INT", "INT", 
"CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", 
"CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", 
"CONTROL", "CONTROL", "CONTROL", "INT", "INT", "INT", 
"INT", "INT", "INT", "INT", "INT", "INT", "INT", 
"INT", "INT", "INT", "INT", "INT", "CONTROL", 
"CONTROL", "CONTROL", "INT", "INT", "INT", "INT", 
"INT", "INT", "CONTROL", "CONTROL", "CONTROL", "CONTROL", 
"CONTROL", "CONTROL", "CONTROL", "CONTROL", "CONTROL", "INT", 
"INT", "INT", "CONTROL", "CONTROL", "CONTROL", "INT", 
"INT", "INT", "CONTROL", "CONTROL", "CONTROL", "INT", 
"INT", "INT"), time = c("0month", "3month", "4month", "0month", 
"3month", "4month", "0month", "3month", "4month", "0month", "3month", 
"4month", "0month", "3month", "4month", "0month", "3month", "4month", 
"0month", "3month", "4month", "0month", "3month", "4month", "0month", 
"3month", "4month", "0month", "3month", "4month", "0month", "3month", 
"4month", "0month", "3month", "4month", "0month", "3month", "4month", 
"0month", "3month", "4month", "0month", "3month", "4month", "0month", 
"3month", "4month", "0month", "3month", "4month", "0month", "3month", 
"4month", "0month", "3month", "4month", "0month", "3month", "4month", 
"0month", "3month", "4month", "0month", "3month", "4month", "0month", 
"3month", "4month", "0month", "3month", "4month", "0month", "3month", 
"4month", "0month", "3month", "4month", "0month", "3month", "4month", 
"0month", "3month", "4month", "0month", "3month", "4month", "0month", 
"3month", "4month", "0month", "3month", "4month", "0month", "3month", 
"4month", "0month", "3month", "4month", "0month", "3month", "4month", 
"0month", "3month", "4month", "0month", "3month", "4month", "0month", 
"3month", "4month", "0month", "3month", "4month", "0month", "3month", 
"4month"), time_cont = c(0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 
0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 
0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 
0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 
0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 
0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4, 0, 3, 4
), ID = c("HF_01", "HF_01", "HF_01", "HF_02", "HF_02", "HF_02", 
"HF_03", "HF_03", "HF_03", "HF_04", "HF_04", "HF_04", "HF_05", 
"HF_05", "HF_05", "HF_06", "HF_06", "HF_06", "HF_07", "HF_07", 
"HF_07", "HF_08", "HF_08", "HF_08", "HF_09", "HF_09", "HF_09", 
"HF_10", "HF_10", "HF_10", "HF_11", "HF_11", "HF_11", "HF_12", 
"HF_12", "HF_12", "HF_13", "HF_13", "HF_13", "HF_14", "HF_14", 
"HF_14", "HF_15", "HF_15", "HF_15", "HF_16", "HF_16", "HF_16", 
"HF_17", "HF_17", "HF_17", "HF_18", "HF_18", "HF_18", "HF_19", 
"HF_19", "HF_19", "HF_20", "HF_20", "HF_20", "HF_21", "HF_21", 
"HF_21", "HF_22", "HF_22", "HF_22", "HF_23", "HF_23", "HF_23", 
"HF_24", "HF_24", "HF_24", "HF_25", "HF_25", "HF_25", "HF_26", 
"HF_26", "HF_26", "HF_27", "HF_27", "HF_27", "HF_28", "HF_28", 
"HF_28", "HF_29", "HF_29", "HF_29", "HF_30", "HF_30", "HF_30", 
"HF_31", "HF_31", "HF_31", "HF_32", "HF_32", "HF_32", "HF_33", 
"HF_33", "HF_33", "HF_34", "HF_34", "HF_34", "HF_36", "HF_36", 
"HF_36", "HF_37", "HF_37", "HF_37", "HF_38", "HF_38", "HF_38", 
"HF_39", "HF_39", "HF_39", "HF_40", "HF_40", "HF_40"), y = c(18.675, 
17.85, 17.175, 19.125, 17.55, 17.25, 19.5, 17.625, 20.325, 21.825, 
19.2, 20.7, 18.825, 18.225, 17.85, 19.125, 15.975, 17.25, 16.425, 
16.35, 14.025, 19.725, 18.6, 18.375, 17.85, 16.35, 16.5, 20.175, 
18.6, NA, 17.7, 18.15, 17.4, 17.625, 16.275, 18.075, 18.375, 
19.725, 17.7, 16.05, 12.675, 14.85, 16.8, 17.25, 16.275, 16.875, 
15, 14.025, 19.8, 19.275, 20.175, 22.275, 20.625, 18.3, 18.9, 
18.9, 19.2, 19.35, 17.7, 17.85, 19.725, 18.825, 17.175, 20.4, 
NA, 18.9, 18.225, 16.8, 17.025, 19.575, 19.875, 20.625, 18, 17.4, 
17.1, 13.725, 15.9, NA, 20.25, 16.8, 18.675, 18.3, 17.1, 16.95, 
17.55, 17.625, 18.825, 19.05, 15, 15.675, 17.4, 15.675, 15.3, 
16.35, 17.55, 16.875, 15.75, 12.6, 15.15, 18.375, 18.375, 17.85, 
19.65, 18.3, 19.875, 23.025, 18, 16.875, 15.525, 16.125, 15.9, 
22.875, 18.825, 18, 16.8, 19.95, 19.05), cov = c(4.095, 3.885, 
3.63, 4.41, 4.29, 4.44, 3.435, 2.235, NA, 4.08, 3.075, NA, 6.7815, 
5.01, 4.305, 4.53, 3.885, 3.735, 4.095, NA, 3.435, 5.115, 4.41, 
3.735, 4.14, 3.63, 3.375, 4.995, NA, NA, 4.665, 5.295, 4.38, 
3.99, 4.44, 4.08, 5.94, 5.325, 6.165, 4.575, 4.44, 3.795, 4.335, 
4.215, 4.185, 3.555, NA, 3.285, 4.98, 3.825, 3.795, 4.635, 4.53, 
4.725, 4.395, 4.185, 4.065, 4.59, NA, NA, 4.2, 3.765, 4.14, 3.27, 
NA, NA, 3.345, 3.3525, 3.51, 3.66, 3.81, 3.75, 4.245, 4.215, 
3.36, 4.395, 3.945, NA, 3.69, 3.9, 3.495, 4.68, 3.255, 4.065, 
4.41, NA, 3.75, 5.43, 4.515, 4.47, 4.47, 4.14, NA, 3.39, NA, 
NA, 4.605, 4.35, 4.32, 5.16, 4.365, 3.84, 4.32, 4.62, 3.69, 5.445, 
2.985, 2.265, 4.17, 3.99, 3.39, 4.185, 3.42, 3.84, 5.235, 4.56, 
4.515)), row.names = c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 
11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 
24L, 25L, 26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L, 
37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 
50L, 51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 
63L, 64L, 65L, 66L, 67L, 68L, 69L, 70L, 71L, 72L, 73L, 74L, 75L, 
76L, 77L, 78L, 79L, 80L, 81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 
89L, 90L, 91L, 92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L, 
101L, 102L, 106L, 107L, 108L, 109L, 110L, 111L, 112L, 113L, 114L, 
115L, 116L, 117L, 118L, 119L, 120L), class = "data.frame")
Solution

  • Your request is pretty easy to fulfill; just take out the facet_wrap() and add color/fill/linetype mappings

    ggplot(pred, aes(x = time_cont, y = predicted,
                 ymin = conf.low, ymax = conf.high)) +
  geom_ribbon(alpha = .1, aes(fill = group)) +
  geom_line(aes(colour=group, linetype = group)) +
  labs(title = "Population-level trajectories")

    enter image description here