Appendix G — Labor Duration Analysis

G.1 Set the Environment

G.1.1 Load Packages

library(actverse) # github.com/danielvartan/actverse
library(askpass)
library(car)
library(checkmate)
library(cli)
library(dplyr)
library(fs)
library(ggplot2)
library(googlesheets4)
library(grDevices)
library(here)
library(hms)
library(lubridate)
library(lubritime) # github.com/danielvartan/lubritime
library(readr)
library(orbis) # github.com/danielvartan/orbis
library(plotr) # github.com/danielvartan/plotr
library(purrr)
library(rutils) # github.com/danielvartan/rutils
library(scaler) # github.com/danielvartan/scaler
library(scales)
library(sodium)
library(stats)
library(stringr)
library(tidyr)

G.1.2 Load Custom Functions

source(here("R", "anonymize_id.R"))
source(here("R", "ga.R"))

get_n <- function(data, group) {
  assert_tibble(data)
  assert_subset("group", names(data))
  assert_set_equal(levels(data$group), c("Poor", "Good"))
  assert_choice(group, c("Poor", "Good"))

  data |>
    rename(.group = group) |>
    filter(.group == group) |>
    nrow()
}

get_labor_duration <- function(data, group) {
  assert_tibble(data)
  assert_subset(c("group", "labor_duration"), names(data))
  assert_set_equal(levels(data$group), c("Poor", "Good"))
  assert_choice(group, c("Poor", "Good"))

  data |>
    rename(.group = group) |>
    filter(.group == group) |>
    pull(labor_duration)
}

# See Maxwell et al. (2018, p. 150)

f_unb <- function(k, n_total, f_test) {
  assert_int(k, lower = 2)
  assert_int(n_total, lower = 2)
  assert_number(f_test)

  ((k - 1) / n_total) %>%
  `*`(((n_total - k - 2) / (n_total - k)) * (f_test - 1)) |>
  sqrt()
}

G.1.3 Set Data Directories

data_dir <- here("data")
processed_data_dir <- here(data_dir, "processed")
actigraphy_processed_data_dir <- here(data_dir, "processed", "actigraphy")
sri_data_dir <- here(data_dir, "sri")
sleep_prop_dir <- here(data_dir, "sleep-prop")
lux_data_dir <- here(data_dir, "lux")

for (i in ls(pattern = "_dir$")) {
  if (!dir_exists(get(i))) dir_create(get(i), recurse = TRUE)
}

G.1.4 Set Keys

gs4_auth(
  cache = here(".secrets"),
  email = TRUE,
  use_oob = FALSE
)

salt <- Sys.getenv("PREGNANCY_SALT") # askpass()

G.2 Download Actigraphy Processed Data from OSF

See the osf-routines.qmd Quarto notebook to download and decrypt the actigraphy processed data from OSF.

G.3 Download and Import the Control Form Data

control_form_data <-
  "13wtDr4fRD1wJSM-qdLOdGSchF8oXU1bG0Jtccu24GhY" |>
  read_sheet(
    sheet = "Dataset",
    col_types = "c",
    na = c("", "NA")
  )

G.4 Download and Import the Field Form Data

field_form_data <-
  "1tY_TT0nPXFsErYQS2VED0-hqTzgHudA9BfhRhYG19fw" |>
  read_sheet(
    sheet = "Dataset",
    col_types = "c",
    na = c("", "NA")
  )

G.5 Download and Import the Primary Sleep Data

"18wKHeqbNfrQDzjyOKIX5MCgjVK0WMho6xvdQQpA7ZLE" |>
  sheet_names() %>%
  magrittr::extract(
    !. %in% c(
      c(
        "Documentation",
        "Codebook",
        "Validation",
        "Template",
        "[TEMPLATE]"
      )
    )
  ) |>
  unique()

sheet_names

primary_sleep_data <- tibble()

for (i in sheet_names) {
  i_data <-
    "18wKHeqbNfrQDzjyOKIX5MCgjVK0WMho6xvdQQpA7ZLE" |>
    read_sheet(
      sheet = i,
      col_types = "c",
      na = c("", "NA")
    ) |>
      drop_na(date)

  if (nrow(i_data) == 0) next

  i_data <-
    i_data |>
    mutate(id = i) |>
    relocate(id)

  primary_sleep_data <-
    primary_sleep_data |>
    bind_rows(i_data)
}

G.6 Download and Import the Secondary Sleep Data

sheet_names <-
  "12xFtuHk3dBB6KqU8qhLS17ZcrZqs0AjGWaQtnUE82Uk" |>
  sheet_names() %>%
  magrittr::extract(
    !. %in% c(
      c(
        "Documentation",
        "Codebook",
        "Validation",
        "Template",
        "[TEMPLATE]"
      )
    )
  ) |>
  unique()

sheet_names

secondary_sleep_data <- tibble()

for (i in sheet_names) {
  i_data <-
    "12xFtuHk3dBB6KqU8qhLS17ZcrZqs0AjGWaQtnUE82Uk" |>
    read_sheet(
      sheet = i,
      col_types = "c",
      na = c("", "NA")
    ) |>
      drop_na(date)

  if (nrow(i_data) == 0) next

  i_data <-
    i_data |>
    mutate(id = i) |>
    relocate(id)

  secondary_sleep_data <-
    secondary_sleep_data |>
    bind_rows(i_data)
}

G.7 Download and Import the Medical Records Data

medical_record_data <-
  "10yvvTOCjZsRKPO8A4yn0NCvWe4Jn_Kq-Zu_CV4miGJI" |>
  read_sheet(
    sheet = "Dataset",
    col_types = "c",
    na = c("", "NA")
  )

G.8 Filter and Tidy the Control Form Data

