Temporal Dynamics of NREM EEG Microstates across
the Night

Lorena R.R. Gianotti, André Minder, Mirjam Studler, Daria Knoch

Department of Social Neuroscience and Social Psychology,
Institute of Psychology, University of Bern, Switzerland

Last edited: 2026-03-25

---

Preparation

# Load packages
library(janitor)
library(ggthemes)
library(readxl)
library(lmerTest)
library(kableExtra)
library(emmeans)
library(hms)
library(sjPlot)
library(lme4)
library(nlme)
library(tidyverse)
library(gridExtra)
library(ggplot2)

# Define parameters
size_big   <- 14
size_small <- 12

# Load microstates data
data_ms <- readr::read_csv("data/cycles/TemporalParameters_5 classes_GrandGrandMeanTemplate.csv",
                           show_col_types = FALSE)

data_ms_clean <- data_ms |>
janitor::clean_names() |>
mutate(subject = parse_number(dataset),
.before = subject,
.keep = "unused") |>
rename(ID = subject)

# Pivot to long format

## duration

data_long_duration <- data_ms_clean |>
select(ID, group, condition, contains("mean_duration"), -contains("all")) |>
pivot_longer(
cols = -c("ID", "group", "condition"),
names_to = "class",
values_to = "duration"
) |>
mutate(condition = stringr::str_sub(condition, 6, -1)) |>
rename(cycle = condition,
stage = group) |>
mutate(
class = stringr::str_sub(class, -1, -1) |> toupper(),
cycle = as.factor(cycle),
class = as.factor(class),
stage = as.factor(stage),
ID    = as.factor(ID),
duration = duration * 1000,
stage = forcats::fct_relevel(stage, c("N1", "N2", "N3")),
class = forcats::fct_relevel(class, c("A", "B", "C", "D", "E"))
)

## occurrence

data_long_occurrence <- data_ms_clean |>
select(ID, group, condition, contains("mean_occurrence"), -contains("all")) |>
pivot_longer(
cols = -c("ID", "group", "condition"),
names_to = "class",
values_to = "occurrence"
) |>
mutate(condition = stringr::str_sub(condition, 6, -1)) |>
rename(cycle = condition,
stage = group) |>
mutate(
class = stringr::str_sub(class, -1, -1) |> toupper(),
cycle = as.factor(cycle),
class = as.factor(class),
stage = as.factor(stage),
ID    = as.factor(ID),
stage = forcats::fct_relevel(stage, c("N1", "N2", "N3")),
class = forcats::fct_relevel(class, c("A", "B", "C", "D", "E"))
)

## coverage

data_long_coverage <- data_ms_clean |>
select(ID, group, condition, contains("coverage")) |>
pivot_longer(
cols = -c("ID", "group", "condition"),
names_to = "class",
values_to = "coverage"
) |>
mutate(condition = stringr::str_sub(condition, 6, -1)) |>
rename(cycle = condition,
stage = group) |>
mutate(
class = stringr::str_sub(class, -1, -1) |> toupper(),
cycle = as.factor(cycle),
class = as.factor(class),
stage = as.factor(stage),
ID    = as.factor(ID),
stage = forcats::fct_relevel(stage, c("N1", "N2", "N3")),
class = forcats::fct_relevel(class, c("A", "B", "C", "D", "E"))
)

Descriptive Statistics

Time segments analysed for each of the 12 combinations of sleep stages and sleep cycles

# count

tab_n <- data_ms_clean |>
count(group, condition, name = "n")

# compute mean and SD (seconds) + round

tab_time <- data_ms_clean |>
group_by(group, condition) |>
summarise(
mean_total_time = mean(total_time, na.rm = TRUE),
sd_total_time   = sd(total_time, na.rm = TRUE),
.groups = "drop"
) |>
mutate(
mean_total_time = round(mean_total_time, 0),
sd_total_time   = round(sd_total_time, 0),
time_summary = paste0(mean_total_time, " (± ", sd_total_time, ")")
) |>
select(group, condition, time_summary)

# merge

table_final <- tab_n |>
left_join(tab_time, by = c("group", "condition")) |>
arrange(group, condition)

# create table

table_final |>
knitr::kable(
format   = "html",
booktabs = TRUE,
caption  = NULL,
col.names = c("Sleep stage", "Sleep cycle", "n",
"Time segments (mean ± SD) in seconds"),
align = c("l","l","c","c"),
escape = TRUE
) |>
kableExtra::kable_styling(
full_width = FALSE,
position = "left",
bootstrap_options = c(),
htmltable_class = "apa-table"
)

Sleep stage	Sleep cycle	n	Time segments (mean ± SD) in seconds
N1	cycle1	48	200 (± 157)
N1	cycle2	37	178 (± 286)
N1	cycle3	34	123 (± 108)
N1	cycle4	31	123 (± 128)
N2	cycle1	54	1446 (± 602)
N2	cycle2	54	2030 (± 1057)
N2	cycle3	53	2348 (± 785)
N2	cycle4	52	2390 (± 632)
N3	cycle1	53	2608 (± 900)
N3	cycle2	51	1627 (± 710)
N3	cycle3	49	1105 (± 710)
N3	cycle4	36	583 (± 538)

Inference Statistics of the Microstate Temporal Characteristics

Duration

Model

duration_model <- lmer(
duration ~ cycle * stage * class +
(1|ID) + (1|cycle:ID) + (1|stage:ID) + (1|class:ID),
data = data_long_duration, REML = TRUE
)

cat('<div class="scroll-y">')

sjPlot::tab_model(duration_model)

	duration
