Quantified Self Experiments / Continous Glucose Monitoring, food testing and diet decision-making

High glucose levels is one of the main reasons of metabolic pathologies, especially in diabetes. Glucose spikes come from food consumption, especially from high-processed simple-carbohydrate rich food. Health specialists recommend to exclude food which high in GI (glyceminc index) from our diet.

These recommendations doesnt always work because every person have different glycemic responses to food and how do you know if them work especially for you? I did a personal experiment to find out how different food affect my glucose levels and expose methodology to compare different foods to make dietary decisions.

Continous glucose measuring

My current methodology:

 

How to test food?

Now we want to test the food, how to do that? Here is a paper which explains how GI calculation is being performed. There is iAUC concept: we take glucose baseline before we eat certain food and then measure Area Under the Curve (AUC).

Generally speaking AUC is a time-weighed average glucose for 2 hours after we eat. First of all we generate a list of food we want to test, then we find out how much carbohydrates in that food per 100g. I often use https://nutritiondata.self.com/ or just google the food. After we know how much carbo per 100g, we should calculate portion with 50g of carbo.

For example, today i've tested carrots. Carrots have 9.6g carbo per 100g, which means i should eat 50*100/9.6 = 520g of carrots. Here is mine :) It took 35 minutes for me to eat half kilo. Rice test was a lot easier (20g carbo per 100g, 250g portion).

After the test i dont eat anything for 2 hours and spend this time in a rest.I make my tests in the morning with single isolated product. Sometimes i make 2 tests per days, with 2h+ pause between then and make sure glucose returned to resting values before next test.

For carrots i see a glucose raise in the Glimp and next step is to calculate iAUC. I export Glimp data into CSV and run my R script (provided below) which calculates iAUC and visualize data.

Here is what i've gotten for carrots.

Dots represents average glucose trends. Bars represents max iAUC in mg/dl for each period

Here we can see that my iAUC for carrots is 32.5 mg/dl or 1.8 mmol/L. iAUC formula based on relative glucose values which overcomes sensors systemic bias (too low or too high absolute values, but accurate trends), so generally speaking iAUC represents relative change in glucose levels.

Here is exact formula:

Is 1.8 mmol/L a high one? To answer that question we can compare foods from personal tests, here is my list.

mmol/L Product
5.46 grape
4.67 beetroot
4.08 boil pumpkin
3.28 potato
3.06 rice
2.83 wholegrain pasta
2.81 banana
2.78 pomegranate
2.58 mango
2.53 millet porridge
2.52 white bread
2.35 buckwheat porridge
2.3 barley porridge
2.19 sweet potato
2.16 watermelon
2.01 tomato
1.62 peach
1.57 quinoa
1.37 apple
1.02 rice noodle
0.89 kefir*
0.87 wholegrain bread
0.67 baked pistachio
0.64 oat porridge
0.62 oat porridge whole
0.54 dark chocolate 100%*
0.51 pistachio
0.19 persimmon
0.11 carob powder
0.1 pear

Here is my iAUC's for different foods in mmol/L per 50g of carbo. Products with "*" have values per 15g of carbo (due to hardness to eat 50g carbo portion). Carrots place is somewhere between peach and tomato. It took about 2 months of wearing CGM (spent 5 sensors) to build this list without rush (i did one period every 2-3 months).

To produce GI we need to measure how much iAUC is for pure glucose and calculate rate = product iAUC / glucose iAUC.

If we just want to compare different foods and predict glucose levels spike, we can just use iAUC as is. So we need our CGM data and a simple list of dates with food we have tested, for example in google spreadsheets.

If i cannot eat portion of 50g carbo then i go down to 15g portion and test it. Also sometimes there is 45g of carbo in my portions, that means i will multiply result by 50/45 = 1.11. But if i've tested only 10g of product, i do not multiply it by 50/10=5 because there are no guarantee for strict linear relationship between glucose levels and amount of carbo.

Diet decisions

Having list of foods and their glycemic responses i can make dietary decisions:

By knowing how much different products influence glucose levels we can decide how much to eat. I know that 50g of mango carbo's (300g mango portion) increase my glucose (iAUC) by 2.58 mmol/L and by knowing that i can decide to eat half portion (150-200g) and be confident that my glucose will stay in range. We can see that 50g carbo from tomatoes is slightly above 2 mmol/L, but there is 5g of carbo per 100g of tomato. So eating even half kilo of tomatoes would not spike my glucose and it's safe to eat them.

Also we should take into account that

Anything else?

There is interesting paper "Glucotypes reveal new patterns of glucose dysregulation" about glucotypes where scientists clustered people by their glucemic response variations into 3 types: low, moderate and severe.

This approach seems to be pretty interesting, huge part of people from severe group and significant part of moderate group developed type 2 diabetes (T2D) in the long run, but most part of low variability group didnt develop T2D. That points me that it may be worth to aim for low variability group.

Scientists made R Code for clustering available on github and provide a webtool for easy glucotyping. Data needs to be formatted into .tsv (tab separated) file with date and glucose values (every 5 minutes) and looking like this:

My R script generates that file in script folder and names it "classify.tsv". Alternative ways is playing with google sheets.

The webtool gives me low glucotype. This period includes some tests which reveal spiking food and i'm not going to eat it anymore (for example rice test at 9 Aug). And even with these spikes i've gotten low glucotype. So it seems my glucose metabolism is fine and my T2D risk is pretty low. On my first period with CGM i got moderate glucotype, then after some foods tests and diet adjusting my results are always stay in a low group.

Data availability & Information

