# Laden der relevanten Pakete
library(lme4) # Mehrebenenmodelle
library(lmerTest) # Signifikanztests für Mehrebenenmodelle
library(performance) # Modellbeurteilung (hier vor allem ICC)
library(parameters) # Modellschätzungen extrahieren und darstellen
library(report) # Einfaches Erstellen von statistischen Berichten
library(tidyverse) # Datenmanagement und Visualisierung: https://www.tidyverse.org/
theme_set(theme_classic(base_size = 12)) # Voreinstellung für Plots


# Daten laden
# Daten aus https://osf.io/7mdbn/
# Aufbereitung nach https://osf.io/7mdbn/files/wtbx2
# Hier Auszug für Mehrebenenanalyse
# Aufbereiteter Datensatz in Blackboard
d <- readRDS("esm.rds") |>
  mutate(pod_duration_total = pod_duration_total / 60) |> # In Stunden
  na.omit()


# Auszug aus Datensatz
d |>
  select(serial, gender, date, time, pod_duration_total, narr_trans) |>
  group_by(serial) |>
  slice_head(n = 3) |>
  head(n = 15)


# Histogramm Dauer
d |>
  ggplot(aes(pod_duration_total)) +
  geom_histogram()


# Histogramm Transportation
d |>
  ggplot(aes(narr_trans)) +
  geom_histogram(binwidth = 0.5)


# Beschreibung der Situationsvariablen
d |>
  select(pod_duration_total, narr_trans) |>
  report_table()


# Naiv Situationen Plot
d |>
  ggplot(aes(pod_duration_total, narr_trans)) +
  geom_point(position = position_jitter(), shape = 1) +
  stat_smooth(method = "lm")


# Naiv Situationen Regression
lm(narr_trans ~ pod_duration_total, data = d) |>
  report_table(metrics = "R2_adj", include_effectsize = FALSE)


# Aggregieren auf Personenebene
d <- d |>
  group_by(serial) |>
  mutate(
    pmc_pod_duration_personmean = mean(pod_duration_total, na.rm = TRUE),
    pmc_narr_trans_personmean = mean(narr_trans, na.rm = TRUE)
  )


# Personen Plot
d |>
  distinct(serial, .keep_all = TRUE) |>
  ggplot(aes(pmc_pod_duration_personmean, pmc_narr_trans_personmean)) +
  geom_point(position = position_jitter(), shape = 1) +
  stat_smooth(method = "lm")


# Personen Regression
lm(pmc_narr_trans_personmean ~ pmc_pod_duration_personmean, data = distinct(d, serial, .keep_all = TRUE)) |>
  report_table(metrics = "R2_adj", include_effectsize = FALSE)


# Zentrieren um Personenmittelwert
d <- d |>
  group_by(serial) |>
  mutate(
    pmc_pod_duration_total = pod_duration_total - pmc_pod_duration_personmean,
    pmc_narr_trans = narr_trans - pmc_narr_trans_personmean
  )
d |>
  select(
    serial,
    pod_duration_total,
    pmc_pod_duration_personmean,
    pmc_pod_duration_total
  ) |>
  group_by(serial) |>
  slice_head(n = 3) |>
  head(n = 15)


# Zentriert um Personenmittelwert Plot
d |>
  ggplot(aes(pmc_pod_duration_total, pmc_narr_trans)) +
  geom_point(position = position_jitter(), shape = 1) +
  stat_smooth(method = "lm")


# Zentriert um Personenmittelwert Regression
lm(pmc_narr_trans ~ pmc_pod_duration_total, data = d) |>
  report_table(metrics = "R2_adj", include_effectsize = FALSE)


# 
m0 <- lmer(narr_trans ~ 1 + (1 | serial), data = d)
m0 |>
  report_table(metrics = "R2", include_effectsize = FALSE)
icc(m0, by_group = TRUE)


# 
lmer(narr_trans ~ pmc_pod_duration_total + pmc_pod_duration_personmean + (1 | serial), data = d) |>
  report_table(metrics = "R2", include_effectsize = FALSE)