Predictors	Estimates	CI	p
(Intercept)	48.86	45.38 – 52.34	<0.001
cycle [2]	10.42	6.20 – 14.64	<0.001
cycle [3]	5.23	0.89 – 9.57	0.018
cycle [4]	2.39	-2.07 – 6.86	0.293
stage [N2]	15.79	11.99 – 19.59	<0.001
stage [N3]	61.22	57.41 – 65.03	<0.001
class [B]	-3.08	-6.88 – 0.71	0.111
class [C]	1.69	-2.10 – 5.49	0.383
class [D]	-5.08	-8.88 – -1.29	0.009
class [E]	-7.83	-11.63 – -4.04	<0.001
cycle [2] × stage [N2]	-8.06	-13.05 – -3.07	0.002
cycle [3] × stage [N2]	-6.22	-11.32 – -1.12	0.017
cycle [4] × stage [N2]	-5.87	-11.08 – -0.66	0.027
cycle [2] × stage [N3]	-19.32	-24.36 – -14.28	<0.001
cycle [3] × stage [N3]	-18.97	-24.12 – -13.82	<0.001
cycle [4] × stage [N3]	-27.07	-32.55 – -21.60	<0.001
cycle [2] × class [B]	-0.88	-6.14 – 4.38	0.742
cycle [3] × class [B]	-0.52	-5.92 – 4.89	0.852
cycle [4] × class [B]	-0.35	-5.92 – 5.22	0.901
cycle [2] × class [C]	15.97	10.70 – 21.23	<0.001
cycle [3] × class [C]	8.32	2.91 – 13.73	0.003
cycle [4] × class [C]	8.38	2.81 – 13.95	0.003
cycle [2] × class [D]	6.57	1.31 – 11.83	0.014
cycle [3] × class [D]	8.60	3.19 – 14.01	0.002
cycle [4] × class [D]	6.87	1.30 – 12.44	0.016
cycle [2] × class [E]	0.90	-4.36 – 6.16	0.737
cycle [3] × class [E]	2.65	-2.76 – 8.05	0.338
cycle [4] × class [E]	1.50	-4.07 – 7.07	0.597
stage [N2] × class [B]	-1.45	-6.20 – 3.30	0.549
stage [N3] × class [B]	-6.43	-11.20 – -1.66	0.008
stage [N2] × class [C]	9.97	5.22 – 14.73	<0.001
stage [N3] × class [C]	29.62	24.85 – 34.39	<0.001
stage [N2] × class [D]	10.64	5.89 – 15.39	<0.001
stage [N3] × class [D]	15.25	10.48 – 20.02	<0.001
stage [N2] × class [E]	5.75	1.00 – 10.50	0.018
stage [N3] × class [E]	1.40	-3.37 – 6.17	0.566
(cycle [2] × stage [N2]) × class [B]	0.49	-6.50 – 7.47	0.892
(cycle [3] × stage [N2]) × class [B]	0.60	-6.51 – 7.72	0.868
(cycle [4] × stage [N2]) × class [B]	1.06	-6.18 – 8.31	0.773
(cycle [2] × stage [N3]) × class [B]	2.66	-4.39 – 9.71	0.460
(cycle [3] × stage [N3]) × class [B]	2.80	-4.39 – 9.98	0.445
(cycle [4] × stage [N3]) × class [B]	2.12	-5.48 – 9.72	0.584
(cycle [2] × stage [N2]) × class [C]	-14.44	-21.43 – -7.45	<0.001
(cycle [3] × stage [N2]) × class [C]	-8.26	-15.37 – -1.14	0.023
(cycle [4] × stage [N2]) × class [C]	-8.60	-15.85 – -1.35	0.020
(cycle [2] × stage [N3]) × class [C]	-19.94	-26.99 – -12.89	<0.001
(cycle [3] × stage [N3]) × class [C]	-15.60	-22.79 – -8.42	<0.001
(cycle [4] × stage [N3]) × class [C]	-21.91	-29.51 – -14.31	<0.001
(cycle [2] × stage [N2]) × class [D]	-6.27	-13.26 – 0.72	0.078
(cycle [3] × stage [N2]) × class [D]	-9.24	-16.36 – -2.13	0.011
(cycle [4] × stage [N2]) × class [D]	-6.45	-13.69 – 0.80	0.081
(cycle [2] × stage [N3]) × class [D]	-8.28	-15.33 – -1.23	0.021
(cycle [3] × stage [N3]) × class [D]	-13.25	-20.43 – -6.07	<0.001
(cycle [4] × stage [N3]) × class [D]	-10.19	-17.79 – -2.59	0.009
(cycle [2] × stage [N2]) × class [E]	-0.71	-7.70 – 6.27	0.841
(cycle [3] × stage [N2]) × class [E]	-1.51	-8.63 – 5.60	0.676
(cycle [4] × stage [N2]) × class [E]	0.49	-6.76 – 7.74	0.894
(cycle [2] × stage [N3]) × class [E]	0.79	-6.26 – 7.83	0.827
(cycle [3] × stage [N3]) × class [E]	0.03	-7.15 – 7.22	0.992
(cycle [4] × stage [N3]) × class [E]	2.41	-5.19 – 10.01	0.534
Random Effects
σ2	74.25
τ00 class:ID	16.90
τ00 cycle:ID	24.38
τ00 stage:ID	20.99
τ00 ID	21.43
ICC	0.53
N ID	54
N cycle	4
N stage	3
N class	5
Observations	2760
Marginal R2 / Conditional R2	0.778 / 0.896

cat('</div>')

Type III Test

anova_duration <- as.data.frame(anova(duration_model))

p_raw <- anova_duration$`Pr(>F)`

anova_duration$`Pr(>F)` <- ifelse(
  p_raw < .001,
  "<b>&lt;0.001</b>",
  ifelse(
    p_raw < .05,
    paste0("<b>", sprintf("%.3f", p_raw), "</b>"),
    sprintf("%.3f", p_raw)
  )
)

cat('<div class="scroll-y">')

sjPlot::tab_df(
  anova_duration,
  title = "Type III Test",
  digits = 2,
  encoding = "UTF-8",
  escape = FALSE
)

Type III Test

Sum.Sq	Mean.Sq	NumDF	DenDF	F.value	Pr..F.
7636.27	2545.42	3	159.85	34.28	<0.001
193400.41	96700.20	2	103.66	1302.30	<0.001
46906.02	11726.51	4	205.19	157.93	<0.001
66072.30	11012.05	6	2311.75	148.30	<0.001
2730.65	227.55	12	2194.83	3.06	<0.001
26314.89	3289.36	8	2188.36	44.30	<0.001
7216.20	300.68	24	2180.24	4.05	<0.001

cat('</div>')

Post-hoc Tests

emm_duration <- emmeans(duration_model, ~ cycle | stage * class)

contrast_means_duration <- as.data.frame(pairs(emm_duration, adjust = "holm")) |>
  dplyr::select(stage, class, contrast, estimate, SE, df, t.ratio, p.value) |>
  dplyr::mutate(
    stage    = factor(stage, levels = c("N1", "N2", "N3")),
    estimate = round(estimate, 2),
    SE       = round(SE, 2),
    df       = round(df, 2),
    t.ratio  = round(t.ratio, 2)
  ) |>
  dplyr::arrange(stage, class, contrast)

for (stg in c("N1", "N2", "N3")) {
  
  tab <- contrast_means_duration |>
    dplyr::filter(stage == stg) |>
    dplyr::select(-stage)
  
  p_raw <- tab$p.value
  
  tab$p.value <- ifelse(
    p_raw < .001,
    "<b>&lt;0.001</b>",
    ifelse(
      p_raw < .05,
      paste0("<b>", sprintf("%.3f", p_raw), "</b>"),
      sprintf("%.3f", p_raw)
    )
  )
  
  colnames(tab) <- c("Class", "Cycle contrast", "Estimate", "SE", "df", "t", "p")
  
  cat('<div class="scroll-y">')
  cat(paste0("<h4>Duration: ", stg, "</h4>"))
  
   tab_out <- sjPlot::tab_df(
    tab,
    title = NULL,
    digits = 2,
    encoding = "UTF-8",
    escape = FALSE
  )

   cat(tab_out$knitr)
   cat("</div><br>")
}

