rglmlme4sjplotgamlss

R how to use sjPlot::tab_model() to put lmer, glmer, and gamlss models into a table

I am looking to make a summary table for a set of linear models. The models are either lmer(), glmer(), or gamlss(). I am trying to print the results of 6 of these models into a single table. However, when I try to do this with sjPlot::tab_model, I receive the following error message, "Error: $ operator is invalid for atomic vectors".

I removed the two gamlss() models to see if the table would print and it did. It seems to me that the issue is that sjPlot::tab_model cannot work with gamlss() models.

Databases

recruitment_data <- structure(list(Site_long = 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, 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, 2L, 2L, 2L, 2L), .Label = c("Hanauma Bay", 
"Waikiki"), class = "factor"), Shelter = 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, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 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), .Label = c("High", 
"Low"), class = "factor"), `Module #` = structure(c(7L, 7L, 7L, 
7L, 7L, 7L, 7L, 7L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 11L, 11L, 
11L, 11L, 11L, 11L, 11L, 11L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 
10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 12L, 12L, 12L, 12L, 12L, 
12L, 12L, 12L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L), .Label = c("111", "112", "113", "114", "115", 
"116", "211", "212", "213", "214", "215", "216"), class = "factor"), 
    TimeStep = c(4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 4L, 4L, 5L, 
    5L, 6L, 6L, 7L, 7L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 4L, 4L, 
    5L, 5L, 6L, 6L, 7L, 7L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 4L, 
    4L, 5L, 5L, 6L, 6L, 7L, 7L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 
    4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 
    7L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 4L, 4L, 5L, 5L, 6L, 6L, 
    7L, 7L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L), Side = structure(c(1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L), .Label = c("N", "S"), class = "factor"), 
    recruits = c(0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 2, 2, 1, 
    0, 0, 3, 0, 0, 1, 1, 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 2, 
    3, 0, 1, 1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1, 0, 0, 1, 0, 
    0, 2, 2, 0, 1, 1, 0, 1, 3, 1, 0, 1, 1, 0, 1, 5, 2, 0, 1, 
    1, 2, 0, 3, 1, 2, 2, 3, 6, 5, 2, 0, 1, 2, 0, 4, 1, 4, 1, 
    0, 0, 4, 0, 1)), row.names = c(NA, -96L), class = c("grouped_df", 
"tbl_df", "tbl", "data.frame"), vars = c("Site_long", "Shelter", 
"Module #", "TimeStep"), drop = TRUE)


