Study 1 analyses
Dani Cosme
2024-08-19
In this report, we reproduce the exploratory Study 1 analyses examining downstream effects of the interventions.
prep data
First, we load the relevant packages and data, and define the plotting aesthetics.
load packages
define functions
# MLM results table function
table_model = function(model_data, eff_size = FALSE, logistic = FALSE, lm = FALSE) {
if (lm == TRUE) {
results = model_data %>%
broom::tidy(conf.int = TRUE) %>%
rename("SE" = std.error,
"t" = statistic,
"p" = p.value) %>%
mutate_at(vars(-contains("term"), -contains("p")), round, 2) %>%
mutate(term = gsub("article_cond", "", term),
term = gsub("\\(Intercept\\)", "control", term),
term = gsub("sharing_type", "sharing type", term),
term = gsub("msg_share_narrow", " (narrow)", term),
term = gsub("msg_rel_self_between", "self-relevance between", term),
term = gsub("msg_rel_social_between", "social relevance between", term),
term = gsub("msg_rel_self_within", "self-relevance within", term),
term = gsub("msg_rel_social_within", "social relevance within", term),
term = gsub("action_current", "current behavior", term),
term = gsub("n_c", "word count", term),
term = gsub("self$", "self - control", term),
term = gsub("other$", "other - control", term),
term = gsub(":", " x ", term),
p = ifelse(p < .001, "< .001",
ifelse(p == 1, "1.000", gsub("0.(.*)", ".\\1", sprintf("%.3f", p)))),
`b [95% CI]` = sprintf("%.2f [%0.2f, %.2f]", estimate, conf.low, conf.high))
} else {
results = model_data %>%
broom.mixed::tidy(conf.int = TRUE) %>%
filter(effect == "fixed") %>%
rename("SE" = std.error,
"t" = statistic,
"p" = p.value) %>%
select(-group, -effect) %>%
mutate_at(vars(-contains("term"), -contains("p")), round, 2) %>%
mutate(term = gsub("article_cond", "", term),
term = gsub("\\(Intercept\\)", "control", term),
term = gsub("sharing_type", "sharing type", term),
term = gsub("msg_share_narrow", " (narrow)", term),
term = gsub("msg_rel_self_between", "self-relevance between", term),
term = gsub("msg_rel_social_between", "social relevance between", term),
term = gsub("msg_rel_self_within", "self-relevance within", term),
term = gsub("msg_rel_social_within", "social relevance within", term),
term = gsub("action_current", "current behavior", term),
term = gsub("n_c", "word count", term),
term = gsub("self$", "self - control", term),
term = gsub("other$", "other - control", term),
term = gsub(":", " x ", term),
p = ifelse(p < .001, "< .001",
ifelse(p == 1, "1.000", gsub("0.(.*)", ".\\1", sprintf("%.3f", p)))),
`b [95% CI]` = sprintf("%.2f [%0.2f, %.2f]", estimate, conf.low, conf.high))
}
if (eff_size == TRUE) {
eff_size = effectsize::effectsize(model_data, type = "d") %>%
data.frame() %>%
rename("term" = Parameter) %>%
mutate(term = gsub("article_cond", "", term),
term = gsub("article_cond", "", term),
term = gsub("\\(Intercept\\)", "control", term),
term = gsub("sharing_type", "sharing type", term),
term = gsub("msg_rel_self_between", "self-relevance between", term),
term = gsub("msg_rel_social_between", "social relevance between", term),
term = gsub("msg_rel_self_within", "self-relevance within", term),
term = gsub("msg_rel_social_within", "social relevance within", term),
`d [95% CI]` = sprintf("%.2f [%0.2f, %.2f]", Std_Coefficient, CI_low, CI_high)) %>%
select(term, `d [95% CI]`)
if (lm == TRUE) {
results %>%
left_join(., eff_size) %>%
select(term, `b [95% CI]`, `d [95% CI]`, t, p) %>%
kable() %>%
kableExtra::kable_styling()
} else {
results %>%
left_join(., eff_size) %>%
select(term, `b [95% CI]`, `d [95% CI]`, df, t, p) %>%
kable() %>%
kableExtra::kable_styling()
}
} else if (logistic == TRUE | lm == TRUE) {
results %>%
select(term, `b [95% CI]`, t, p) %>%
kable() %>%
kableExtra::kable_styling()
} else {
results %>%
select(term, `b [95% CI]`, df, t, p) %>%
kable() %>%
kableExtra::kable_styling()
}
}
# Run bayesian mediation model
run_brm_model = function(model_name, model_formula, y_var, data) {
if (file.exists(sprintf("models/model_%s.RDS", model_name))) {
assign(get("model_name"), readRDS(sprintf("models/model_%s.RDS", model_name)))
} else {
assign(get("model_name"),
brm(
model_formula,
data = data,
cores = 4,
thin = 4,
seed = seed,
control = list(adapt_delta = .99, max_treedepth = 15)
))
saveRDS(eval(parse(text = model_name)), sprintf("models/model_%s.RDS", model_name))
return(eval(parse(text = model_name)))
}
}
# Get path estimates from bayesian mediation models
create_paths = function(model, x_var, y_var) {
paths = posterior_samples(model) %>%
mutate(a1 = get(sprintf("b_msgrelself_article_cond%s", x_var)),
a2 = get(sprintf("b_msgrelsocial_article_cond%s", x_var)),
b1 = get(sprintf("b_%s_msg_rel_self", y_var)),
b2 = get(sprintf("b_%s_msg_rel_social", y_var)),
c_prime = get(sprintf("b_%s_article_cond%s", y_var, x_var)),
a1b1 = a1 * b1,
a2b2 = a2 * b2,
c = c_prime + a1b1 + a2b2,
cor1 = get(sprintf("cor_SID__msgrelself_Intercept__%s_msg_rel_self", y_var)),
cor2 = get(sprintf("cor_SID__msgrelsocial_Intercept__%s_msg_rel_social", y_var)),
sd_a1 = sd_SID__msgrelself_Intercept,
sd_b1 = get(sprintf("sd_SID__%s_msg_rel_self", y_var)),
sd_a2 = sd_SID__msgrelsocial_Intercept,
sd_b2 = get(sprintf("sd_SID__%s_msg_rel_social", y_var)),
cov_a1b1 = cor1*sd_a1*sd_b1,
cov_a2b2 = cor2*sd_a2*sd_b2,
a1b1_cov_a1b1 = a1b1 + cov_a1b1,
a2b2_cov_a2b2 = a2b2 + cov_a2b2,
model = x_var,
outcome = y_var)
return(paths)
}
create_paths_words = function(model, x_var, y_var) {
y_var = gsub("_", "", y_var)
paths = posterior_samples(model) %>%
mutate(a1 = get(sprintf("b_nc_article_cond%s", x_var)),
b1 = get(sprintf("b_%s_n_c", y_var)),
c_prime = get(sprintf("b_%s_article_cond%s", y_var, x_var)),
a1b1 = a1 * b1,
c = c_prime + a1b1,
cor1 = get(sprintf("cor_SID__nc_article_cond%s__%s_n_c", x_var, y_var)),
sd_a1 = get(sprintf("sd_SID__nc_article_cond%s", x_var)),
sd_b1 = get(sprintf("sd_SID__%s_n_c", y_var)),
cov_a1b1 = cor1*sd_a1*sd_b1,
a1b1_cov_a1b1 = a1b1 + cov_a1b1,
model = x_var,
outcome = y_var)
return(paths)
}
get_paths = function(model, x_var, y_var) {
create_paths(model, x_var, y_var) %>%
select(a1:a2b2_cov_a2b2, -contains("sd"), -contains("cor"), -starts_with("cov")) %>%
gather(path, value) %>%
article_cond_by(path) %>%
summarize(median = median(value),
`Mdn [95% CI]` = sprintf("%.2f [%.2f, %.2f]", median(value), quantile(value, probs = .025), quantile(value, probs = .975))) %>%
mutate(path = factor(path, levels = c("a1", "b1", "a1b1", "a1b1_cov_a1b1", "a2", "b2", "a2b2", "a2b2_cov_a2b2", "c", "c_prime")),
variable = ifelse(grepl("1", path), "self-relevance",
ifelse(grepl("2", path), "social relevance", ""))) %>%
arrange(path) %>%
select(variable, path, `Mdn [95% CI]`) %>%
kable() %>%
kableExtra::kable_styling()
}
get_paths_words = function(model, x_var, y_var) {
create_paths_words(model, x_var, y_var) %>%
select(a1:a1b1_cov_a1b1, -contains("sd"), -contains("cor"), -starts_with("cov")) %>%
gather(path, value) %>%
article_cond_by(path) %>%
summarize(median = median(value),
`Mdn [95% CI]` = sprintf("%.2f [%.2f, %.2f]", median(value), quantile(value, probs = .025), quantile(value, probs = .975))) %>%
mutate(path = factor(path, levels = c("a1", "b1", "a1b1", "a1b1_cov_a1b1", "c", "c_prime")),
variable = ifelse(grepl("1", path), "self-relevance",
ifelse(grepl("2", path), "social relevance", ""))) %>%
arrange(path) %>%
select(variable, path, `Mdn [95% CI]`) %>%
kable() %>%
kableExtra::kable_styling()
}
percent_mediated = function(model, x_var, y_var) {
create_paths(model, x_var, y_var) %>%
select(a1b1_cov_a1b1, a2b2_cov_a2b2, c) %>%
gather(path, value) %>%
article_cond_by(path) %>%
summarize(median = median(value)) %>%
select(path, median) %>%
spread(path, median) %>%
mutate(self = round((a1b1_cov_a1b1 / c) * 100, 0),
social = round((a2b2_cov_a2b2 / c) * 100, 0),
total = self + social) %>%
select(self, social, total) %>%
kable(caption = "percent mediated") %>%
kableExtra::kable_styling()
}
percent_mediated_words = function(model, x_var, y_var) {
create_paths_words(model, x_var, y_var) %>%
select(a1b1_cov_a1b1, c) %>%
gather(path, value) %>%
article_cond_by(path) %>%
summarize(median = median(value)) %>%
select(path, median) %>%
spread(path, median) %>%
mutate(word_count = round((a1b1_cov_a1b1 / c) * 100, 0)) %>%
select(word_count) %>%
kable(caption = "percent mediated") %>%
kableExtra::kable_styling()
}define aesthetics
palette_condition = c("self" = "#ee9b00",
"control" = "#0a9396",
"other" = "#005f73")
palette_dv = c("self-relevance" = "#ee9b00",
"social relevance" = "#005f73",
"broadcast sharing" = "#5F0F40",
"narrowcast sharing" = "#D295BF")
palette_sharing = c("broadcast sharing" = "#5F0F40",
"narrowcast sharing" = "#D295BF")
plot_aes = theme_minimal() +
theme(legend.position = "top",
legend.text = element_text(size = 12),
text = element_text(size = 16, family = "Futura Medium"),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.text = element_text(color = "black"),
axis.line = element_line(colour = "black"),
axis.ticks.y = element_blank())descriptives
exploratory downstream effects
petitions
Do the interventions affect climate petition engagement?
categories = unique(filter(exploratory_dvs, grepl("petition", scale_name))$scale_name)
for (category in categories){
cat(paste0('\n\n### ', category, '{.tabset}\n\n'))
data_loop = exploratory_dvs %>%
filter(scale_name == !!category)
if (category == "petition_link_clicks"){
model = lm(value ~ article_cond, data = data_loop)
print(table_model(model, lm = TRUE))
} else {
model = lmer(value ~ article_cond + (1 | SID) + (1 | item), data = data_loop)
print(table_model(model))
}
predicted = ggeffects::ggpredict(model, terms = "article_cond") %>%
data.frame()
print(
predicted %>%
mutate(x = factor(x, levels = c("self", "control", "other"))) %>%
ggplot(aes(x, predicted, color = x)) +
geom_line(aes(group = 1), size = 1, color = "black") +
geom_pointrange(aes(ymin = conf.low, ymax = conf.high), size = 1, linewidth = 1) +
scale_color_manual(values = palette_condition) +
labs(x = "", y = "predicted\n") +
plot_aes +
theme(legend.position = "none")
)
}petition_share_broad
| term | b [95% CI] | df | t | p |
|---|---|---|---|---|
| control | -0.09 [-0.17, -0.02] | 918.6 | -2.52 | .012 |
| other - control | 0.21 [0.08, 0.35] | 918.6 | 3.21 | .001 |
| self - control | 0.18 [0.05, 0.31] | 918.6 | 2.68 | .007 |
petition_share_narrow
| term | b [95% CI] | df | t | p |
|---|---|---|---|---|
| control | -0.07 [-0.14, 0.00] | 908.23 | -1.86 | .063 |
| other - control | 0.18 [0.05, 0.31] | 908.18 | 2.74 | .006 |
| self - control | 0.11 [-0.02, 0.23] | 908.18 | 1.61 | .108 |
actions
Do the interventions affect climate actions?
across categories
model = lmer(value ~ article_cond + (1 | SID),
data = filter(exploratory_dvs, grepl("action_env_impact", scale_name)))
table_model(model)| term | b [95% CI] | df | t | p |
|---|---|---|---|---|
| control | -0.06 [-0.11, -0.01] | 1610.02 | -2.27 | .023 |
| other - control | 0.17 [0.08, 0.26] | 1610.01 | 3.77 | < .001 |
| self - control | 0.07 [-0.02, 0.16] | 1610.01 | 1.60 | .110 |
by category
model = lmer(value ~ article_cond * category + (1 | SID),
data = filter(exploratory_dvs, grepl("action_env_impact", scale_name)))
cis = confint(emmeans::emmeans(model, specs = revpairwise ~ article_cond | category,
type = "response", adjust = "none", pbkrtest.limit = 19367))$contrasts %>%
data.frame()
emmeans::emmeans(model, specs = revpairwise ~ article_cond | category,
type = "response", adjust = "none", pbkrtest.limit = 19367)$contrasts %>%
data.frame() %>%
left_join(., cis) %>%
mutate(p = ifelse(p.value < .001, "< .001",
ifelse(p.value == 1, "1.000", gsub("0.(.*)", ".\\1", sprintf("%.3f", p.value)))),
`b [95% CI]` = sprintf("%.2f [%0.2f, %.2f]", estimate, lower.CL, upper.CL),
df = round(df, 2),
t.ratio = round(t.ratio, 2)) %>%
rename("t" = t.ratio) %>%
filter(!contrast == "self - other") %>%
select(category, contrast, `b [95% CI]`, df, t, p) %>%
kable() %>%
kableExtra::kable_styling()| category | contrast | b [95% CI] | df | t | p |
|---|---|---|---|---|---|
| collective | other - control | 0.15 [0.06, 0.25] | 2152.58 | 3.19 | .001 |
| collective | self - control | 0.09 [-0.01, 0.18] | 2152.58 | 1.78 | .075 |
| conversations | other - control | 0.13 [0.01, 0.25] | 5252.29 | 2.10 | .036 |
| conversations | self - control | 0.12 [-0.00, 0.24] | 5252.29 | 1.93 | .053 |
| diet | other - control | 0.19 [0.09, 0.29] | 2448.90 | 3.81 | < .001 |
| diet | self - control | 0.07 [-0.03, 0.17] | 2448.90 | 1.45 | .147 |
| energy | other - control | 0.18 [0.06, 0.30] | 5252.29 | 3.00 | .003 |
| energy | self - control | 0.08 [-0.04, 0.20] | 5252.29 | 1.35 | .178 |
| recycle | other - control | 0.18 [0.06, 0.30] | 5252.29 | 2.92 | .004 |
| recycle | self - control | 0.11 [-0.01, 0.23] | 5252.29 | 1.79 | .074 |
| transit | other - control | 0.18 [0.07, 0.28] | 3085.08 | 3.34 | < .001 |
| transit | self - control | -0.01 [-0.11, 0.10] | 3085.08 | -0.10 | .923 |
person-level outcomes
Do the interventions affect person-level outcomes?
categories = unique(filter(exploratory_dvs, grepl("knowledge|self_efficacy", scale_name))$scale_name)
for (category in categories){
cat(paste0('\n\n### ', category, '{.tabset}\n\n'))
data_loop = exploratory_dvs %>%
filter(scale_name == !!category)
model = lm(value ~ article_cond, data = data_loop)
print(table_model(model, lm = TRUE))
predicted = ggeffects::ggpredict(model, terms = "article_cond") %>%
data.frame()
print(
predicted %>%
mutate(x = factor(x, levels = c("self", "control", "other"))) %>%
ggplot(aes(x, predicted, color = x)) +
geom_line(aes(group = 1), size = 1, color = "black") +
geom_pointrange(aes(ymin = conf.low, ymax = conf.high), size = 1, linewidth = 1) +
scale_color_manual(values = palette_condition) +
labs(x = "", y = "predicted\n") +
plot_aes +
theme(legend.position = "none")
)
}climate_knowledge
| term | b [95% CI] | t | p |
|---|---|---|---|
| control | -0.10 [-0.17, -0.03] | -2.87 | .004 |
| other - control | 0.22 [0.10, 0.34] | 3.52 | < .001 |
| self - control | 0.20 [0.08, 0.32] | 3.26 | .001 |
self_efficacy
| term | b [95% CI] | t | p |
|---|---|---|---|
| control | -0.04 [-0.11, 0.03] | -1.18 | .238 |
| other - control | 0.13 [0.01, 0.25] | 2.12 | .034 |
| self - control | 0.04 [-0.08, 0.16] | 0.65 | .514 |
cite packages
## - Angelo Canty, B. D. Ripley (2024). _boot: Bootstrap R (S-Plus) Functions_. R package version 1.3-30. A. C. Davison, D. V. Hinkley (1997). _Bootstrap Methods and Their Applications_. Cambridge University Press, Cambridge. ISBN 0-521-57391-2, <doi:10.1017/CBO9780511802843>.
## - Bates D, Mächler M, Bolker B, Walker S (2015). "Fitting Linear Mixed-Effects Models Using lme4." _Journal of Statistical Software_, *67*(1), 1-48. doi:10.18637/jss.v067.i01 <https://doi.org/10.18637/jss.v067.i01>.
## - Bates D, Maechler M, Jagan M (2024). _Matrix: Sparse and Dense Matrix Classes and Methods_. R package version 1.7-0, <https://CRAN.R-project.org/package=Matrix>.
## - Bolker B, Robinson D (2024). _broom.mixed: Tidying Methods for Mixed Models_. R package version 0.2.9.5, <https://CRAN.R-project.org/package=broom.mixed>.
## - Bürkner P (2017). "brms: An R Package for Bayesian Multilevel Models Using Stan." _Journal of Statistical Software_, *80*(1), 1-28. doi:10.18637/jss.v080.i01 <https://doi.org/10.18637/jss.v080.i01>. Bürkner P (2018). "Advanced Bayesian Multilevel Modeling with the R Package brms." _The R Journal_, *10*(1), 395-411. doi:10.32614/RJ-2018-017 <https://doi.org/10.32614/RJ-2018-017>. Bürkner P (2021). "Bayesian Item Response Modeling in R with brms and Stan." _Journal of Statistical Software_, *100*(5), 1-54. doi:10.18637/jss.v100.i05 <https://doi.org/10.18637/jss.v100.i05>.
## - Eddelbuettel D, Francois R, Allaire J, Ushey K, Kou Q, Russell N, Ucar I, Bates D, Chambers J (2024). _Rcpp: Seamless R and C++ Integration_. R package version 1.0.13, <https://CRAN.R-project.org/package=Rcpp>. Eddelbuettel D, François R (2011). "Rcpp: Seamless R and C++ Integration." _Journal of Statistical Software_, *40*(8), 1-18. doi:10.18637/jss.v040.i08 <https://doi.org/10.18637/jss.v040.i08>. Eddelbuettel D (2013). _Seamless R and C++ Integration with Rcpp_. Springer, New York. doi:10.1007/978-1-4614-6868-4 <https://doi.org/10.1007/978-1-4614-6868-4>, ISBN 978-1-4614-6867-7. Eddelbuettel D, Balamuta J (2018). "Extending R with C++: A Brief Introduction to Rcpp." _The American Statistician_, *72*(1), 28-36. doi:10.1080/00031305.2017.1375990 <https://doi.org/10.1080/00031305.2017.1375990>.
## - Grolemund G, Wickham H (2011). "Dates and Times Made Easy with lubridate." _Journal of Statistical Software_, *40*(3), 1-25. <https://www.jstatsoft.org/v40/i03/>.
## - Kassambara A (2023). _ggpubr: 'ggplot2' Based Publication Ready Plots_. R package version 0.6.0, <https://CRAN.R-project.org/package=ggpubr>.
## - Kay M (2023). _tidybayes: Tidy Data and Geoms for Bayesian Models_. doi:10.5281/zenodo.1308151 <https://doi.org/10.5281/zenodo.1308151>, R package version 3.0.6, <http://mjskay.github.io/tidybayes/>.
## - Kleiman E (2021). _EMAtools: Data Management Tools for Real-Time Monitoring/Ecological Momentary Assessment Data_. R package version 0.1.4, <https://CRAN.R-project.org/package=EMAtools>.
## - Kuznetsova A, Brockhoff PB, Christensen RHB (2017). "lmerTest Package: Tests in Linear Mixed Effects Models." _Journal of Statistical Software_, *82*(13), 1-26. doi:10.18637/jss.v082.i13 <https://doi.org/10.18637/jss.v082.i13>.
## - Lenth R (2024). _emmeans: Estimated Marginal Means, aka Least-Squares Means_. R package version 1.10.2, <https://CRAN.R-project.org/package=emmeans>.
## - Makowski D, Lüdecke D, Patil I, Thériault R, Ben-Shachar M, Wiernik B (2023). "Automated Results Reporting as a Practical Tool to Improve Reproducibility and Methodological Best Practices Adoption." _CRAN_. <https://easystats.github.io/report/>.
## - Müller K, Wickham H (2023). _tibble: Simple Data Frames_. R package version 3.2.1, <https://CRAN.R-project.org/package=tibble>.
## - R Core Team (2024). _R: A Language and Environment for Statistical Computing_. R Foundation for Statistical Computing, Vienna, Austria. <https://www.R-project.org/>.
## - Rinker TW, Kurkiewicz D (2018). _pacman: Package Management for R_. version 0.5.0, <http://github.com/trinker/pacman>.
## - Wickham H (2016). _ggplot2: Elegant Graphics for Data Analysis_. Springer-Verlag New York. ISBN 978-3-319-24277-4, <https://ggplot2.tidyverse.org>.
## - Wickham H (2023). _forcats: Tools for Working with Categorical Variables (Factors)_. R package version 1.0.0, <https://CRAN.R-project.org/package=forcats>.
## - Wickham H (2023). _stringr: Simple, Consistent Wrappers for Common String Operations_. R package version 1.5.1, <https://CRAN.R-project.org/package=stringr>.
## - Wickham H, Averick M, Bryan J, Chang W, McGowan LD, François R, Grolemund G, Hayes A, Henry L, Hester J, Kuhn M, Pedersen TL, Miller E, Bache SM, Müller K, Ooms J, Robinson D, Seidel DP, Spinu V, Takahashi K, Vaughan D, Wilke C, Woo K, Yutani H (2019). "Welcome to the tidyverse." _Journal of Open Source Software_, *4*(43), 1686. doi:10.21105/joss.01686 <https://doi.org/10.21105/joss.01686>.
## - Wickham H, Bryan J, Barrett M, Teucher A (2024). _usethis: Automate Package and Project Setup_. R package version 2.2.3, <https://CRAN.R-project.org/package=usethis>.
## - Wickham H, François R, D'Agostino McGowan L (2024). _emo: Easily Insert 'Emoji'_. R package version 0.0.0.9000, commit 3f03b11491ce3d6fc5601e210927eff73bf8e350, <https://github.com/hadley/emo>.
## - Wickham H, François R, Henry L, Müller K, Vaughan D (2023). _dplyr: A Grammar of Data Manipulation_. R package version 1.1.4, <https://CRAN.R-project.org/package=dplyr>.
## - Wickham H, Henry L (2023). _purrr: Functional Programming Tools_. R package version 1.0.2, <https://CRAN.R-project.org/package=purrr>.
## - Wickham H, Hester J, Bryan J (2024). _readr: Read Rectangular Text Data_. R package version 2.1.5, <https://CRAN.R-project.org/package=readr>.
## - Wickham H, Hester J, Chang W, Bryan J (2022). _devtools: Tools to Make Developing R Packages Easier_. R package version 2.4.5, <https://CRAN.R-project.org/package=devtools>.
## - Wickham H, Vaughan D, Girlich M (2024). _tidyr: Tidy Messy Data_. R package version 1.3.1, <https://CRAN.R-project.org/package=tidyr>.
## - Xie Y (2024). _knitr: A General-Purpose Package for Dynamic Report Generation in R_. R package version 1.47, <https://yihui.org/knitr/>. Xie Y (2015). _Dynamic Documents with R and knitr_, 2nd edition. Chapman and Hall/CRC, Boca Raton, Florida. ISBN 978-1498716963, <https://yihui.org/knitr/>. Xie Y (2014). "knitr: A Comprehensive Tool for Reproducible Research in R." In Stodden V, Leisch F, Peng RD (eds.), _Implementing Reproducible Computational Research_. Chapman and Hall/CRC. ISBN 978-1466561595.
## - Zhu H (2024). _kableExtra: Construct Complex Table with 'kable' and Pipe Syntax_. R package version 1.4.0, <https://CRAN.R-project.org/package=kableExtra>.