Duration: N1

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	-10.42	2.15	1521.67	-4.84	<0.001
A	cycle1 - cycle3	-5.23	2.21	1588.34	-2.36	0.073
A	cycle1 - cycle4	-2.39	2.28	1636.51	-1.05	0.485
A	cycle2 - cycle3	5.19	2.33	1725.57	2.22	0.079
A	cycle2 - cycle4	8.03	2.39	1760.47	3.36	0.004
A	cycle3 - cycle4	2.83	2.42	1775.76	1.17	0.485
B	cycle1 - cycle2	-9.54	2.15	1521.67	-4.43	<0.001
B	cycle1 - cycle3	-4.71	2.21	1588.34	-2.13	0.134
B	cycle1 - cycle4	-2.04	2.28	1636.51	-0.90	0.541
B	cycle2 - cycle3	4.82	2.33	1725.57	2.07	0.134
B	cycle2 - cycle4	7.50	2.39	1760.47	3.13	0.009
B	cycle3 - cycle4	2.67	2.42	1775.76	1.10	0.541
C	cycle1 - cycle2	-26.39	2.15	1521.67	-12.27	<0.001
C	cycle1 - cycle3	-13.55	2.21	1588.34	-6.12	<0.001
C	cycle1 - cycle4	-10.77	2.28	1636.51	-4.73	<0.001
C	cycle2 - cycle3	12.83	2.33	1725.57	5.50	<0.001
C	cycle2 - cycle4	15.62	2.39	1760.47	6.53	<0.001
C	cycle3 - cycle4	2.78	2.42	1775.76	1.15	0.251
D	cycle1 - cycle2	-16.99	2.15	1521.67	-7.90	<0.001
D	cycle1 - cycle3	-13.83	2.21	1588.34	-6.25	<0.001
D	cycle1 - cycle4	-9.26	2.28	1636.51	-4.07	<0.001
D	cycle2 - cycle3	3.16	2.33	1725.57	1.36	0.175
D	cycle2 - cycle4	7.73	2.39	1760.47	3.23	0.004
D	cycle3 - cycle4	4.57	2.42	1775.76	1.88	0.120
E	cycle1 - cycle2	-11.32	2.15	1521.67	-5.26	<0.001
E	cycle1 - cycle3	-7.87	2.21	1588.34	-3.56	0.002
E	cycle1 - cycle4	-3.90	2.28	1636.51	-1.71	0.262
E	cycle2 - cycle3	3.45	2.33	1725.57	1.48	0.262
E	cycle2 - cycle4	7.42	2.39	1760.47	3.10	0.008
E	cycle3 - cycle4	3.98	2.42	1775.76	1.64	0.262

Duration: N2

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	-2.36	1.91	1188.83	-1.24	0.602
A	cycle1 - cycle3	0.99	1.92	1189.40	0.52	0.605
A	cycle1 - cycle4	3.48	1.93	1189.42	1.80	0.361
A	cycle2 - cycle3	3.36	1.92	1189.40	1.75	0.361
A	cycle2 - cycle4	5.84	1.93	1189.42	3.02	0.015
A	cycle3 - cycle4	2.48	1.94	1188.87	1.28	0.602
B	cycle1 - cycle2	-1.96	1.91	1188.83	-1.03	0.913
B	cycle1 - cycle3	0.91	1.92	1189.40	0.47	0.913
B	cycle1 - cycle4	2.77	1.93	1189.42	1.43	0.678
B	cycle2 - cycle3	2.87	1.92	1189.40	1.49	0.678
B	cycle2 - cycle4	4.73	1.93	1189.42	2.45	0.087
B	cycle3 - cycle4	1.86	1.94	1188.87	0.96	0.913
C	cycle1 - cycle2	-3.89	1.91	1188.83	-2.03	0.169
C	cycle1 - cycle3	0.93	1.92	1189.40	0.48	0.630
C	cycle1 - cycle4	3.70	1.93	1189.42	1.91	0.169
C	cycle2 - cycle3	4.81	1.92	1189.40	2.50	0.062
C	cycle2 - cycle4	7.59	1.93	1189.42	3.93	<0.001
C	cycle3 - cycle4	2.77	1.94	1188.87	1.43	0.306
D	cycle1 - cycle2	-2.66	1.91	1188.83	-1.39	0.492
D	cycle1 - cycle3	1.64	1.92	1189.40	0.85	0.789
D	cycle1 - cycle4	3.06	1.93	1189.42	1.58	0.457
D	cycle2 - cycle3	4.30	1.92	1189.40	2.24	0.128
D	cycle2 - cycle4	5.72	1.93	1189.42	2.96	0.019
D	cycle3 - cycle4	1.42	1.94	1188.87	0.73	0.789
E	cycle1 - cycle2	-2.55	1.91	1188.83	-1.33	0.913
E	cycle1 - cycle3	-0.14	1.92	1189.40	-0.07	1.000
E	cycle1 - cycle4	1.48	1.93	1189.42	0.77	1.000
E	cycle2 - cycle3	2.41	1.92	1189.40	1.25	0.913
E	cycle2 - cycle4	4.03	1.93	1189.42	2.09	0.223
E	cycle3 - cycle4	1.62	1.94	1188.87	0.84	1.000

Duration: N3

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	8.90	1.95	1236.88	4.58	<0.001
A	cycle1 - cycle3	13.74	1.97	1249.62	6.98	<0.001
A	cycle1 - cycle4	24.68	2.14	1453.67	11.55	<0.001
A	cycle2 - cycle3	4.84	1.98	1275.17	2.44	0.015
A	cycle2 - cycle4	15.78	2.15	1477.57	7.33	<0.001
A	cycle3 - cycle4	10.94	2.17	1486.24	5.04	<0.001
B	cycle1 - cycle2	7.13	1.95	1236.88	3.66	<0.001
B	cycle1 - cycle3	11.46	1.97	1249.62	5.82	<0.001
B	cycle1 - cycle4	22.91	2.14	1453.67	10.73	<0.001
B	cycle2 - cycle3	4.33	1.98	1275.17	2.18	0.029
B	cycle2 - cycle4	15.78	2.15	1477.57	7.33	<0.001
B	cycle3 - cycle4	11.45	2.17	1486.24	5.27	<0.001
C	cycle1 - cycle2	12.88	1.95	1236.88	6.62	<0.001
C	cycle1 - cycle3	21.02	1.97	1249.62	10.68	<0.001
C	cycle1 - cycle4	38.22	2.14	1453.67	17.89	<0.001
C	cycle2 - cycle3	8.14	1.98	1275.17	4.10	<0.001
C	cycle2 - cycle4	25.34	2.15	1477.57	11.77	<0.001
C	cycle3 - cycle4	17.20	2.17	1486.24	7.92	<0.001
D	cycle1 - cycle2	10.61	1.95	1236.88	5.46	<0.001
D	cycle1 - cycle3	18.39	1.97	1249.62	9.35	<0.001
D	cycle1 - cycle4	28.00	2.14	1453.67	13.11	<0.001
D	cycle2 - cycle3	7.78	1.98	1275.17	3.92	<0.001
D	cycle2 - cycle4	17.39	2.15	1477.57	8.08	<0.001
D	cycle3 - cycle4	9.61	2.17	1486.24	4.43	<0.001
E	cycle1 - cycle2	7.22	1.95	1236.88	3.71	<0.001
E	cycle1 - cycle3	11.06	1.97	1249.62	5.62	<0.001
E	cycle1 - cycle4	20.76	2.14	1453.67	9.72	<0.001
E	cycle2 - cycle3	3.84	1.98	1275.17	1.94	0.053
E	cycle2 - cycle4	13.55	2.15	1477.57	6.29	<0.001
E	cycle3 - cycle4	9.70	2.17	1486.24	4.47	<0.001

Plot

duration_model_emm <- emmeans(duration_model, ~ cycle * stage | class)
duration_model_emm <- as.data.frame(duration_model_emm)

duration_emm_plot <- duration_model_emm |> 
  select(class, cycle, stage, emmean, lower.CL, upper.CL) |> 
  slice(1:60) |> 
  ggplot(aes(group = stage, color = stage, y = emmean, x = cycle)) +
  geom_line(fill = cycle, alpha = 1) +
  geom_point(stat="identity", size = 3, alpha = 1) +
  geom_errorbar(aes(ymin=lower.CL, ymax=upper.CL), width=.2, alpha = 1) +
  ggtitle("Microstate Classes") +
  facet_wrap( ~ class, ncol = 5) +
  theme_bw() +
  theme(
    plot.title = element_text(size = size_big, hjust= 0.5),
    axis.text = element_text(size = size_small, family = "sans"),
    axis.title = element_text(size = size_big, family = "sans"),
    strip.text.x = element_text(size = size_big, family = "sans"),
    strip.background = element_blank(),
    legend.title = element_blank(),
    legend.text = element_text(size = size_small),
    legend.position="right",
  ) +
  labs(color = "Sleep Stages") +
  scale_y_continuous(name = "Duration (ms)", limits = c(0,150), expand = c(0,0)) +
  scale_x_discrete(name= "Sleep Cycles") +
  scale_color_manual(values = c("#1b9e77","#d95f02","#7570b3")) #color-blind friendly colors