control_form_data <-
  control_form_data |>
  `names<-`(paste0("X", seq_along(control_form_data))) |>
  dplyr::select(X5, X6, X12, X13, X14, X39, X40) |>
  rename(
    cpf = X5,
    id = X6,
    birth_center_delivery = X12,
    birth = X13,
    labor_start = X14,
    ultrasound_ga = X39,
    ultrasound_date = X40
  ) |>
  mutate(
    birth_center_delivery = case_when(
      tolower(trimws(birth_center_delivery)) == "sim" ~ TRUE,
      tolower(trimws(birth_center_delivery)) == "não" ~ FALSE,
      TRUE ~ NA
    ),
    ultrasound_date = ymd(ultrasound_date),
    ga_start = if_else(
      is.na(ultrasound_ga),
      as.Date(NA),
      ga_start(
        ultrasound = ultrasound_date,
        ga =
          str_extract(ultrasound_ga, "^\\d{1,2}") |>
          as.integer() |>
          weeks() |>
          add(
            str_extract(ultrasound_ga, "\\d(?= di)") |>
            as.integer() |>
            days()
          )
      )
    ),
    birth = as_datetime(birth, tz = "America/Sao_Paulo"),
    labor_start = as_datetime(labor_start, tz = "America/Sao_Paulo"),
    labor_duration =
      interval(
        start = labor_start,
        end = birth
      ) |>
        lubritime::int_duration()
  ) |>
  dplyr::select(-starts_with("cpf")) |>
  dplyr::select(
    id, birth_center_delivery, ga_start, labor_start, labor_duration
  )

control_form_data |> glimpse()

G.9 Filter and Tidy the Field Form Data

field_form_data <-
  field_form_data |>
  `names<-`(paste0("X", seq_along(field_form_data))) |>
  dplyr::select(
    X1, X4, X5, X11, X91, X92, X95, X96, X97, X99, X101, X102, X106, X119,
    X118, X125, X126
  ) |>
  rename(
    timestamp = X1,
    birth_date = X4,
    cpf = X5,
    weight_before = X96,
    weight_current = X95,
    height = X97,
    ethnicity = X106,
    education = X118,
    family_income = X99,
    health_plan = X101,
    solo = X102,
    exercise = X119,
    work = X92,
    study = X91,
    gestations = X11,
    deliveries = X125,
    abortions = X126
  ) |>
  mutate(
    id = anonymize_id(cpf, salt),
    timestamp = mdy_hms(timestamp, tz = "America/Sao_Paulo"),
    birth_date = dmy(birth_date),
    age = scaler:::age(birth_date, timestamp),
    ethnicity = factor(
      ethnicity,
      levels = c(
        "Indígena",
        "Preta",
        "Parda",
        "Amarela",
        "Branca"
      ),
      ordered = TRUE
    ),
    education = factor(
      education,
      levels = c(
        "Não frequentou a escola",
        "Fundamental incompleto",
        "Fundamental completo",
        "Ensino médio incompleto",
        "Ensino médio completo",
        "Ensino superior incompleto",
        "Ensino superior completo",
        "Mestrado incompleto",
        "Mestrado completo",
        "Doutorado incompleto",
        "Doutorado completo",
        "Pós-doutorado incompleto",
        "Pós-doutorado completo"
      ),
      ordered = TRUE
    ),
    gestations = case_match(
      gestations,
      c("Sim", "0") ~ "1",
      "Não" ~ "0",
      .default = gestations
    ),
    deliveries = case_when(
      is.na(deliveries) ~ "0",
      TRUE ~ deliveries
    ),
    abortions = case_when(
      is.na(abortions) ~ "0",
      TRUE ~ abortions
    )
  ) |>
  mutate(
    across(
      .cols = all_of(
        c(
          "weight_before", "weight_current", "height",
          "family_income", "gestations", "deliveries",
          "abortions"
        )
      ),
      .fns = as.numeric
    )
  ) |>
  mutate(
    across(
      .cols = all_of(
        c("health_plan", "solo", "exercise", "work", "study")
      ),
      .fns = function(x) {
        case_match(
          x,
          "Sim" ~ TRUE,
          "Não" ~ FALSE
        )
      }
    )
  ) |>
  mutate(
    gestations = case_when(
      id == "ae5fc5cd" ~ 3,
      TRUE ~ gestations
    ),
  ) |>
  mutate(
    bmi_before = weight_before / ((height / 100)^2),
    bmi_current = weight_current / ((height / 100)^2)
  ) |>
  dplyr::select(-cpf) |>
  relocate(
    id, timestamp, birth_date, age, weight_before, weight_current, height, bmi_before, bmi_current, family_income
  )

field_form_data |> glimpse()

G.10 Filter and Tidy the Primary Sleep Data

primary_sleep_data <-
  primary_sleep_data |>
  mutate(
    date = as.Date(date),
    crossover = as.logical(crossover),
    across(
      .cols = all_of(c("his", "hfs", "fps", "tts", "waso")),
      .fns = hms::parse_hms
    ),
    awakenings = as.numeric(awakenings),
    tts_fps = as.numeric(tts) / as.numeric(fps)
  ) |>
  relocate(tts_fps, .after = waso)

primary_sleep_data |> glimpse()

G.11 Filter and Tidy the Secondary Sleep Data

secondary_sleep_data <-
  secondary_sleep_data |>
  mutate(
    date = as.Date(date),
    crossover = as.logical(crossover),
    across(
      .cols = all_of(c("his", "hfs", "fps", "tts", "waso")),
      .fns = hms::parse_hms
    ),
    awakenings = as.numeric(awakenings),
    tts_fps = as.numeric(tts) / as.numeric(fps)
  ) |>
  relocate(tts_fps, .after = waso)

secondary_sleep_data |> glimpse()

G.12 Filter and Tidy the Medical Record Data

