Hall et al. (2018) Data Preparation

Companion script for cgm.qmd (Parts 5-6) — CGM Tutorial, JSM Short Course

Author

Irina Gaynanova

Published

July 31, 2026

Note

This page is published for reference only. The code below is shown but not executed, since it needs the raw Hall et al. (2018) files, which aren’t hosted in this repository (see Downloading the raw data below). hall_data.Rda - the output of running this script - is included in the repo, and is all cgm.qmd needs.

This script needs a few packages not required by cgm.qmd itself (DBI, RSQLite, for the clinical SQLite database):

install.packages("dplyr")
install.packages("tidyr")
install.packages("DBI")
install.packages("RSQLite")
install.packages("biosensors.usc")

Purpose

cgm.qmd (Parts 5-6) uses the full Hall et al. (2018) dataset (S1, S5 and S6 supplementary files - see below for download links, and see cgm.qmd for citation and per-file descriptions). Reading and cleaning these raw files (a ~105K-row raw CGM file, a SQLite clinical database, and building the non-parametric density/quantile representations for all 57 subjects) is the slow part of that tutorial, and doesn’t need to be redone every time the tutorial is edited or re-rendered.

This script does that wrangling once and saves every derived object cgm.qmd needs to hall_data.Rda. cgm.qmd then just load()s that file. hall_data.Rda is already included in this repository, so you do not need to run this script (or download the raw files below) to follow the tutorial. Re-run this script only if you want to reproduce that step yourself, or when the wrangling logic changes.

Downloading the raw data

The raw files are not included in this repository (they’re sizeable and are already publicly hosted). Hall et al. (2018), PLOS Biology, doi:10.1371/journal.pbio.2005143, CC BY 4.0 licensed. Download the three supplementary files below into a local HallFullData/ folder (not tracked by git - see .gitignore):

Supplementary file	Contents	Download	Save as
S1 Data	Raw CGM recordings, all 57 subjects	https://doi.org/10.1371/journal.pbio.2005143.s010	`HallFullData/pbio.2005143.s010`
S5 Data	Clinical variables and glucotype	https://doi.org/10.1371/journal.pbio.2005143.s014	`HallFullData/pbio.2005143.s014.db`
S6 Data	Glucose response to standardized meals	https://doi.org/10.1371/journal.pbio.2005143.s015	`HallFullData/pbio.2005143.s015.tsv`

All paths below are relative to this project’s root (same convention as cgm.qmd), so raw files are expected under HallFullData/.

Part 5 inputs: meal-aligned glucose curves

Reads “S6 Data” (pbio.2005143.s015.tsv) and builds meal-aligned curves for all subjects/meals, then narrows to one Cereal Flakes (“CF”) curve per subject (first occurrence) for the FPCA/function-on-scalar illustrations in cgm.qmd, and pulls each of those subjects’ A1C.

hall_meals_full = read.delim("HallFullData/pbio.2005143.s015.tsv", stringsAsFactors = FALSE)
hall_meals_full$time = as.POSIXct(hall_meals_full$time, format = "%Y-%m-%d %H:%M:%S")

# "Low" is Dexcom's convention for glucose below the sensor's detection floor (<=40 mg/dL)
hall_meals_full$GlucoseValue[hall_meals_full$GlucoseValue == "Low"] = "40"
hall_meals_full$GlucoseValue = as.numeric(hall_meals_full$GlucoseValue)

# "PB 1", "CF 2", etc. -> meal type + replicate number
hall_meals_full$meal_type = trimws(sub("[0-9]+$", "", hall_meals_full$Meal))
hall_meals_full$replicate = as.integer(sub("\\D+", "", hall_meals_full$Meal))

hall_meals_full = hall_meals_full %>%
  dplyr::rename(id = userID, gl = GlucoseValue) %>%
  # a small number of (id, Meal, time) combinations have duplicate CGM readings; average them
  dplyr::group_by(id, Meal, meal_type, replicate, time) %>%
  dplyr::summarize(gl = mean(gl), .groups = "drop")

head(hall_meals_full)

Meal onset isn’t given directly, but each window starts 30 minutes before the meal, so onset = first sampled time + 30 minutes. We use that to compute time relative to meal onset, and drop the small number of meals with a substantially incomplete window (missing CGM data during that period).