duration_emm_plot

ggsave("duration_plot.png",
        plot = duration_emm_plot,
        width = 9,
        height = 7, 
        dpi = 1200,
        bg = "transparent",
        path = "R:/MicrostateAnalysis/sleepData/MS_SleepMarkers/scripts/graphs/sleep/cycles")

Occurence

Model

occurrence_model <- lmer(
occurrence ~ cycle * stage * class +
(1|ID) + (1|cycle:ID) + (1|stage:ID) + (1|class:ID),
data = data_long_occurrence, REML = TRUE
)

cat('<div class="scroll-y">')

sjPlot::tab_model(occurrence_model)

	occurrence
Predictors	Estimates	CI	p
(Intercept)	4.65	4.52 – 4.79	<0.001
cycle [2]	-1.23	-1.38 – -1.09	<0.001
cycle [3]	-1.02	-1.17 – -0.87	<0.001
cycle [4]	-0.81	-0.96 – -0.66	<0.001
stage [N2]	-1.80	-1.93 – -1.67	<0.001
stage [N3]	-3.02	-3.15 – -2.89	<0.001
class [B]	-0.15	-0.29 – -0.01	0.036
class [C]	0.28	0.14 – 0.42	<0.001
class [D]	-0.63	-0.77 – -0.48	<0.001
class [E]	-0.87	-1.01 – -0.73	<0.001
cycle [2] × stage [N2]	1.12	0.96 – 1.29	<0.001
cycle [3] × stage [N2]	1.03	0.86 – 1.20	<0.001
cycle [4] × stage [N2]	0.85	0.68 – 1.02	<0.001
cycle [2] × stage [N3]	1.44	1.27 – 1.61	<0.001
cycle [3] × stage [N3]	1.31	1.14 – 1.48	<0.001
cycle [4] × stage [N3]	1.40	1.21 – 1.58	<0.001
cycle [2] × class [B]	-0.09	-0.27 – 0.08	0.288
cycle [3] × class [B]	-0.10	-0.28 – 0.08	0.263
cycle [4] × class [B]	-0.19	-0.37 – -0.00	0.047
cycle [2] × class [C]	0.42	0.25 – 0.59	<0.001
cycle [3] × class [C]	0.40	0.22 – 0.58	<0.001
cycle [4] × class [C]	0.42	0.23 – 0.60	<0.001
cycle [2] × class [D]	0.70	0.53 – 0.88	<0.001
cycle [3] × class [D]	0.86	0.68 – 1.04	<0.001
cycle [4] × class [D]	0.79	0.61 – 0.98	<0.001
cycle [2] × class [E]	0.30	0.13 – 0.47	0.001
cycle [3] × class [E]	0.35	0.17 – 0.53	<0.001
cycle [4] × class [E]	0.29	0.11 – 0.48	0.002
stage [N2] × class [B]	-0.11	-0.27 – 0.05	0.168
stage [N3] × class [B]	-0.05	-0.21 – 0.11	0.530
stage [N2] × class [C]	0.32	0.17 – 0.48	<0.001
stage [N3] × class [C]	0.20	0.04 – 0.36	0.014
stage [N2] × class [D]	0.94	0.78 – 1.10	<0.001
stage [N3] × class [D]	0.82	0.67 – 0.98	<0.001
stage [N2] × class [E]	0.87	0.71 – 1.02	<0.001
stage [N3] × class [E]	0.81	0.65 – 0.97	<0.001
(cycle [2] × stage [N2]) × class [B]	0.08	-0.15 – 0.32	0.479
(cycle [3] × stage [N2]) × class [B]	0.09	-0.15 – 0.32	0.467
(cycle [4] × stage [N2]) × class [B]	0.21	-0.03 – 0.45	0.093
(cycle [2] × stage [N3]) × class [B]	0.10	-0.14 – 0.33	0.420
(cycle [3] × stage [N3]) × class [B]	0.12	-0.12 – 0.36	0.332
(cycle [4] × stage [N3]) × class [B]	0.15	-0.10 – 0.40	0.243
(cycle [2] × stage [N2]) × class [C]	-0.41	-0.65 – -0.18	<0.001
(cycle [3] × stage [N2]) × class [C]	-0.34	-0.57 – -0.10	0.005
(cycle [4] × stage [N2]) × class [C]	-0.28	-0.52 – -0.04	0.024
(cycle [2] × stage [N3]) × class [C]	-0.45	-0.68 – -0.22	<0.001
(cycle [3] × stage [N3]) × class [C]	-0.34	-0.58 – -0.10	0.005
(cycle [4] × stage [N3]) × class [C]	-0.37	-0.62 – -0.12	0.004
(cycle [2] × stage [N2]) × class [D]	-0.75	-0.99 – -0.52	<0.001
(cycle [3] × stage [N2]) × class [D]	-0.87	-1.11 – -0.64	<0.001
(cycle [4] × stage [N2]) × class [D]	-0.66	-0.90 – -0.42	<0.001
(cycle [2] × stage [N3]) × class [D]	-0.75	-0.98 – -0.51	<0.001
(cycle [3] × stage [N3]) × class [D]	-0.88	-1.12 – -0.64	<0.001
(cycle [4] × stage [N3]) × class [D]	-0.74	-1.00 – -0.49	<0.001
(cycle [2] × stage [N2]) × class [E]	-0.27	-0.50 – -0.03	0.024
(cycle [3] × stage [N2]) × class [E]	-0.26	-0.49 – -0.02	0.032
(cycle [4] × stage [N2]) × class [E]	-0.13	-0.37 – 0.11	0.294
(cycle [2] × stage [N3]) × class [E]	-0.31	-0.54 – -0.07	0.010
(cycle [3] × stage [N3]) × class [E]	-0.26	-0.49 – -0.02	0.035
(cycle [4] × stage [N3]) × class [E]	-0.24	-0.49 – 0.02	0.067
Random Effects
σ2	0.08
τ00 class:ID	0.05
τ00 cycle:ID	0.03
τ00 stage:ID	0.03
τ00 ID	0.03
ICC	0.64
N ID	54
N cycle	4
N stage	3
N class	5
Observations	2760
Marginal R2 / Conditional R2	0.760 / 0.914

cat('</div>')

Type III Test

anova_occurrence <- as.data.frame(anova(occurrence_model))

p_raw <- anova_occurrence$`Pr(>F)`

anova_occurrence$`Pr(>F)` <- ifelse(
  p_raw < .001,
  "<b>&lt;0.001</b>",
  ifelse(
    p_raw < .05,
    paste0("<b>", sprintf("%.3f", p_raw), "</b>"),
    sprintf("%.3f", p_raw)
  )
)

cat('<div class="scroll-y">')

sjPlot::tab_df(
  anova_occurrence,
  title = "Type III Test",
  digits = 2,
  encoding = "UTF-8",
  escape = FALSE
)

Type III Test