medical_record_data <-
  medical_record_data |>
  mutate(
    across(
      .cols = all_of(c("partograph_begin", "baby_birth")),
      .fns = ~ ymd_hms(.x, tz = "America/Sao_Paulo")
    ),
    across(
      .cols = -all_of(c("id", "partograph_begin", "baby_birth")),
      .fns = as.numeric
    ),
    across(
      .cols = starts_with("apgar_"),
      .fns = as.integer
    )
  )

medical_record_data |> glimpse()

medical_record_data <-
  medical_record_data |>
  mutate(
    across(
      .cols = all_of(c("partograph_begin", "baby_birth")),
      .fns = ~ ymd_hms(.x, tz = "America/Sao_Paulo")
    ),
    across(
      .cols = -all_of(c("id", "partograph_begin", "baby_birth")),
      .fns = as.numeric
    ),
    across(
      .cols = matches("^begin_dilation$|^apgar_"),
      .fns = as.integer
    )
  ) |>
  mutate(
    partogram_duration =
      interval(
        start = partograph_begin,
        end = baby_birth
      ) |>
        lubritime::int_duration()
  )

medical_record_data

G.13 Import, Tidy and Transform the Processed Actigraphy Data

G.13.1 Get Paths to Processed Actigraphy Data Files

files <-
  dir_ls(
    path = actigraphy_processed_data_dir,
    type = "file",
    regexp = ".txt$"
  )

G.13.2 Import, Tidy, Transform, and Write Data Files

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  i_id <- str_remove(basename(i), "_actigraphy-processed-data.txt")

  i_ga_start <-
    control_form_data |>
    filter(id == i_id) |>
    pull(ga_start)

  i_data <-
    i |>
    read_acttrust(tz = "America/Sao_Paulo") |>
    mutate(
      ga =
        i_ga_start |>
        ga_point(timestamp) |>
        as.numeric() |>
        as.duration() |>
        add(
          as_hms(timestamp) |>
            as.numeric() |>
            as.duration()
        ),
      ga_week =
        ga |>
        as.numeric() |>
        divide_by(
          ddays(7) |>
            as.numeric()
        ) |>
        as.integer()
    ) |>
    suppressMessages() |>
    suppressWarnings()

  i_data |>
    write_rds(
      path(
        actigraphy_processed_data_dir,
        paste0(i_id, "_actigraphy-processed-data.rds")
      )
    )

  cli_progress_update()
}

cli_progress_done()

G.14 Import and Group the Validated Actigraphy Data

G.14.1 Get Paths to the Data Files

files <-
  dir_ls(
    path = actigraphy_processed_data_dir,
    type = "file",
    regexp = "_actigraphy-processed-data.rds$"
  )

G.14.2 Import and Group the Data

actigraphy_data <-
  files |>
  map_dfr(read_rds)

actigraphy_data |> glimpse()

G.15 Import, Group and Transform SRI Data

See the SRI Analysis Report Quarto notebook report to process and write the Sleep Irregularity Index (SRI) data.

G.15.1 SRI by Gestational Age Week

G.15.1.1 Get Paths to the Data Files

files <-
  dir_ls(
    path = sri_data_dir,
    type = "file",
    regexp = "_sri-data.rds$"
  )

G.15.1.2 Import, Group and Transform the Data

sri_data <- tibble()

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  sri_data <-
    i |>
    read_rds() |>
    mutate(
      time =
        time |>
        as.numeric() |>
        divide_by(60) |>
        as.integer() |>
        dminutes() |>
        as.numeric() |>
        as_hms()
    ) |>
    bind_rows(sri_data)

  cli_progress_update()
}

cli_progress_done()

sri_data <-
  sri_data |>
  mutate(ga_week = as.factor(ga_week)) |>
  arrange(ga_week)

sri_data |> glimpse()

G.15.2 SRI Mean Data by Gestational Age Week

G.15.2.1 Get Paths to Data Files

files <-
  dir_ls(
    path = sri_data_dir,
    type = "file",
    regexp = "_sri-mean-data.rds$"
  )

G.15.2.2 Import and Group Data

sri_mean_data <- tibble()

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  sri_mean_data <-
    i |>
    read_rds(i) |>
    bind_rows(sri_mean_data)

  cli_progress_update()
}

cli_progress_done()

G.15.2.3 Transform SRI Mean Data

sri_mean_data <-
  sri_mean_data |>
  mutate(ga_week = as.factor(ga_week)) |>
  arrange(ga_week)

sri_mean_data |> glimpse()

G.16 Get SRI Data Unique IDs

sri_data_ids <- sri_data |> pull(id) |> unique()

sri_data_ids

G.17 Import, Group and Transform Sleep Proportions Data

See the SRI Analysis Report Quarto notebook report to process and write the Sleep Proportions Index data.

G.17.1 Sleep Proportions by Gestational Age Week

G.17.1.1 Get Paths to Data Files

files <-
  dir_ls(
    path = sleep_prop_dir,
    type = "file",
    regexp = "_sleep-prop-data.rds$"
  )

G.17.1.2 Group and Transform Data

sleep_prop_data <- tibble()

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  sleep_prop_data <-
    i |>
    read_rds() |>
    mutate(
      time =
        time |>
        as.numeric() |>
        divide_by(60) |>
        as.integer() |>
        dminutes() |>
        as.numeric() |>
        as_hms()
    ) |>
    bind_rows(sleep_prop_data)

  cli_progress_update()
}

cli_progress_done()

sleep_prop_data <-
  sleep_prop_data |>
  mutate(ga_week = as.factor(ga_week)) |>
  arrange(ga_week)

sleep_prop_data |> glimpse()

G.17.2 Sleep Proportions Mean Data by Gestational Age Week

G.17.2.1 Get Paths to Mean Data Files

files <-
  dir_ls(
    path = sleep_prop_dir,
    type = "file",
    regexp = "_sleep-prop-mean-data.rds$"
  )

G.17.2.2 Import and Group Mean Data

