Aerobic Threshold Estimation using DFA a1
The demarcation of training zones has long been a topic amongst athletes, coaches, scientists, etc. While there are a ton of different methods for estimating the heart rate and power associated with the anaerobic threshold (e.g. critical power, functional threshold power, etc.) there is rather less about the aerobic threshold (e.g., VT1, AET, LT1). This is unfortunate because much of the literature on endurance performance suggests that it is important to spend most of your time training below this threshold. I don't believe that the exact demarcation of this zone is very important, but although people that do a lot of endurance training might think they know what their aerobic threshold power/heart rate are, they might be likely to over-estimate it. In any case if you are interested in long distance performance it is probably the most important performance parameter to track.
The developer of HRV4Training: Marco Altini implemented the estimation of the HRV derived DFA a1 described in these papers (1, 2) in the app HRVLogger.
I did a preliminary test with a 6 hour steady state ride where I rode at around 200 watts (with two food stops) (here). Marco's app helpfully allows you to export a csv file of all the features computed (including alpha1), which you can then merge with other data as desired. I set the app to record 2 minute smoothed values, aligned the ride file (a Garming .fit file), and computed an equivalent rolling 2 minute estimate of power. The aforementioned papers suggest that AET is power/heart rate where alpha1 = .75. I fit a generalized additive model to power and heart rate (seperately) as a function of this alpha1 parameter, which you can see below. According to this hard threshold I spent around 73% of the ride below AET. This is relatively close to "intensity" (NP / FTP) in this case (which was .7). You can find the code at the bottom of the page.
This seems interesting so far. Being an ignoramous as far as physiology is concerned I cannot comment on the validity of the feature or the threshold, but it seems better than using zones determined by max heart rate, which is what Seiler suggests. In my case this gives a heart rate of around 140, and my power at that heart rate is around 175 or so in contrast to the DFA a1 estimate (from this ride anyway). Although it is certainly the case that if I do rides at this heart rate/power I am more recovered and more ready for hard intervals I can easily sustain above those numbers for upwards of 8 hours. Does that mean they are wrong? Seiler doesn't seem particularly keen on the max HR based estimation of the aerobic threshold. The DFA a1 estimate puts the threshold around 215w and 154bpm.
library(data.table)
library(ggplot2)
library(lubridate)
library(cycleRtools)
library(mgcv)
dat = fread("~/Dropbox/Apps/Heart Rate Variability Logger/2021-3-13_Features_.csv")
dat = dat[, list(date, alpha1, heart_rate)]
garmin = read_fit("~/hrv/6425268604_ACTIVITY.fit")
setDT(garmin)
garmin[, date := as_datetime(timestamp) - hours(8)]
garmin = garmin[date >= min(dat$date) & date <= max(dat$date)]
garmin[, c("start", "end") := .(date - 2 * 60, date)]
garmin[, power := .SD[.SD, on = .(date >= start, date <= end),
by = .EACHI, lapply(.SD, mean, na.rm = TRUE), .SDcols = c("power.W")][, (1L:2L) := NULL]]
test = merge(dat, garmin[, list(date, power)], by = "date")
test[, below_aet := ifelse(alpha1 >= .75, TRUE, FALSE)]
fit_power = gam(power ~ alpha1, data = test)
fit_heart_rate = gam(heart_rate ~ alpha1, data = test)
predict_power = as.data.table(do.call("cbind", predict(fit_power, newdata = test, se.fit = TRUE)))
predict_heart_rate = as.data.table(do.call("cbind", predict(fit_heart_rate, newdata = test, se.fit = TRUE)))
preds = rbind(predict_heart_rate, predict_power)
preds[, lower := fit + qnorm(.025) * se.fit]
preds[, upper := fit + qnorm(.975) * se.fit]
aet = data.table('power' = round(predict(fit_power, newdata = data.table(alpha1 = .75))),
'heart_rate' = round(predict(fit_heart_rate, newdata = data.table(alpha1 = .75))))
test = melt(test, id.vars = c("date", "alpha1", "below_aet"))
test = cbind(test, preds)
ggplot(test, aes(alpha1, value)) + geom_point() +
geom_vline(aes(xintercept = .75), linetype = "dashed") +
geom_line(aes(alpha1, fit)) +
geom_ribbon(aes(ymin = lower, ymax = upper), alpha = .25) +
facet_wrap(~ variable, scales = "free_y") +
labs(title = paste0('aet power:hr ', aet$power, ':', aet$heart_rate))
ggsave("~/hrv/bd_to_coast_aet.png", width = 10, height = 5)