Sum.Sq	Mean.Sq	NumDF	DenDF	F.value	Pr..F.
8.55	2.85	3	159.65	34.95	<0.001
184.38	92.19	2	97.84	1130.62	<0.001
34.12	8.53	4	216.07	104.60	<0.001
102.09	17.02	6	2306.10	208.68	<0.001
7.91	0.66	12	2181.90	8.09	<0.001
38.69	4.84	8	2178.57	59.32	<0.001
11.79	0.49	24	2174.56	6.03	<0.001

cat('</div>')

Post-hoc Tests

emm_occurrence <- emmeans(occurrence_model, ~ cycle | stage * class)

contrast_means_occurrence <- as.data.frame(pairs(emm_occurrence, adjust = "holm")) |>
  dplyr::select(stage, class, contrast, estimate, SE, df, t.ratio, p.value) |>
  dplyr::mutate(
    stage    = factor(stage, levels = c("N1", "N2", "N3")),
    estimate = round(estimate, 2),
    SE       = round(SE, 2),
    df       = round(df, 2),
    t.ratio  = round(t.ratio, 2)
  ) |>
  dplyr::arrange(stage, class, contrast)

for (stg in c("N1", "N2", "N3")) {

  tab <- contrast_means_occurrence |>
    dplyr::filter(stage == stg) |>
    dplyr::select(-stage)

  p_raw <- tab$p.value

  tab$p.value <- ifelse(
    p_raw < .001,
    "<b>&lt;0.001</b>",
    ifelse(
      p_raw < .05,
      paste0("<b>", sprintf("%.3f", p_raw), "</b>"),
      sprintf("%.3f", p_raw)
    )
  )

  colnames(tab) <- c("Class", "Cycle contrast", "Estimate", "SE", "df", "t", "p")

  cat('<div class="scroll-y">')
  cat(paste0("<h4>Occurrence: ", stg, "</h4>"))

  tab_out <- sjPlot::tab_df(
    tab,
    title = NULL,
    digits = 2,
    encoding = "UTF-8",
    escape = FALSE
  )

  cat(tab_out$knitr)
  cat("</div><br>")
}

Occurrence: N1

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	1.23	0.07	1401.41	17.04	<0.001
A	cycle1 - cycle3	1.02	0.07	1467.62	13.68	<0.001
A	cycle1 - cycle4	0.81	0.08	1517.22	10.59	<0.001
A	cycle2 - cycle3	-0.21	0.08	1607.29	-2.73	0.013
A	cycle2 - cycle4	-0.42	0.08	1644.31	-5.24	<0.001
A	cycle3 - cycle4	-0.21	0.08	1662.87	-2.55	0.013
B	cycle1 - cycle2	1.33	0.07	1401.41	18.35	<0.001
B	cycle1 - cycle3	1.12	0.07	1467.62	15.05	<0.001
B	cycle1 - cycle4	1.00	0.08	1517.22	13.03	<0.001
B	cycle2 - cycle3	-0.21	0.08	1607.29	-2.62	0.018
B	cycle2 - cycle4	-0.33	0.08	1644.31	-4.09	<0.001
B	cycle3 - cycle4	-0.12	0.08	1662.87	-1.51	0.131
C	cycle1 - cycle2	0.82	0.07	1401.41	11.25	<0.001
C	cycle1 - cycle3	0.62	0.07	1467.62	8.35	<0.001
C	cycle1 - cycle4	0.40	0.08	1517.22	5.18	<0.001
C	cycle2 - cycle3	-0.19	0.08	1607.29	-2.44	0.015
C	cycle2 - cycle4	-0.42	0.08	1644.31	-5.19	<0.001
C	cycle3 - cycle4	-0.23	0.08	1662.87	-2.77	0.011
D	cycle1 - cycle2	0.53	0.07	1401.41	7.32	<0.001
D	cycle1 - cycle3	0.16	0.07	1467.62	2.14	0.098
D	cycle1 - cycle4	0.02	0.08	1517.22	0.27	0.784
D	cycle2 - cycle3	-0.37	0.08	1607.29	-4.72	<0.001
D	cycle2 - cycle4	-0.51	0.08	1644.31	-6.32	<0.001
D	cycle3 - cycle4	-0.14	0.08	1662.87	-1.70	0.180
E	cycle1 - cycle2	0.94	0.07	1401.41	12.91	<0.001
E	cycle1 - cycle3	0.67	0.07	1467.62	9.02	<0.001
E	cycle1 - cycle4	0.52	0.08	1517.22	6.79	<0.001
E	cycle2 - cycle3	-0.26	0.08	1607.29	-3.34	0.002
E	cycle2 - cycle4	-0.41	0.08	1644.31	-5.14	<0.001
E	cycle3 - cycle4	-0.15	0.08	1662.87	-1.86	0.063

Occurrence: N2

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	0.11	0.06	1075.54	1.74	0.328
A	cycle1 - cycle3	-0.01	0.06	1075.57	-0.12	1.000
A	cycle1 - cycle4	-0.04	0.07	1075.28	-0.58	1.000
A	cycle2 - cycle3	-0.12	0.06	1075.57	-1.85	0.322
A	cycle2 - cycle4	-0.15	0.07	1075.28	-2.30	0.129
A	cycle3 - cycle4	-0.03	0.07	1075.25	-0.46	1.000
B	cycle1 - cycle2	0.12	0.06	1075.54	1.91	0.283
B	cycle1 - cycle3	0.01	0.06	1075.57	0.11	0.989
B	cycle1 - cycle4	-0.06	0.07	1075.28	-0.87	0.989
B	cycle2 - cycle3	-0.12	0.06	1075.57	-1.79	0.297
B	cycle2 - cycle4	-0.18	0.07	1075.28	-2.75	0.036
B	cycle3 - cycle4	-0.06	0.07	1075.25	-0.98	0.989
C	cycle1 - cycle2	0.11	0.06	1075.54	1.65	0.299
C	cycle1 - cycle3	-0.07	0.06	1075.57	-1.06	0.299
C	cycle1 - cycle4	-0.18	0.07	1075.28	-2.70	0.035
C	cycle2 - cycle3	-0.17	0.06	1075.57	-2.69	0.035
C	cycle2 - cycle4	-0.28	0.07	1075.28	-4.33	<0.001
C	cycle3 - cycle4	-0.11	0.07	1075.25	-1.65	0.299
D	cycle1 - cycle2	0.16	0.06	1075.54	2.51	0.049
D	cycle1 - cycle3	0.00	0.06	1075.57	0.06	0.951
D	cycle1 - cycle4	-0.17	0.07	1075.28	-2.53	0.049
D	cycle2 - cycle3	-0.16	0.06	1075.57	-2.43	0.049
D	cycle2 - cycle4	-0.33	0.07	1075.28	-5.02	<0.001
D	cycle3 - cycle4	-0.17	0.07	1075.25	-2.59	0.049
E	cycle1 - cycle2	0.08	0.06	1075.54	1.23	0.340
E	cycle1 - cycle3	-0.10	0.06	1075.57	-1.50	0.340
E	cycle1 - cycle4	-0.20	0.07	1075.28	-3.08	0.011
E	cycle2 - cycle3	-0.18	0.06	1075.57	-2.71	0.027
E	cycle2 - cycle4	-0.28	0.07	1075.28	-4.29	<0.001
E	cycle3 - cycle4	-0.10	0.07	1075.25	-1.59	0.340