sleep_prop_mean_data <- tibble()

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  sleep_prop_mean_data <-
    i |>
    read_rds(i) |>
    bind_rows(sleep_prop_mean_data)

  cli_progress_update()
}

cli_progress_done()

G.17.2.3 Transform Mean Data

sleep_prop_mean_data <-
  sleep_prop_mean_data |>
  mutate(ga_week = as.factor(ga_week)) |>
  arrange(ga_week)

sleep_prop_mean_data |> glimpse()

G.18 Import, Group and Transform Illuminance Levels Data

See the SRI Analysis Report Quarto notebook report to process and write the Illuminance Data data.

G.18.1 Illuminance Levels by Gestational Age Week

G.18.1.1 Get Paths to Data Files

files <-
  dir_ls(
    path = lux_data_dir,
    type = "file",
    regexp = "_lux-data.rds$"
  )

G.18.1.2 Group and Transform Data

lux_data <- tibble()

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  lux_data <-
    i |>
    read_rds() |>
    mutate(
      time =
        time |>
        as.numeric() |>
        divide_by(60) |>
        as.integer() |>
        dminutes() |>
        as.numeric() |>
        as_hms()
    ) |>
    bind_rows(lux_data)

  cli_progress_update()
}

cli_progress_done()

lux_data <-
  lux_data |>
  mutate(ga_week = as.factor(ga_week)) |>
  arrange(ga_week)

lux_data |> glimpse()

G.18.2 Illuminance Levels Mean Data by Gestational Age Week

G.18.2.1 Get Paths to Data Files

files <-
  dir_ls(
    path = lux_data_dir,
    type = "file",
    regexp = "_lux-mean-data.rds$"
  )

G.18.2.2 Import and Group Mean Data

lux_mean_data <- tibble()

cli_progress_bar(total = length(files), clear = FALSE)

for (i in files) {
  lux_mean_data <-
    i |>
    read_rds(i) |>
    bind_rows(lux_mean_data)

  cli_progress_update()
}

cli_progress_done()

G.18.2.3 Transform Mean Data

lux_mean_data <-
  lux_mean_data |>
  mutate(ga_week = as.factor(ga_week)) |>
  arrange(ga_week)

lux_mean_data |> glimpse()

G.19 Define Valid Cases

valid_cases <-
  field_form_data |>
  left_join(
    control_form_data,
    by = "id"
  ) |>
  left_join(
    sleep_data,
    by = "id"
  ) |>
  left_join(
    medical_record_data,
    by = "id"
  ) |>
  filter(birth_center_delivery == TRUE) |>
  drop_na(labor_duration, baby_weight, fps_mean, waso_mean) |>
  pull(id)

valid_cases

valid_cases |> glimpse()

G.20 Create Analysis Data

data <-
  control_form_data |>
  left_join(
    sleep_data,
    by = "id"
  ) |>
  left_join(
    medical_record_data,
    by = "id"
  ) |>
  left_join(
    field_form_data,
    by = "id"
  ) |>
  mutate(
    waso_prop = (as.numeric(waso_mean)) / as.numeric(fps_mean)
  ) |>
  filter(id %in% valid_cases)

data |> glimpse()

G.21 Quantitative Variables

G.21.1 Days from Actigraphy Processed Data

ids <- model_data |> pull(id) |> unique()

files <-
  dir_ls(
    path = actigraphy_processed_data_dir,
    type = "file",
    regexp = paste0(ids, "_actigraphy-processed-data.rds$", collapse = "|")
  )

days_by_id <- tibble()

for (i in files) {
  i_id <- str_remove(basename(i), "_actigraphy-processed-data.rds")
  i_data <- i |> read_rds()
  i_start <- i_data |> pull(timestamp) |> first()
  i_end <- i_data |> pull(timestamp) |> last()

  days_by_id <-
    tibble(
      id = i_id,
      days = (i_end - i_start) |> as.numeric()
    ) |>
    bind_rows(days_by_id)
}

days_by_id <- days_by_id |> arrange(id)

days_by_id |> print(n = Inf)

days_by_id |> pull(days) |> summary()

G.22 Partogram Duration Model

form <- formula(partogram_duration ~ begin_dilation + baby_head_perimeter + baby_abdominal_perimeter)

# baby_weight
# baby_length
# baby_head_perimeter
# baby_thoracic_perimeter
# baby_abdominal_perimeter

model <-
  parsnip::linear_reg() |>
  parsnip::set_engine("lm")

fit <-
  model |>
  parsnip::fit(
    formula = form,
    data = medical_record_data
  )

fit |>
  broom::tidy() |>
  print(n = Inf)

fit |>
  broom::glance() |>
  pivot_longer(cols = everything()) |>
  print(n = Inf)

fit_augment <- fit |> broom::augment(new_data = medical_record_data)

fit_error <-
  fit_augment |>
  transmute(
    .pred = .pred,
    .resid = .resid,
    partogram_duration = partogram_duration |> as.numeric(),
    error = .resid / partogram_duration
  )

fit_augment <-
  fit |>
  broom::augment(
    new_data = tibble(
      begin_dilation = 1
    )
  )

medical_record_data |>
  ggplot2::ggplot(
    ggplot2::aes(
      x = partogram_duration,
      y = begin_dilation
    )
  ) +
  ggplot2::geom_point() +
  ggplot2::geom_smooth(
    method = "lm",
    color = "#FC2913"
  ) +
  ggplot2::labs(
    x = "Duração do Partograma (h)",
    y = "Dilatação do Colo do Útero (cm)"
  ) +
  ggplot2::scale_x_time()

G.22.1 Correlation Matrices

data |>
  plotr:::plot_ggally(
    cols = c(
      "labor_duration", "age", "family_income", "bmi_before"
    ),
    labels = c(
      "Duração do parto", "Idade", "Renda familiar", "IMC"
    )
  )