survival_data <- structure(list(Date = structure(c(17288, 17288, 17288, 17288, 
17288, 17288, 17292, 17299, 17299, 17304, 17306, 17386, 17386, 
17386, 17386, 17386, 17386, 17387, 17387, 17387, 17389, 17389, 
17389, 17390, 17398, 17404, 17475, 17475, 17477, 17480, 17482, 
17484, 17484, 17484, 17484, 17484, 17484, 17489, 17575, 17575, 
17575, 17575, 17575, 17582, 17586, 17594, 17600, 17601, 17603, 
17603), class = "Date"), Year = 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, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("17", 
"18"), class = "factor"), Site = structure(c(1L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 
1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L), .Label = c("HAN", 
"WAI"), class = "factor"), Treatment = c("CLO", "CLO", "CLO", 
"OPE", "OPE", "OPE", "CLO", "CLO", "OPE", "OPE", "CLO", "CLO", 
"CLO", "CLO", "OPE", "OPE", "OPE", "CLO", "OPE", "OPE", "CLO", 
"CLO", "OPE", "OPE", "CLO", "CLO", "CLO", "OPE", "CLO", "OPE", 
"CLO", "CLO", "CLO", "CLO", "OPE", "OPE", "OPE", "OPE", "CLO", 
"CLO", "CLO", "OPE", "OPE", "OPE", "CLO", "OPE", "OPE", "CLO", 
"CLO", "OPE"), `Module #` = c(212L, 214L, 216L, 211L, 213L, 215L, 
116L, 114L, 115L, 113L, 112L, 212L, 214L, 216L, 211L, 213L, 215L, 
116L, 111L, 115L, 112L, 114L, 115L, 113L, 114L, 114L, 112L, 115L, 
116L, 113L, 114L, 212L, 214L, 216L, 211L, 213L, 215L, 111L, 212L, 
214L, 216L, 213L, 215L, 211L, 116L, 115L, 113L, 114L, 112L, 111L
), n.x = c(1L, 1L, 2L, 2L, 3L, 3L, 6L, 7L, 5L, 4L, 2L, 2L, 2L, 
1L, 2L, 3L, 5L, 10L, 1L, 4L, 10L, 13L, 6L, 5L, 2L, 2L, 8L, 6L, 
10L, 8L, 12L, 2L, 6L, 2L, 2L, 3L, 5L, 2L, 2L, 5L, 1L, 8L, 9L, 
2L, 10L, 15L, 10L, 16L, 12L, 4L), n.y = c(1, 1, 0, 2, 3, 3, 6, 
7, 4, 4, 2, 2, 2, 1, 2, 3, 5, 9, 1, 3, 10, 11, 5, 5, 1, 1, 7, 
6, 7, 8, 11, 2, 5, 2, 2, 3, 5, 2, 2, 5, 1, 7, 7, 0, 8, 14, 9, 
9, 9, 4), `%_Survival` = c(100, 100, 0, 100, 100, 100, 100, 100, 
80, 100, 100, 100, 100, 100, 100, 100, 100, 90, 100, 75, 100, 
84.6153846153846, 83.3333333333333, 100, 50, 50, 87.5, 100, 70, 
100, 91.6666666666667, 100, 83.3333333333333, 100, 100, 100, 
100, 100, 100, 100, 100, 87.5, 77.7777777777778, 0, 80, 93.3333333333333, 
90, 56.25, 75, 100), Quarter = c(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, 3, 3, 3, 3, 3, 
3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4), Mortality = c(0, 
0, 2, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 2, 
1, 0, 1, 1, 1, 0, 3, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 2, 
2, 2, 1, 1, 7, 3, 0), Survival = c(100, 100, 0, 100, 100, 100, 
100, 100, 80, 100, 100, 100, 100, 100, 100, 100, 100, 90, 100, 
75, 100, 84.6153846153846, 83.3333333333333, 100, 50, 50, 87.5, 
100, 70, 100, 91.6666666666667, 100, 83.3333333333333, 100, 100, 
100, 100, 100, 100, 100, 100, 87.5, 77.7777777777778, 0, 80, 
93.3333333333333, 90, 56.25, 75, 100), Module = c(212L, 214L, 
216L, 211L, 213L, 215L, 116L, 114L, 115L, 113L, 112L, 212L, 214L, 
216L, 211L, 213L, 215L, 116L, 111L, 115L, 112L, 114L, 115L, 113L, 
114L, 114L, 112L, 115L, 116L, 113L, 114L, 212L, 214L, 216L, 211L, 
213L, 215L, 111L, 212L, 214L, 216L, 213L, 215L, 211L, 116L, 115L, 
113L, 114L, 112L, 111L), Survival_prop = c(1, 1, 0, 1, 1, 1, 
1, 1, 0.8, 1, 1, 1, 1, 1, 1, 1, 1, 0.9, 1, 0.75, 1, 0.846153846153846, 
0.833333333333333, 1, 0.5, 0.5, 0.875, 1, 0.7, 1, 0.916666666666667, 
1, 0.833333333333333, 1, 1, 1, 1, 1, 1, 1, 1, 0.875, 0.777777777777778, 
0, 0.8, 0.933333333333333, 0.9, 0.5625, 0.75, 1), Date_new = structure(c(17378, 
17378, 17378, 17378, 17378, 17378, 17382, 17389, 17389, 17394, 
17396, 17476, 17476, 17476, 17476, 17476, 17476, 17477, 17477, 
17477, 17479, 17479, 17479, 17480, 17488, 17494, 17565, 17565, 
17567, 17570, 17572, 17574, 17574, 17574, 17574, 17574, 17574, 
17579, 17665, 17665, 17665, 17665, 17665, 17672, 17676, 17684, 
17690, 17691, 17693, 17693), class = "Date"), Site_long = structure(c(1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 
1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
2L), .Label = c("Hanauma Bay", "Waikiki"), class = "factor"), 
    Treatment_long = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 1L, 
    1L, 2L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 
    2L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 
    2L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 2L), .Label = c("Closed", 
    "Open"), class = "factor"), Shelter = c("Low", "Low", "Low", 
    "High", "High", "High", "Low", "Low", "High", "High", "Low", 
    "Low", "Low", "Low", "High", "High", "High", "Low", "High", 
    "High", "Low", "Low", "High", "High", "Low", "Low", "Low", 
    "High", "Low", "High", "Low", "Low", "Low", "Low", "High", 
    "High", "High", "High", "Low", "Low", "Low", "High", "High", 
    "High", "Low", "High", "High", "Low", "Low", "High")), row.names = c(NA, 
-50L), vars = c("Date", "Year", "Site", "Treatment", "Module #"
), labels = structure(list(Date = structure(c(17288, 17288, 17288, 
17288, 17288, 17288, 17292, 17299, 17299, 17304, 17306, 17386, 
17386, 17386, 17386, 17386, 17386, 17387, 17387, 17387, 17389, 
17389, 17389, 17390, 17398, 17404, 17475, 17475, 17477, 17480, 
17482, 17484, 17484, 17484, 17484, 17484, 17484, 17489, 17575, 
17575, 17575, 17575, 17575, 17582, 17586, 17594, 17600, 17601, 
17603, 17603), class = "Date"), Year = c(17, 17, 17, 17, 17, 
17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 
17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 
17, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18), Site = structure(c(1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 
1L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
2L), .Label = c("HAN", "WAI"), class = "factor"), Treatment = c("CLO", 
"CLO", "CLO", "OPE", "OPE", "OPE", "CLO", "CLO", "OPE", "OPE", 
"CLO", "CLO", "CLO", "CLO", "OPE", "OPE", "OPE", "CLO", "OPE", 
"OPE", "CLO", "CLO", "OPE", "OPE", "CLO", "CLO", "CLO", "OPE", 
"CLO", "OPE", "CLO", "CLO", "CLO", "CLO", "OPE", "OPE", "OPE", 
"OPE", "CLO", "CLO", "CLO", "OPE", "OPE", "OPE", "CLO", "OPE", 
"OPE", "CLO", "CLO", "OPE"), `Module #` = c(212L, 214L, 216L, 
211L, 213L, 215L, 116L, 114L, 115L, 113L, 112L, 212L, 214L, 216L, 
211L, 213L, 215L, 116L, 111L, 115L, 112L, 114L, 115L, 113L, 114L, 
114L, 112L, 115L, 116L, 113L, 114L, 212L, 214L, 216L, 211L, 213L, 
215L, 111L, 212L, 214L, 216L, 213L, 215L, 211L, 116L, 115L, 113L, 
114L, 112L, 111L)), row.names = c(NA, -50L), class = "data.frame", vars = c("Date", 
"Year", "Site", "Treatment", "Module #"), drop = TRUE), indices = list(
    0L, 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), drop = TRUE, group_sizes = 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), biggest_group_size = 1L, class = c("grouped_df", "tbl_df", 
"tbl", "data.frame"))

Analyses

library(lme4)
library(gamlss)

recruitment_glmer_n3 <- glmer(recruits ~ Site_long*Shelter + (1|module_recruit), data = n3, family = poisson, na.action = "na.fail")
summary(recruitment_glmer_n3)

survival_gamlss <- gamlss(Survival_prop ~ Site_long*Treatment_long + (1|module_survival), data = survival_results_long_2, family = BEINF())
summary(survival_gamlss)

Table Code

library(sjPlot)

tab_model(recruitment_glmer_n3, survival_gamlss)

Is there a way to use sjPlot::tab_model to get html table outputs with gamlss model objects or is there another package that you would recommend for making publication-quality tables for linear model objects of gamlss as well as lmer and glmer? Thank you!

Solution

  • I have added (better) support for gamlss-models, so now you should be able to plot both models in one table with tab_model() (at least it worked for me). However, you need to update sjPlot and insight from GitHub first.

    Then, it should work.