Occurrence: N3

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	-0.20	0.07	1121.82	-3.10	0.004
A	cycle1 - cycle3	-0.29	0.07	1134.18	-4.31	<0.001
A	cycle1 - cycle4	-0.58	0.07	1334.75	-8.10	<0.001
A	cycle2 - cycle3	-0.08	0.07	1158.64	-1.23	0.219
A	cycle2 - cycle4	-0.38	0.07	1358.23	-5.23	<0.001
A	cycle3 - cycle4	-0.30	0.07	1368.01	-4.06	<0.001
B	cycle1 - cycle2	-0.21	0.07	1121.82	-3.13	0.004
B	cycle1 - cycle3	-0.30	0.07	1134.18	-4.54	<0.001
B	cycle1 - cycle4	-0.55	0.07	1334.75	-7.59	<0.001
B	cycle2 - cycle3	-0.10	0.07	1158.64	-1.43	0.152
B	cycle2 - cycle4	-0.34	0.07	1358.23	-4.70	<0.001
B	cycle3 - cycle4	-0.25	0.07	1368.01	-3.35	0.002
C	cycle1 - cycle2	-0.17	0.07	1121.82	-2.63	0.017
C	cycle1 - cycle3	-0.34	0.07	1134.18	-5.19	<0.001
C	cycle1 - cycle4	-0.63	0.07	1334.75	-8.75	<0.001
C	cycle2 - cycle3	-0.17	0.07	1158.64	-2.57	0.017
C	cycle2 - cycle4	-0.46	0.07	1358.23	-6.30	<0.001
C	cycle3 - cycle4	-0.29	0.07	1368.01	-3.90	<0.001
D	cycle1 - cycle2	-0.16	0.07	1121.82	-2.44	0.029
D	cycle1 - cycle3	-0.27	0.07	1134.18	-4.06	<0.001
D	cycle1 - cycle4	-0.63	0.07	1334.75	-8.76	<0.001
D	cycle2 - cycle3	-0.11	0.07	1158.64	-1.63	0.104
D	cycle2 - cycle4	-0.47	0.07	1358.23	-6.48	<0.001
D	cycle3 - cycle4	-0.36	0.07	1368.01	-4.94	<0.001
E	cycle1 - cycle2	-0.20	0.07	1121.82	-2.98	0.006
E	cycle1 - cycle3	-0.38	0.07	1134.18	-5.68	<0.001
E	cycle1 - cycle4	-0.64	0.07	1334.75	-8.88	<0.001
E	cycle2 - cycle3	-0.18	0.07	1158.64	-2.72	0.007
E	cycle2 - cycle4	-0.44	0.07	1358.23	-6.12	<0.001
E	cycle3 - cycle4	-0.26	0.07	1368.01	-3.58	0.001

Plot

occurrence_model_emm <- emmeans(occurrence_model, ~ cycle * stage | class)
occurrence_model_emm <- as.data.frame(occurrence_model_emm)

occurrence_model_emm_plot <- occurrence_model_emm |> 
  select(class, cycle, stage, emmean, lower.CL, upper.CL) |> 
  slice(1:60) |> 
  ggplot(aes(group = stage, color = stage, y = emmean, x = cycle)) +
  geom_line(fill = cycle, alpha = 1) +
  geom_point(stat="identity", size = 3, alpha = 1) +
  geom_errorbar(aes(ymin=lower.CL, ymax=upper.CL), width=.2, alpha = 1) +
  ggtitle("Microstate Classes") +
  facet_wrap( ~ class, ncol = 5) +
  theme_bw() +
  theme(
    plot.title = element_text(size = size_big, hjust= 0.5),
    axis.text = element_text(size = size_small, family = "sans"),
    axis.title = element_text(size = size_big, family = "sans"),
    strip.text.x = element_text(size = size_big, family = "sans"),
    strip.background = element_blank(),
    legend.title = element_blank(),
    legend.text = element_text(size = size_small),
    legend.position="right",
  ) +
  labs(color = "Sleep Stages") +
  scale_y_continuous(name = "Occurence (per second)", limits = c(0,6), expand = c(0,0)) +
  scale_x_discrete(name= "Sleep Cycles") +
  scale_color_manual(values = c("#1b9e77","#d95f02","#7570b3")) #color-blind friendly colors

occurrence_model_emm_plot

ggsave("occurrence_plot.png",
        plot = occurrence_model_emm_plot,
        width = 9,
        height = 7,         
        dpi = 1200,
        bg = "transparent",
        path = "R:/MicrostateAnalysis/sleepData/MS_SleepMarkers/scripts/graphs/sleep/cycles")

Coverage

Model

coverage_model <- lmer(
coverage ~ cycle * stage * class +
(1|ID) + (1|cycle:ID) + (1|stage:ID) + (1|class:ID),
data = data_long_coverage, REML = TRUE
)

cat('<div class="scroll-y">')

sjPlot::tab_model(coverage_model)

	coverage
Predictors	Estimates	CI	p
(Intercept)	22.24	21.09 – 23.39	<0.001
cycle [2]	-2.98	-3.79 – -2.18	<0.001
cycle [3]	-2.94	-3.76 – -2.11	<0.001
cycle [4]	-2.63	-3.49 – -1.78	<0.001
stage [N2]	-3.85	-4.57 – -3.12	<0.001
stage [N3]	-4.34	-5.07 – -3.61	<0.001
class [B]	-1.95	-3.57 – -0.32	0.019
class [C]	2.17	0.54 – 3.79	0.009
class [D]	-4.66	-6.28 – -3.04	<0.001
class [E]	-6.76	-8.39 – -5.14	<0.001
cycle [2] × stage [N2]	2.90	1.83 – 3.97	<0.001
cycle [3] × stage [N2]	2.75	1.66 – 3.83	<0.001
cycle [4] × stage [N2]	1.94	0.83 – 3.04	0.001
cycle [2] × stage [N3]	3.29	2.21 – 4.36	<0.001
cycle [3] × stage [N3]	3.08	1.98 – 4.17	<0.001
cycle [4] × stage [N3]	2.91	1.75 – 4.07	<0.001
cycle [2] × class [B]	-0.59	-1.73 – 0.54	0.307
cycle [3] × class [B]	-0.48	-1.65 – 0.69	0.426
cycle [4] × class [B]	-0.80	-2.01 – 0.40	0.192
cycle [2] × class [C]	7.42	6.29 – 8.56	<0.001
cycle [3] × class [C]	5.35	4.18 – 6.52	<0.001
cycle [4] × class [C]	5.58	4.37 – 6.78	<0.001
cycle [2] × class [D]	5.96	4.83 – 7.10	<0.001
cycle [3] × class [D]	7.22	6.05 – 8.39	<0.001
cycle [4] × class [D]	6.32	5.12 – 7.53	<0.001
cycle [2] × class [E]	2.12	0.98 – 3.26	<0.001
cycle [3] × class [E]	2.60	1.43 – 3.77	<0.001
cycle [4] × class [E]	2.07	0.87 – 3.28	0.001
stage [N2] × class [B]	-0.85	-1.87 – 0.18	0.105
stage [N3] × class [B]	-1.57	-2.60 – -0.55	0.003
stage [N2] × class [C]	5.49	4.46 – 6.51	<0.001
stage [N3] × class [C]	9.37	8.34 – 10.40	<0.001
stage [N2] × class [D]	8.38	7.35 – 9.40	<0.001
stage [N3] × class [D]	8.72	7.69 – 9.75	<0.001
stage [N2] × class [E]	6.23	5.20 – 7.25	<0.001
stage [N3] × class [E]	5.19	4.16 – 6.22	<0.001
(cycle [2] × stage [N2]) × class [B]	0.41	-1.10 – 1.91	0.596
(cycle [3] × stage [N2]) × class [B]	0.41	-1.12 – 1.95	0.597
(cycle [4] × stage [N2]) × class [B]	1.12	-0.44 – 2.69	0.161
(cycle [2] × stage [N3]) × class [B]	0.91	-0.61 – 2.43	0.239
(cycle [3] × stage [N3]) × class [B]	1.02	-0.53 – 2.57	0.198
(cycle [4] × stage [N3]) × class [B]	0.90	-0.74 – 2.54	0.281
(cycle [2] × stage [N2]) × class [C]	-6.88	-8.39 – -5.38	<0.001
(cycle [3] × stage [N2]) × class [C]	-4.91	-6.45 – -3.38	<0.001
(cycle [4] × stage [N2]) × class [C]	-4.82	-6.38 – -3.25	<0.001
(cycle [2] × stage [N3]) × class [C]	-8.60	-10.12 – -7.08	<0.001
(cycle [3] × stage [N3]) × class [C]	-6.37	-7.92 – -4.82	<0.001
(cycle [4] × stage [N3]) × class [C]	-8.17	-9.81 – -6.53	<0.001
(cycle [2] × stage [N2]) × class [D]	-6.20	-7.70 – -4.69	<0.001
(cycle [3] × stage [N2]) × class [D]	-7.50	-9.04 – -5.97	<0.001
(cycle [4] × stage [N2]) × class [D]	-5.43	-6.99 – -3.86	<0.001
(cycle [2] × stage [N3]) × class [D]	-6.81	-8.33 – -5.28	<0.001
(cycle [3] × stage [N3]) × class [D]	-8.55	-10.10 – -7.00	<0.001
(cycle [4] × stage [N3]) × class [D]	-6.48	-8.12 – -4.84	<0.001
(cycle [2] × stage [N2]) × class [E]	-1.83	-3.33 – -0.32	0.018
(cycle [3] × stage [N2]) × class [E]	-1.73	-3.27 – -0.20	0.027
(cycle [4] × stage [N2]) × class [E]	-0.56	-2.12 – 1.01	0.486
(cycle [2] × stage [N3]) × class [E]	-1.95	-3.48 – -0.43	0.012
(cycle [3] × stage [N3]) × class [E]	-1.49	-3.04 – 0.06	0.060
(cycle [4] × stage [N3]) × class [E]	-0.81	-2.46 – 0.83	0.331
Random Effects
σ2	3.44
τ00 class:ID	14.60
τ00 cycle:ID	0.00
τ00 stage:ID	0.00
τ00 ID	0.00
N ID	54
N cycle	4
N stage	3
N class	5
Observations	2760
Marginal R2 / Conditional R2	0.839 / NA