data |>
  plotr:::plot_ggally(
    cols = c(
      "labor_duration", "gestations", "deliveries", "abortions"
    ),
    labels = c(
      "Duração do parto", "# Gestações", "# Partos", "# Abortos"
    )
  )

data |>
  plotr:::plot_ggally(
    cols = c(
      "labor_duration", "baby_weight", "baby_length",
      "baby_head_perimeter",  "baby_thoracic_perimeter",
      "baby_abdominal_perimeter"
    ),
    labels = c(
      "Duração do parto", "Peso do bebê", "Comprimento do bebê",
      "Per. da cabeça do bebê", "Per. torácico do bebê", "Per. abdominal do bebê"
    )
  )

data |>
  plotr:::plot_ggally(
    cols = c("labor_duration", "fps_mean", "waso_mean", "awakenings_mean"),
    labels = c("Duração do parto", "FPS", "WASO", "Desp.")
  )

data |>
  mutate(
    his_mean = his_mean |> lubritime::link_to_timeline(),
    hfs_mean = hfs_mean |> lubritime::link_to_timeline()
  ) |>
  plotr:::plot_ggally(
    cols = c(
      "his_mean", "hfs_mean", "fps_mean", "tts_mean", "tts_fps_mean",
       "waso_mean", "awakenings_mean"
    ),
    labels = c(
      "HIS", "HFS", "FPS", "TTS", "TTS/FPS", "WASO", "Desp."
    )
  )

data |>
  plotr:::plot_ggally(
    cols = c(
      "fps_mean", "tts_mean", "waso_mean", "awakenings_mean",
      "age", "family_income", "bmi_before"
    ),
    labels = c(
      "FPS", "TTS", "WASO", "Desp.", "Idade", "Renda", "IMC"
    )
  )

grDevices::dev.off()

G.22.2 Actograms

actigraphy_data <-
  here("data") |>
  dir_ls(type = "file") |>
  str_subset(paste0("f1e30a89")) |>
  read_acttrust(tz = "America/Sao_Paulo")

actigraphy_data |>
  actogram(
    days = -1,
    latitude = brazil_state_latitude("sp"),
    longitude = brazil_state_longitude("sp")
    # locale = "pt_BR.UTF-8",
    # x_label = "2 Dias — Horas",
    # y_label = "Dias",
    # labels = c(
    #   "1" = "Sono",
    #   "2" = "Despertar",
    #   "4" = "Offwrist",
    #   "base" = "Atividade (PIM)",
    #   "lp" = "Fase clara",
    #   "dp" = "Fase escura"
    # ),
    # colors = c(
    #   "1" = "#410085",
    #   "2" = "#FFB426",
    #   "4" = "#FC2913",
    #   "base" = "#000040",
    #   "lp" = "#FFFFFF",
    #   "dp" = "#DBD7D3"
    # )
  )

G.22.3 Sleep Irregularity Index (SRI)

See the SRI Analysis Report for statistics and visualizations related to Sleep Irregularity Index (SRI).

G.22.4 Sleep Proportions

See the SRI Analysis Report for statistics and visualizations related to Sleep Proportions.

G.22.5 Illuminance Levels

See the SRI Analysis Report for statistics and visualizations related to Illuminance Levels.

G.22.6 Groups

group_data <-
  control_form_data |>
  left_join(
    sleep_data,
    by = "id"
  ) |>
  left_join(
    medical_record_data,
    by = "id"
  ) |>
  filter(id %in% valid_cases) |>
  dplyr::select(labor_duration, fps_mean, waso_mean, baby_weight) |>
  drop_na() |>
  mutate(
    waso_prop = (as.numeric(waso_mean)) / as.numeric(fps_mean),
    group_fps =
      case_when(
        fps_mean < hms::parse_hm("07:00") ~ "<7h (Poor)",
        fps_mean >= hms::parse_hm("07:00") ~ ">7h (Good)"
      ) |>
      factor(levels = c("<7h (Poor)", ">7h (Good)")),
    group_waso =
      case_when(
        waso_prop > 0.1 ~ ">10% (Poor)",
        waso_prop < 0.1 ~ "<10% (Good)"
      ) |>
      factor(levels = c(">10% (Poor)", "<10% (Good)"))
  ) |>
  drop_na()

row_list <- list(
  c("group_waso", ">10% (Poor)"),
  c("group_waso", "<10% (Good)"),
  c("group_fps", "<7h (Poor)"),
  c("group_fps", ">7h (Good)")
)

out <- tibble()

for (i in row_list) {
  out <-
    out |>
    bind_rows(
    tibble(
      n =
        group_data |>
        filter(!!as.symbol(i[1]) == i[2]) |>
        nrow(),
      mean =
        group_data |>
        filter(!!as.symbol(i[1]) == i[2]) |>
        pull(labor_duration) |>
        as.numeric() |>
        mean(na.rm = TRUE),
      sd =
        group_data |>
        filter(!!as.symbol(i[1]) == i[2]) |>
        pull(labor_duration) |>
        as.numeric() |>
        stats::sd(na.rm = TRUE),
      labor_duration = paste0(
        mean |>
          hms::as_hms() |>
          lubritime::round_time() |>
          as.character(),
        " ± ",
        sd |>
          hms::as_hms() |>
          lubritime::round_time() |>
          as.character()
      )
    ) |>
      mutate(group = i[2]) |>
      dplyr::select(group, everything())
  )
}

out |>
  dplyr::select(n, labor_duration) |>
  as.data.frame() %>%
  `row.names<-`(
    c(
      "WASO >10% (Poor)",
      "WASO <10% (Good)",
      "FPS <7h (Poor)",
      "FPS >7h (Good)"
    )
  )

Table G.1: Source: Created by the author

Diferenças na duração do trabalho de parto em função do tempo acordado após o início do sono (Wake After Sleep Onset (WASO)) e da duração do sono (Fase Principal do Sono (FSP)). NA = Não disponível (Not Available). (\(n = 8\)).

