add rmarkdowns

source/rmarkdown/analyses_diversity/_child_model_observed_richness.Rmd

@@ -0,0 +1,77 @@
+### `r stringr::str_to_sentence("{{group}} - {{primerset}} - {{unit}}")`
+
+```{r}
+data_subset <- combined %>%
+  filter(
+    group == "{{group}}",
+    primerset == "{{primerset}}",
+    unit == "{{unit}}"
+  )
+```
+
+
+```{r m-richness-{{group}}-{{primerset}}-{{unit}}}
+form_richness <- formula(
+  observed ~
+    log(total_count)
+  + landgebruik
+  + diepte
+  + landgebruik:diepte
+  + (1 | cmon_plot_id)
+  )
+
+cat("Fitting Poisson model")
+
+m_richness <- glmmTMB(
+  formula = form_richness,
+  data = data_subset,
+  family = poisson())
+
+test <- check_overdispersion(m_richness)
+test <- test$p_value < 0.05
+
+if (test) {
+  cat("Overdispersion detected: fitting Negative binomial model instead.")
+  m_richness <- glmmTMB(
+    formula = form_richness,
+    data = data_subset,
+    family = nbinom2())
+}
+```
+
+
+```{r m-richness-checks-{{group}}-{{primerset}}-{{unit}}}
+performance::check_model(m_richness)
+```
+
+
+```{r m-richness-summary-{{group}}-{{primerset}}-{{unit}}}
+summary(m_richness) |> print(digits = 2)
+```
+
+```{r m-richness-anova-{{group}}-{{primerset}}-{{unit}}}
+car::Anova(m_richness) |> print(digits = 2)
+```
+
+```{r m-richness-preds-total-count-{{group}}-{{primerset}}-{{unit}}}
+marginaleffects::plot_predictions(
+  m_richness,
+  condition = c("total_count"),
+  re.form = NA,
+  vcov = TRUE,
+  type = "response") +
+  labs(y = "Voorspeld aantal taxa") +
+  scale_x_log10()
+```
+
+```{r m-richness-preds-landgebruikxdiepte-{{group}}-{{primerset}}-{{unit}}}
+marginaleffects::plot_predictions(
+  m_richness,
+  condition = c("landgebruik", "diepte"),
+  re.form = NA,
+  vcov = TRUE,
+  type = "response") +
+  labs(y = "Voorspeld aantal taxa")
+```
+
+

source/rmarkdown/analyses_diversity/_child_model_shannon.Rmd

@@ -0,0 +1,78 @@
+### `r stringr::str_to_sentence("{{group}} - {{primerset}} - {{unit}}")`
+
+```{r}
+data_subset <- combined %>%
+  filter(
+    group == "{{group}}",
+    primerset == "{{primerset}}",
+    unit == "{{unit}}"
+  )
+```
+
+
+```{r m-shannon-{{group}}-{{primerset}}-{{unit}}}
+form_shannon <- formula(
+  shannon ~
+    log(total_count)
+  + landgebruik
+  + diepte
+  + landgebruik:diepte
+  + (1 | cmon_plot_id)
+  )
+
+zeroes <- min(data_subset$shannon) == 0
+
+if (!zeroes) {
+  cat("Fitting Gamma model")
+
+  m_shannon <- glmmTMB(
+    formula = form_shannon,
+    data = data_subset,
+    family = Gamma(link = "log"))
+} else {
+  cat("Zeroes detected (only one taxum observed): Fitting zero-inflated Gamma model")
+
+  m_shannon <- glmmTMB(
+    formula = form_shannon,
+    ziformula = ~ 1,
+    data = data_subset,
+    family = ziGamma(link = "log"))
+}
+```
+
+
+```{r m-shannon-checks-{{group}}-{{primerset}}-{{unit}}}
+performance::check_model(m_shannon)
+```
+
+
+```{r m-shannon-summary-{{group}}-{{primerset}}-{{unit}}}
+summary(m_shannon) |> print(digits = 2)
+```
+
+```{r m-shannon-anova-{{group}}-{{primerset}}-{{unit}}}
+car::Anova(m_shannon) |> print(digits = 2)
+```
+
+```{r m-shannon-preds-total-count-{{group}}-{{primerset}}-{{unit}}}
+marginaleffects::plot_predictions(
+  m_shannon,
+  condition = c("total_count"),
+  re.form = NA,
+  vcov = TRUE,
+  type = "response") +
+  labs(y = "Shannon index") +
+  scale_x_log10()
+```
+
+```{r m-shannon-preds-landgebruikxdiepte-{{group}}-{{primerset}}-{{unit}}}
+marginaleffects::plot_predictions(
+  m_shannon,
+  condition = c("landgebruik", "diepte"),
+  re.form = NA,
+  vcov = TRUE,
+  type = "response") +
+  labs(y = "Shannon index")
+```
+
+

source/rmarkdown/analyses_diversity/_child_model_simpson.Rmd

@@ -0,0 +1,66 @@
+### `r stringr::str_to_sentence("{{group}} - {{primerset}} - {{unit}}")`
+
+```{r}
+data_subset <- combined %>%
+  filter(
+    group == "{{group}}",
+    primerset == "{{primerset}}",
+    unit == "{{unit}}"
+  )
+```
+
+
+```{r m-simpson-{{group}}-{{primerset}}-{{unit}}}
+form_simpson <- formula(
+  simpson ~
+    log(total_count)
+  + landgebruik
+  + diepte
+  + landgebruik:diepte
+  + (1 | cmon_plot_id)
+  )
+
+cat("Fitting ordinal beta model")
+
+m_simpson <- glmmTMB(
+  formula = form_simpson,
+  data = data_subset,
+  family = ordbeta())
+```
+
+
+```{r m-simpson-checks-{{group}}-{{primerset}}-{{unit}}, error = TRUE}
+performance::check_model(m_simpson)
+```
+
+
+```{r m-simpson-summary-{{group}}-{{primerset}}-{{unit}}}
+summary(m_simpson) |> print(digits = 2)
+```
+
+```{r m-simpson-anova-{{group}}-{{primerset}}-{{unit}}, error = TRUE}
+car::Anova(m_simpson) |> print(digits = 2)
+```
+
+```{r m-simpson-preds-total-count-{{group}}-{{primerset}}-{{unit}}, error = TRUE}
+marginaleffects::plot_predictions(
+  m_simpson,
+  condition = c("total_count"),
+  re.form = NA,
+  vcov = TRUE,
+  type = "response") +
+  labs(y = "Simpson index") +
+  scale_x_log10()
+```
+
+```{r m-simpson-preds-landgebruikxdiepte-{{group}}-{{primerset}}-{{unit}}, error = TRUE}
+marginaleffects::plot_predictions(
+  m_simpson,
+  condition = c("landgebruik", "diepte"),
+  re.form = NA,
+  vcov = TRUE,
+  type = "response") +
+  labs(y = "Simpson index")
+```
+
+