cat('</div>')

Type III Test

anova_coverage <- as.data.frame(anova(coverage_model))

p_raw <- anova_coverage$`Pr(>F)`

anova_coverage$`Pr(>F)` <- ifelse(
  p_raw < .001,
  "<b>&lt;0.001</b>",
  ifelse(
    p_raw < .05,
    paste0("<b>", sprintf("%.3f", p_raw), "</b>"),
    sprintf("%.3f", p_raw)
  )
)

cat('<div class="scroll-y">')

sjPlot::tab_df(
  anova_coverage,
  title = "Type III Test",
  digits = 2,
  encoding = "UTF-8",
  escape = FALSE
)

Type III Test

Sum.Sq	Mean.Sq	NumDF	DenDF	F.value	Pr..F.
0.00	0.00	3	2168.75	0.00	1.000
0.00	0.00	2	2165.43	0.00	1.000
1007.75	251.94	4	15.40	73.15	<0.001
0.00	0.00	6	2161.67	0.00	1.000
539.12	44.93	12	2168.75	13.04	<0.001
2740.42	342.55	8	2165.43	99.46	<0.001
1424.88	59.37	24	2161.67	17.24	<0.001

cat('</div>')

Post-hoc Tests

emm_coverage <- emmeans(coverage_model, ~ cycle | stage * class)

contrast_means_coverage <- as.data.frame(pairs(emm_coverage, adjust = "holm")) |>
  dplyr::select(stage, class, contrast, estimate, SE, df, t.ratio, p.value) |>
  dplyr::mutate(
    stage    = factor(stage, levels = c("N1", "N2", "N3")),
    estimate = round(estimate, 2),
    SE       = round(SE, 2),
    df       = round(df, 2),
    t.ratio  = round(t.ratio, 2)
  ) |>
  dplyr::arrange(stage, class, contrast)

for (stg in c("N1", "N2", "N3")) {

  tab <- contrast_means_coverage |>
    dplyr::filter(stage == stg) |>
    dplyr::select(-stage)

  p_raw <- tab$p.value

  tab$p.value <- ifelse(
    p_raw < .001,
    "<b>&lt;0.001</b>",
    ifelse(
      p_raw < .05,
      paste0("<b>", sprintf("%.3f", p_raw), "</b>"),
      sprintf("%.3f", p_raw)
    )
  )

  colnames(tab) <- c("Class", "Cycle contrast", "Estimate", "SE", "df", "t", "p")

  cat('<div class="scroll-y">')
  cat(paste0("<h4>Occurrence: ", stg, "</h4>"))

  tab_out <- sjPlot::tab_df(
    tab,
    title = NULL,
    digits = 2,
    encoding = "UTF-8",
    escape = FALSE
  )

  cat(tab_out$knitr)
  cat("</div><br>")
}

Occurrence: N1

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	2.98	0.41	2371.80	7.28	<0.001
A	cycle1 - cycle3	2.94	0.42	2391.30	6.96	<0.001
A	cycle1 - cycle4	2.63	0.43	2397.23	6.06	<0.001
A	cycle2 - cycle3	-0.05	0.45	2399.67	-0.10	1.000
A	cycle2 - cycle4	-0.35	0.46	2401.62	-0.76	1.000
A	cycle3 - cycle4	-0.30	0.47	2390.40	-0.65	1.000
B	cycle1 - cycle2	3.58	0.41	2371.80	8.72	<0.001
B	cycle1 - cycle3	3.41	0.42	2391.30	8.09	<0.001
B	cycle1 - cycle4	3.44	0.43	2397.23	7.90	<0.001
B	cycle2 - cycle3	-0.16	0.45	2399.67	-0.36	1.000
B	cycle2 - cycle4	-0.14	0.46	2401.62	-0.30	1.000
B	cycle3 - cycle4	0.02	0.47	2390.40	0.05	1.000
C	cycle1 - cycle2	-4.44	0.41	2371.80	-10.83	<0.001
C	cycle1 - cycle3	-2.41	0.42	2391.30	-5.71	<0.001
C	cycle1 - cycle4	-2.94	0.43	2397.23	-6.77	<0.001
C	cycle2 - cycle3	2.03	0.45	2399.67	4.52	<0.001
C	cycle2 - cycle4	1.50	0.46	2401.62	3.25	0.002
C	cycle3 - cycle4	-0.53	0.47	2390.40	-1.14	0.256
D	cycle1 - cycle2	-2.98	0.41	2371.80	-7.28	<0.001
D	cycle1 - cycle3	-4.28	0.42	2391.30	-10.14	<0.001
D	cycle1 - cycle4	-3.69	0.43	2397.23	-8.49	<0.001
D	cycle2 - cycle3	-1.30	0.45	2399.67	-2.89	0.012
D	cycle2 - cycle4	-0.71	0.46	2401.62	-1.54	0.248
D	cycle3 - cycle4	0.59	0.47	2390.40	1.26	0.248
E	cycle1 - cycle2	0.86	0.41	2371.80	2.10	0.213
E	cycle1 - cycle3	0.34	0.42	2391.30	0.81	1.000
E	cycle1 - cycle4	0.56	0.43	2397.23	1.29	0.980
E	cycle2 - cycle3	-0.52	0.45	2399.67	-1.16	0.980
E	cycle2 - cycle4	-0.30	0.46	2401.62	-0.65	1.000
E	cycle3 - cycle4	0.22	0.47	2390.40	0.48	1.000