G.22.6.1 FSP

group_data |>
  plot_ggally(
    cols = c("labor_duration", "baby_weight", "fps_mean", "group_fps"),
    mapping = ggplot2::aes(color = group_fps)
  )

group_data |>
  mutate(
    labor_duration =
      labor_duration |>
      as.numeric() |>
      hms::as_hms(),
    group_fps = case_match(
      group_fps,
      "<7h (Poor)" ~ "<7h (Ruim)",
      ">7h (Good)" ~ ">7h (Bom)"
    )
  ) |>
  ggplot2::ggplot(
    ggplot2::aes(
      x = labor_duration,
      y = baby_weight,
      color = group_fps
    )
  ) +
  ggplot2::geom_point() +
  ggplot2::stat_smooth(method = lm) +
  ggplot2::labs(
    x = "Duração do parto",
    y = "Peso do bebê (g)",
    color = "FPS"
  ) +
  ggplot2::scale_x_continuous(labels = labels_hms) +
  viridis::scale_color_viridis(end = 0.5, discrete = TRUE)

Figure G.1: Fonte: Elaborado pela autora.

Relação da variável labor_duration (Duração do parto) com a covariante baby_weight (Peso do bebê) (gramas) por grupos com duração de sono baixa/ruim (<7h) (cor roxa) e com duração de sono boa/adequada (>7h) (cor verde).
A duração do sono foi medida por meio da média da Fase Principal de Sono (FSP) nos últimos sete dias que antecederam o parto. As linhas representam regressões lineares. A área cinza representa o intervalo de confiança da estimativa.

G.22.6.2 WASO

group_data |>
  plot_ggally(
    cols = c("labor_duration", "baby_weight", "waso_prop", "group_waso"),
    mapping = ggplot2::aes(color = group_waso)
  )

group_data |>
  mutate(
    labor_duration =
      labor_duration |>
      as.numeric() |>
      hms::as_hms(),
    group_waso = case_match(
      group_waso,
      ">10% (Poor)" ~ ">10% (Ruim)",
      "<10% (Good)" ~ "<10% (Bom)"
    )
  ) |>
  ggplot2::ggplot(
    ggplot2::aes(
      x = labor_duration,
      y = baby_weight,
      color = group_waso
    )
  ) +
  ggplot2::geom_point() +
  ggplot2::stat_smooth(method = lm) +
  ggplot2::labs(
    x = "Duração do parto",
    y = "Peso do bebê (g)",
    color = "WASO"
  ) +
  ggplot2::scale_x_continuous(labels = labels_hms) +
  viridis::scale_color_viridis(end = 0.5, discrete = TRUE, direction = -1) +
  ggplot2::guides(color = ggplot2::guide_legend(reverse = TRUE))

Figure G.2: Fonte: Elaborado pela autora.

Relação da variável labor_duration (Duração do parto) com a covariante baby_weight (Peso do bebê) (em gramas) por grupos com baixa qualidade de sono (>10% de WASO/FPS) (cor roxa) e com qualidade de sono boa/adequada (<10% de WASO/FPS) (cor verde).
Os dados se referem aos últimos sete dias que antecederam o parto. FPS = Fase Principal do Sono, representando a duração do período entre o horário de início do sono e o horário final do sono. WASO = Wake After Sleep Onset ou tempo acordado após o início do sono, usado como proxy para medir a qualidade do sono. As linhas representam regressões lineares. A área em cinza representa o intervalo de confiança da estimativa.

G.23 Hypothesis Testing (Sleep and Labor Duration)

Variável dependente: labor_duration

\[ \begin{cases} \text{H}_{0}: \overline{X}_{\text{Ruim}} <= \overline{X}_{\text{Bom}} \\ \text{H}_{a}: \overline{X}_{\text{Ruim}} > \overline{X}_{\text{Bom}} \end{cases} \]

See Ellis (2010) to learn more about Hedges’ g and other effect size measures.

\[ \text{g de Hedges} = \cfrac{\overline{\text{X}}_{\text{Ruim}} - \overline{\text{X}}_{\text{Bom}}}{\text{SD}_{\text{pooled}}} \]

model_data <-
  control_form_data |>
  left_join(
    sleep_data,
    by = "id"
  ) |>
  left_join(
    medical_record_data,
    by = "id"
  ) |>
  filter(id %in% valid_cases) |>
  dplyr::select(labor_duration, fps_mean, waso_mean, baby_weight) |>
  drop_na()

model_data

G.23.1 FPS (Fase Principal do Sono) (Main Sleep Phase)

TipConsensus Statement

• Adults should sleep 7 or more hours per night on a regular basis to promote optimal health.
• Sleeping less than 7 hours per night on a regular basis is associated with adverse health outcomes, including weight gain and obesity, diabetes, hypertension, heart disease and stroke, depression, and increased risk of death. Sleeping less than 7 hours per night is also associated with impaired immune function, increased pain, impaired performance, increased errors, and greater risk of accidents.
• Sleeping more than 9 hours per night on a regular basis may be appropriate for young adults, individuals recovering from sleep debt, and individuals with illnesses. For others, it is uncertain whether sleeping more than 9 hours per night is associated with health risk.
• People concerned they are sleeping too little or too much should consult their healthcare provider.

(Watson et al., 2015)

model_data_fps <-
  model_data |>
  dplyr::select(labor_duration, fps_mean, baby_weight) |>
  mutate(
    group =
      case_when(
        fps_mean < hms::parse_hm("07:00") ~ "Poor",
        fps_mean >= hms::parse_hm("07:00") ~ "Good"
      ) |>
      factor(levels = c("Poor", "Good"))
  ) |>
  mutate(
    labor_duration = as.numeric(labor_duration)
  ) |>
  relocate(group) |>
  dplyr::select(group, labor_duration, baby_weight) |>
  drop_na()