meal_curves = hall_meals_full %>%
  dplyr::group_by(id, Meal) %>%
  dplyr::mutate(onset = min(time) + 30 * 60,
                rel_time_hr = as.numeric(difftime(time, onset, units = "hours")),
                window_min = as.numeric(difftime(max(time), min(time), units = "mins"))) %>%
  dplyr::ungroup() %>%
  dplyr::filter(window_min >= 170) %>%
  dplyr::select(-window_min)

dplyr::n_distinct(meal_curves$id)
dplyr::n_distinct(paste(meal_curves$id, meal_curves$Meal))

Restrict to Cereal Flakes (“CF”), earliest occurrence per subject, so every subject contributes exactly one curve, then pivot to a subject-by-timepoint matrix for FPCA:

cf_first = meal_curves %>%
  dplyr::filter(meal_type == "CF") %>%
  dplyr::group_by(id) %>%
  dplyr::filter(onset == min(onset)) %>%
  dplyr::ungroup() %>%
  dplyr::mutate(rel_min = round(rel_time_hr * 60))

dplyr::n_distinct(cf_first$id)

fpca_wide = cf_first %>%
  dplyr::select(id, rel_min, gl) %>%
  tidyr::pivot_wider(names_from = rel_min, values_from = gl,
                      names_prefix = "t_")

subject_ids = fpca_wide$id
fpca_mat = fpca_wide %>%
  dplyr::select(-id) %>%
  as.matrix()

# Drop rows with all NA (no observed curve), keeping subject_ids in sync
keep_rows = apply(fpca_mat, 1, function(x) any(!is.na(x)))
fpca_mat = fpca_mat[keep_rows, , drop = FALSE]
subject_ids = subject_ids[keep_rows]
dim(fpca_mat)

Pull A1C (from “S5 Data”, the clinical SQLite database) for these subjects, for the function-on-scalar regression in cgm.qmd:

con = DBI::dbConnect(RSQLite::SQLite(), "HallFullData/pbio.2005143.s014.db")
clinical = DBI::dbGetQuery(con, "SELECT userID, A1C FROM clinical")
DBI::dbDisconnect(con)

subject_covs = tibble::tibble(id = subject_ids) %>%
  dplyr::left_join(clinical %>% dplyr::rename(id = userID), by = "id")

sum(!is.na(subject_covs$A1C))

Part 6 inputs: full-cohort distributional data

Reads “S1 Data” (pbio.2005143.s010, raw CGM for all 57 subjects) and “S5 Data” (clinical covariates for all 57 subjects), ignoring meal timing entirely.

hall_full = read.delim("HallFullData/pbio.2005143.s010", stringsAsFactors = FALSE)
hall_full$GlucoseValue[hall_full$GlucoseValue == "Low"] = "40"
hall_full$gl = as.numeric(hall_full$GlucoseValue)
hall_full = hall_full %>%
  dplyr::rename(id = subjectId, time = DisplayTime) %>%
  dplyr::select(id, time, gl)

dplyr::n_distinct(hall_full$id)

con = DBI::dbConnect(RSQLite::SQLite(), "HallFullData/pbio.2005143.s014.db")
clinical_all = DBI::dbGetQuery(con, "SELECT userID AS id, A1C, diagnosis, glucotype FROM clinical")
DBI::dbDisconnect(con)

biosensors.usc::load_data() expects its inputs as CSV files - raw long-format data (id, time, value) and a covariates file (id plus any covariates) - and then converts them into a non-parametric density estimate ($density) and an empirical quantile estimate ($quantiles) per subject. This per-subject density/quantile estimation over all 57 subjects is the single slowest step in the whole pipeline, which is exactly why it belongs here rather than in the tutorial itself.

raw_path = "HallFullData/hall_raw_biosensors.csv"
covs_path = "HallFullData/hall_covs_biosensors.csv"

hall_full %>%
  dplyr::rename(value = gl) %>%
  write.csv(raw_path, row.names = FALSE)

write.csv(clinical_all, covs_path, row.names = FALSE)

hall_biosensor = load_data(raw_path, covs_path)
names(hall_biosensor)

Save everything

save(meal_curves, fpca_mat, subject_ids, subject_covs,
     hall_full, clinical_all, hall_biosensor,
     file = "hall_data.Rda")

file.info("hall_data.Rda")$size

cgm.qmd loads this single file at the start of Part 5 (and again, defensively, at the start of Part 6) via load("hall_data.Rda"), and focuses on fitting/plotting (fpca.face, pffr, clustering, frechetreg_univar2wass, optimal_thresholds) rather than re-deriving these objects.