Welcome for questions, suggestions and critics in comments below.

Original unmodified (exported) Glimp raw data for oura is here.

R Code:

library(mise); mise()
RNGkind(sample.kind = "Rejection")
options(scipen=999)
describe.data <- function(df, digits = 2) {
  for(i in 1:length(colnames(df))) {
    x <- df[,colnames(df)[i]]
    if(is.numeric(x)) {
      print(paste0(colnames(df)[i], "(", length(x) , "): ",round(mean(x),digits),"±",round(sd(x),digits)," ", round(median(x),digits),"[",round(quantile(x,.25),digits),",",round(quantile(x,.75), digits),"] ", round(min(x),digits),"-",round(max(x),digits))) 
    }
  }
}
cache.filename <- 'GlicemiaMisurazioni.csv'; cache.folder <- ''; cache.file <- paste0(cache.folder, cache.filename)
if(file.exists(cache.file)) { load(cache.file) } else {
  url <- 'https://blog.kto.to/uploads/'
  gimp <- read.csv(paste0(url, cache.filename), sep = ';', header = F, skipNul = T)
  save(gimp, file = cache.file)
}
data.csv <- gimp
data.csv <- data.csv[,c("V3","V6","V7")]
colnames(data.csv) <- c("datetime","raw","calibrated")
data.csv$datetime <- as.POSIXct(data.csv$datetime, origin = "1970-01-01")
date.start <- as.POSIXct("2022-08-06", timezone = "UTC")
date.end <- as.POSIXct("2022-08-20", timezone = "UTC")
data.csv <-  data.csv[data.csv$datetime >= date.start & data.csv$datetime <= date.end,]
# check for NA data
library(mice)
md.pattern(data.csv, rotate.names = T)
summary(data.csv)
df <- data.csv[,c("datetime","calibrated")]
describe.data(df)
# visualize missing data
library(naniar)
vis_miss(df, sort_miss = TRUE)
# delete.outliers <- F
# cols <- c("calibrated")
# for(i in 1:length(cols)) {
#   df.out <- analyze.outliers(df[, cols[i]], cols[i])
#   if(delete.outliers) { 
#     df <- df[df.out, ];
#   }
# }
library(ggplot2); library(scales)
ggplot(df, aes(x = datetime, y = calibrated)) + 
  geom_point() + scale_x_datetime(breaks = "1 hour") + 
  theme(axis.text.x = element_text(angle=90))
period.after <- 120 # in minutes
period.pre <- 30 # in minutes
df.f <- df
df.f$iAUC <- NA
n <- nrow(df.f)
for(i in 1:n) {
  point <- df.f[i,]
  ds <- point$datetime; de <- point$datetime + period.after*60; dpre <- point$datetime - period.pre*60
  # get points for next 2 hours (postprandial), excluding start one
  points.after <- df.f[df.f$datetime >= ds & df.f$datetime <= de,]
  # get points for previous 30 mins to calc a baseline
  points.pre <- df.f[df.f$datetime >= dpre & df.f$datetime <= ds,]
  # get baseline as mean from previous 30 minutes
  baseline <- mean(points.pre$calibrated)
  
  points.after$iAUC <- NULL
  date.prev <- ds
  g.prev <- baseline
  n.a <- nrow(points.after)
  iAUC <- 0
  if(n.a > 3) {
    for(j in 1:n.a) {
      point.j <- points.after[j,]
      time.diff <- as.numeric(point.j$datetime - date.prev)
      point.iAUC <- 0
      if(points.after[j,]$calibrated > baseline) {
        # calc iAUC https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2769714/
        # Si = dt*((Gi - G0) + (Gi-1 - G0))/2; iAUC = sum(S)/120
        point.iAUC <- time.diff*((point.j$calibrated-baseline)+(g.prev-baseline))/2
        iAUC <- iAUC + point.iAUC
        #print(paste0(point.j$datetime, " (",time.diff,") > ", points.after[j,]$iAUC, " bl: ", baseline, " cur.p: ", point.j$calibrated, " g.p: ", g.prev))
      }
      g.prev <- point.j$calibrated
      date.prev <- point.j$datetime
    }
    df.f[i,]$iAUC <- round(iAUC / 120,1)
  }
}
df.f <- na.omit(df.f)
library(lubridate)
df.p <- aggregate(cbind(iAUC = df.f$iAUC, glucose = df.f$calibrated), by = list(datetime = floor_date(df.f$datetime, "15 minute")), FUN = max)
q <- quantile(df.p$iAUC,.9)
df.p$label <- ""; df.p[df.p$iAUC > q, "label"] <- df.p[df.p$iAUC > q, "iAUC"]
ggplot(df.p, aes(x = datetime, y = iAUC, label=label)) + 
  geom_bar(stat = "identity") + 
  geom_point(aes(x=datetime, y=glucose/2)) + 
  geom_text(hjust=0, vjust=0) + 
  theme(axis.text.x = element_text(angle=45))
df.tsv <- data.frame("GlucoseDisplayTime" = df$datetime, "GlucoseValue" = df$calibrated)
df.tsv <- aggregate(cbind(GlucoseValue = df.tsv$GlucoseValue), by = list(GlucoseDisplayTime = floor_date(df.tsv$GlucoseDisplayTime, "5 minute")), FUN = mean)
df.tsv$GlucoseValue <- round(df.tsv$GlucoseValue)
write.table(df.tsv, file="classify.tsv", quote = F,sep = "\t", row.names = F)

Statistical analysis

RStudio version 2022.07.1 and R version 4.2.1.