model_data_fps

model_data_fps |> report::report()

G.23.1.1 Welch’s T-Test

G.23.1.1.1 Power Analysis (a posteriori)

See the a priori power analysis at https://danielvartan.github.io/sample-size.

Desired power: \((1 - \beta) \geq 0.8\)

# A posteriori

pwr_analysis_fps_anconva <- pwrss::pwrss.t.2means(
  mu1 = model_data_fps |> get_labor_duration("Poor") |> mean(na.rm = TRUE),
  mu2 = model_data_fps |> get_labor_duration("Good") |> mean(na.rm = TRUE),
  sd1 = model_data_fps |> get_labor_duration("Poor") |> stats::sd(na.rm = TRUE),
  sd2 = model_data_fps |> get_labor_duration("Good") |> stats::sd(na.rm = TRUE),
  paired = FALSE,
  n2 = model_data_fps |> get_n("Good"),
  kappa = get_n(model_data_fps, "Poor") / get_n(model_data_fps, "Good"),
  power = NULL,
  alpha = 0.05,
  welch.df = TRUE,
  alternative = "superior"
)

G.23.1.1.2 Test

model_t_test_fps <- stats::t.test(
  x = model_data_fps |> get_labor_duration("Poor"),
  y = model_data_fps |> get_labor_duration("Good"),
  alternative = "greater"
)

model_t_test_fps

report::report(model_t_test_fps)

G.23.1.1.3 Effect Size

tibble(
  fps_mean_poor =
    model_data_fps |>
    get_labor_duration("Poor") |>
    mean(na.rm = TRUE) |>
    hms::as_hms(),
  fps_mean_good =
    model_data_fps |>
    get_labor_duration("Good") |>
    mean(na.rm = TRUE) |>
    hms::as_hms(),
  diff = (fps_mean_poor - fps_mean_good) |> hms::as_hms()
)

effectsize::cohens_d(
  x = model_data_fps |> get_labor_duration("Poor"),
  y = model_data_fps |> get_labor_duration("Good"),
  pooled_sd = TRUE,
  mu = 0,
  paired = FALSE,
  adjust = TRUE, # For small samples
  ci = 0.95,
  alternative = "greater",
  verbose = TRUE
) |>
  effectsize::interpret_hedges_g()

G.23.1.2 ANCOVA

See Rutherford (2011, pp. 237–261) to learn more about ANOVA-ANCOVA assumptions.

G.23.1.2.1 Test

model_aov_fps <- stats::aov(
  formula = labor_duration ~ baby_weight * group,
  data = model_data_fps
)

model_aov_fps

model_aov_fps |> summary()

model_aov_fps |> report::report()

G.23.1.2.2 Effect Size

Veja Cohen (1992) para saber mais.

Tamanho de efeito para duração do parto e FSP encontrado por Lee & Gay (2004):

k <- 3
num_df <- 2
den_df <- 125
n_total_tst <- den_df + k
f_test <- 5.54

f_unb_tst <- f_unb(k, n_total_tst, f_test)

f_unb_tst

Tamanho de efeito a posteriori:

effect_size_aov_fps <-
  model_aov_fps |>
  effectsize::cohens_f(
    partial = FALSE,
    generalized = FALSE,
    method = "eta",
    model2 = NULL,
    ci = 0.95,
    alternative = "greater",
    verbose = TRUE,
  ) |>
  effectsize::interpret_eta_squared("cohen1992")

effect_size_aov_fps

G.23.1.2.3 Power Analysis (a posteriori)

Veja a análise de poder a priori em: https://danielvartan.github.io/sample-size.

Poder desejado: \((1 - \beta) \geq 0.8\)

# A posteriori

pwr_analysis_fps_anconva <- pwrss::pwrss.f.ancova(
  f2 = effect_size_aov_fps$Cohens_f[3]^2,
  n.way = 1,
  n.levels = 2,
  n.covariates = 1,
  # Not really true, since each group has a different number of observations.
  n = model_data_fps |>
    get_labor_duration("Good") |>
    length(),
  alpha = 0.05,
  verbose = TRUE
)

pwr_analysis_fps_anconva

G.23.2 WASO (Wake After Sleep Onset)

As an estimate of sleep disruption, WASO is reported as the percentage of minutes awake divided by minutes in bed after falling asleep. During a typical 7- to 8-hour sleep period, WASO of 15% represents more than an hour of wake time after falling asleep. WASO between 5% and 10% is typical in healthy, nonpregnant women5 and greater than 15% was considered severe sleep disruption for this study.(Lee & Gay, 2004[p. 2042]).

model_data_waso <-
  model_data |>
  dplyr::select(labor_duration, waso_mean, fps_mean, baby_weight) |>
  mutate(
    waso_prop = (as.numeric(waso_mean)) / as.numeric(fps_mean)
  ) |>
  mutate(
    group =
      case_when(
        waso_prop > 0.1 ~ "Poor",
        waso_prop < 0.1 ~ "Good"
      ) |>
      factor(levels = c("Poor", "Good"))
  ) |>
  mutate(
    labor_duration = as.numeric(labor_duration)
  ) |>
  relocate(group) |>
  dplyr::select(group, labor_duration, baby_weight) |>
  drop_na()

model_data_waso

model_data_waso |> report::report()

G.23.2.1 Welch’s T-Test

G.23.2.1.1 Power Analysis (a posteriori)

See the a priori power analysis at https://danielvartan.github.io/sample-size.

Desired power: \((1 - \beta) \geq 0.8\)