Occurrence: N2

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	0.08	0.36	2328.39	0.23	1.000
A	cycle1 - cycle3	0.19	0.36	2330.07	0.53	1.000
A	cycle1 - cycle4	0.70	0.36	2331.54	1.93	0.320
A	cycle2 - cycle3	0.11	0.36	2330.07	0.30	1.000
A	cycle2 - cycle4	0.61	0.36	2331.54	1.70	0.444
A	cycle3 - cycle4	0.51	0.36	2329.88	1.40	0.646
B	cycle1 - cycle2	0.27	0.36	2328.39	0.75	1.000
B	cycle1 - cycle3	0.25	0.36	2330.07	0.70	1.000
B	cycle1 - cycle4	0.38	0.36	2331.54	1.05	1.000
B	cycle2 - cycle3	-0.02	0.36	2330.07	-0.05	1.000
B	cycle2 - cycle4	0.11	0.36	2331.54	0.31	1.000
B	cycle3 - cycle4	0.13	0.36	2329.88	0.36	1.000
C	cycle1 - cycle2	-0.45	0.36	2328.39	-1.27	1.000
C	cycle1 - cycle3	-0.25	0.36	2330.07	-0.69	1.000
C	cycle1 - cycle4	-0.06	0.36	2331.54	-0.17	1.000
C	cycle2 - cycle3	0.21	0.36	2330.07	0.58	1.000
C	cycle2 - cycle4	0.39	0.36	2331.54	1.09	1.000
C	cycle3 - cycle4	0.18	0.36	2329.88	0.51	1.000
D	cycle1 - cycle2	0.31	0.36	2328.39	0.88	1.000
D	cycle1 - cycle3	0.48	0.36	2330.07	1.33	0.774
D	cycle1 - cycle4	-0.20	0.36	2331.54	-0.55	1.000
D	cycle2 - cycle3	0.16	0.36	2330.07	0.45	1.000
D	cycle2 - cycle4	-0.51	0.36	2331.54	-1.42	0.774
D	cycle3 - cycle4	-0.68	0.36	2329.88	-1.87	0.373
E	cycle1 - cycle2	-0.21	0.36	2328.39	-0.59	1.000
E	cycle1 - cycle3	-0.67	0.36	2330.07	-1.87	0.305
E	cycle1 - cycle4	-0.82	0.36	2331.54	-2.26	0.142
E	cycle2 - cycle3	-0.46	0.36	2330.07	-1.28	0.598
E	cycle2 - cycle4	-0.61	0.36	2331.54	-1.68	0.375
E	cycle3 - cycle4	-0.14	0.36	2329.88	-0.40	1.000

Occurrence: N3

Class	Cycle.contrast	Estimate	SE	df	t	p
A	cycle1 - cycle2	-0.31	0.36	2336.95	-0.84	1.000
A	cycle1 - cycle3	-0.14	0.37	2340.13	-0.38	1.000
A	cycle1 - cycle4	-0.28	0.41	2369.63	-0.69	1.000
A	cycle2 - cycle3	0.17	0.37	2345.04	0.44	1.000
A	cycle2 - cycle4	0.03	0.41	2375.47	0.07	1.000
A	cycle3 - cycle4	-0.14	0.41	2375.93	-0.34	1.000
B	cycle1 - cycle2	-0.63	0.36	2336.95	-1.72	0.428
B	cycle1 - cycle3	-0.68	0.37	2340.13	-1.85	0.385
B	cycle1 - cycle4	-0.38	0.41	2369.63	-0.93	1.000
B	cycle2 - cycle3	-0.06	0.37	2345.04	-0.15	1.000
B	cycle2 - cycle4	0.25	0.41	2375.47	0.61	1.000
B	cycle3 - cycle4	0.30	0.41	2375.93	0.74	1.000
C	cycle1 - cycle2	0.87	0.36	2336.95	2.38	0.051
C	cycle1 - cycle3	0.88	0.37	2340.13	2.39	0.051
C	cycle1 - cycle4	2.31	0.41	2369.63	5.71	<0.001
C	cycle2 - cycle3	0.01	0.37	2345.04	0.03	0.974
C	cycle2 - cycle4	1.45	0.41	2375.47	3.53	0.002
C	cycle3 - cycle4	1.43	0.41	2375.93	3.47	0.002
D	cycle1 - cycle2	0.54	0.36	2336.95	1.47	0.325
D	cycle1 - cycle3	1.19	0.37	2340.13	3.23	0.008
D	cycle1 - cycle4	-0.12	0.41	2369.63	-0.30	0.764
D	cycle2 - cycle3	0.66	0.37	2345.04	1.76	0.316
D	cycle2 - cycle4	-0.66	0.41	2375.47	-1.61	0.325
D	cycle3 - cycle4	-1.31	0.41	2375.93	-3.18	0.008
E	cycle1 - cycle2	-0.47	0.36	2336.95	-1.29	0.392
E	cycle1 - cycle3	-1.25	0.37	2340.13	-3.38	0.004
E	cycle1 - cycle4	-1.54	0.41	2369.63	-3.79	<0.001
E	cycle2 - cycle3	-0.78	0.37	2345.04	-2.08	0.112
E	cycle2 - cycle4	-1.06	0.41	2375.47	-2.60	0.037
E	cycle3 - cycle4	-0.29	0.41	2375.93	-0.70	0.486

Plot

coverage_model_emm <- emmeans(coverage_model, ~ cycle * stage | class)
coverage_model_emm <- as.data.frame(coverage_model_emm)

coverage_model_emm_plot <- coverage_model_emm |> 
  select(class, cycle, stage, emmean, lower.CL, upper.CL) |> 
  slice(1:60) |> 
  ggplot(aes(group = stage, color = stage, y = emmean, x = cycle)) +
  geom_line(fill = cycle, alpha = 1) +
  geom_point(stat="identity", size = 3, alpha = 1) +
  geom_errorbar(aes(ymin=lower.CL, ymax=upper.CL), width=.2, alpha = 1) +
  ggtitle("Microstate Classes") +
  facet_wrap(~ class, ncol = 5) +
  theme_bw() +
   theme(
    plot.title = element_text(size = size_big, hjust= 0.5),
    axis.text = element_text(size = size_small, family = "sans"),
    axis.title = element_text(size = size_big, family = "sans"),
    strip.text.x = element_text(size = size_big, family = "sans"),
    strip.background = element_blank(),
    legend.title = element_blank(),
    legend.text = element_text(size = size_small),
    legend.position="right",
  ) +
  labs(color = "Sleep Stages") +
  scale_y_continuous(name = "Coverage (%)", limits = c(0,35), expand = c(0,0)) +
  scale_x_discrete(name= "Sleep Cycles") +
  scale_color_manual(values = c("#1b9e77","#d95f02","#7570b3")) #color-blind friendly colors

coverage_model_emm_plot

ggsave("coverage_plot.png",
        plot = coverage_model_emm_plot,
        width = 9,
        height = 7,
        dpi = 1200,
        bg = "transparent",
        path = "R:/MicrostateAnalysis/sleepData/MS_SleepMarkers/scripts/graphs/sleep/cycles")