pwr_analysis_waso_t_test <- pwrss::pwrss.t.2means(
  mu1 = model_data_waso |>
    get_labor_duration("Poor") |>
    mean(na.rm = TRUE),
  mu2 = model_data_waso |>
    get_labor_duration("Good") |>
    mean(na.rm = TRUE),
  sd1 = model_data_waso |>
    get_labor_duration("Poor") |>
    stats::sd(na.rm = TRUE),
  sd2 = model_data_waso |>
    get_labor_duration("Good") |>
    stats::sd(na.rm = TRUE),
  paired = FALSE,
  n2 = model_data_waso |> get_n("Good"),
  kappa = get_n(model_data_waso, "Poor") / get_n(model_data_waso, "Good"),
  power = NULL,
  alpha = 0.05,
  welch.df = TRUE,
  alternative = "superior"
)

G.23.2.1.2 Test

model_t_test_waso <- stats::t.test(
  x = model_data_waso |> get_labor_duration("Poor"),
  y = model_data_waso |> get_labor_duration("Good"),
  alternative = "greater"
)

model_t_test_waso

report::report(model_t_test_waso)

G.23.2.1.3 Effect Size

tibble(
  waso_mean_poor =
    model_data_waso |>
    get_labor_duration("Poor") |>
    mean(na.rm = TRUE) |>
    hms::as_hms(),
  waso_mean_good =
    model_data_waso |>
    get_labor_duration("Good") |>
    mean(na.rm = TRUE) |>
    hms::as_hms(),
  diff = (waso_mean_poor - waso_mean_good) |> hms::as_hms()
)

effectsize::cohens_d(
  x = model_data_waso |> get_labor_duration("Poor"),
  y = model_data_waso |> get_labor_duration("Good"),
  pooled_sd = TRUE,
  mu = 0,
  paired = FALSE,
  adjust = TRUE, # For small samples
  ci = 0.95,
  alternative = "greater",
  verbose = TRUE
) |>
  effectsize::interpret_hedges_g()

G.23.2.2 ANCOVA

G.23.2.2.1 Test

model_aov_waso <- stats::aov(
  formula = labor_duration ~ baby_weight * group,
  data = model_data_waso
)

model_aov_waso

model_aov_waso |> summary()

model_aov_waso |> report::report()

G.23.2.2.2 Effect Size

See Cohen (1992) to learn more.

Effect size for labor duration and WASO found by Lee & Gay (2004):

k <- 3
num_df <- 2
den_df <- 125
n_total_waso <- den_df + k
f_test <- 6.93

f_unb_waso <- f_unb(k, n_total_waso, f_test)

f_unb_waso

Tamanho de efeito a posteriori:

effect_size_aov_waso <-
  model_aov_waso |>
  effectsize::cohens_f(
    partial = FALSE,
    method = "eta",
    model2 = NULL,
    ci = 0.95,
    alternative = "greater",
    verbose = TRUE,
  ) |>
  effectsize::interpret_eta_squared("cohen1992")

effect_size_aov_waso

G.23.2.2.3 Power Analysis (a posteriori)

Veja a análise de poder a priori em: https://danielvartan.github.io/sample-size.

Poder desejado: \((1 - \beta) \geq 0.8\)

pwr_analysis_waso_ancova <- pwrss::pwrss.f.ancova(
  f2 = effect_size_aov_waso$Cohens_f[3]^2,
  n.way = 1,
  n.levels = 2,
  n.covariates = 1,
  # Not really true, since each group has a different number of observations.
  n = model_data_waso |>
    get_labor_duration("Poor") |>
    length(),
  alpha = 0.05,
  verbose = FALSE
)

pwr_analysis_waso_ancova

G.24 ANCOVA Assumptions

Assumption 1

Predictor is known (ANOVA). Each condition contains a random sample of the population of such scores (Rutherford, 2011[p. 237]).

Assumption 2

Normality (ANOVA). The scores in each condition are distributed normally (Rutherford, 2011[p. 237]).

Assumption 3

Independence (ANOVA). The scores in each condition are independent of each other (Rutherford, 2011[p. 237]).

Assumption 4

Common variance (homoscedasticity) (ANOVA). The variances of the scores in each experimental condition are homogeneous (Rutherford, 2011[p. 237]).

Assumption 5

Independence of the covariates (homoscedasticity) (ANCOVA). The covariate is independent of the treatments (Rutherford, 2011[p. 237]).

Assumption 6

Linear mean (ANCOVA). In each treatment group the relationship between the covariate and the dependent variable is linear (the covariate and dependent variable are expressed at the first power only) (Rutherford, 2011[p. 240]).

Assumption 7

Homogeneity of the coefficients (ANCOVA). The regression coefficients of the dependent variable on the covariate in each treatment group are homogeneous (Rutherford, 2011[p. 240]).

Cohen, J. (1992). A power primer. Psychological Bulletin, 112(1), 155–159. https://doi.org/10.1037/0033-2909.112.1.155

Ellis, P. D. (Ed.). (2010). The essential guide to effect sizes: Statistical power, meta-analysis, and the interpretation of research results (a). Cambridge University Press.

Lee, K. A., & Gay, C. L. (2004). Sleep in late pregnancy predicts length of labor and type of delivery. American Journal of Obstetrics and Gynecology, 191(6), 2041–2046. https://doi.org/10.1016/j.ajog.2004.05.086

Rutherford, A. (2011). ANOVA and ANCOVA: A GLM approach (2nd ed.). Wiley. https://doi.org/10.1002/9781118491683

Watson, N. F., Badr, M. S., Belenky, G., Bliwise, D. L., Buxton, O. M., Buysse, D., Dinges, D. F., Gangwisch, J., Grandner, M. A., Kushida, C., Malhotra, R. K., Martin, J. L., Patel, S. R., Quan, S. F., & Tasali, E. (2015). Recommended amount of sleep for a healthy adult: A joint consensus statement of the American Academy of Sleep Medicine and Sleep Research Society. Journal of Clinical Sleep Medicine, 11(6), 591–592. https://doi.org/10.5664/jcsm.4758