ORIGINAL RESEARCH article
A New Detection Method Defining the Aerobic Threshold for Endurance Exercise and Training Prescription Based on Fractal Correlation Properties of Heart Rate Variability
- 1College of Medicine, University of Central Florida, Orlando, FL, United States
- 2Lattice Training Ltd., Chesterfield, United Kingdom
- 3School of Health Sciences, College of Education, Health and Human Development, University of Canterbury, Christchurch, New Zealand
- 4Center for Sports and Physical Education, Julius Maximilians University of Wuerzburg, Wuerzburg, Germany
- 5Department of Performance, Neuroscience, Therapy and Health, Faculty of Health Sciences, MSH Medical School Hamburg, University of Applied Sciences and Medical University, Hamburg, Germany
The short-term scaling exponent alpha1 of detrended fluctuation analysis (DFA a1), a nonlinear index of heart rate variability (HRV) based on fractal correlation properties, has been shown to steadily change with increasing exercise intensity. To date, no study has specifically examined using the behavior of this index as a method for defining a low intensity exercise zone. The aim of this report is to compare both oxygen intake (VO2) and heart rate (HR) reached at the first ventilatory threshold (VT1), a well-established delimiter of low intensity exercise, to those derived from a predefined DFA a1 transitional value. Gas exchange and HRV data were obtained from 15 participants during an incremental treadmill run. Comparison of both VO2 and HR reached at VT1 defined by gas exchange (VT1 GAS) was made to those parameters derived from analysis of DFA a1 reaching a value of 0.75 (HRVT). Based on Bland Altman analysis, linear regression, intraclass correlation (ICC) and t testing, there was strong agreement between VT1 GAS and HRVT as measured by both HR and VO2. Mean VT1 GAS was reached at 39.8 ml/kg/min with a HR of 152 bpm compared to mean HRVT which was reached at 40.1 ml/kg/min with a HR of 154 bpm. Strong linear relationships were seen between test modalities, with Pearson’s r values of 0.99 (p < 0.001) and.97 (p < 0.001) for VO2 and HR comparisons, respectively. Intraclass correlation between VT1 GAS and HRVT was 0.99 for VO2 and 0.96 for HR. In addition, comparison of VT1 GAS and HRVT showed no differences by t testing, also supporting the method validity. In conclusion, it appears that reaching a DFA a1 value of 0.75 on an incremental treadmill test is closely associated with crossing the first ventilatory threshold. As training intensity below the first ventilatory threshold is felt to have great importance for endurance sport, utilization of DFA a1 activity may provide guidance for a valid low training zone.
Training zone identification is part of the foundation for exercise intensity distribution study and implementation (Stöggl and Sperlich, 2019. Traditionally, the upper limit of the low intensity range (zone 1 in a 3 zone model) for intensity distribution for endurance exercise and training prescription has been represented by the first ventilatory (VT1) or lactate threshold (LT1) (Seiler and Kjerland, 2006; Esteve-Lanao et al., 2007; Mann et al., 2013; Pallarés et al., 2016). Although there may be different schools of thought on what type of distribution is “optimal” (polarized vs. pyramidal or threshold) both models are defined by having the major portion of training in zone 1. In addition, several training approaches for endurance athletes recommend spending large amounts of exercise time in a low intensity zone (Muñoz et al., 2014; Stöggl and Sperlich, 2014, 2019; Bourgois et al., 2019; Casado et al., 2019). Gold standard methods to obtain VT1 or LT1 revolve around either formal gas exchange testing or invasive blood lactate sampling. These procedures can be costly, require special test equipment, trained operators, ongoing calibration and verification. Even if these methods are utilized, there is disagreement on their accuracy as both visual (Yeh et al., 1983; Meyer et al., 1996) and automated (Ekkekakis et al., 2008) gas exchange analysis can be subject to substantial error. In addition, LT1 assessment can vary depending on the chosen concept of determination (Newell et al., 2007; Faude et al., 2009; Jamnick et al., 2018). Training guided by erroneous intensity targets could lead to potential adverse consequences such as prolonged cardiac parasympathetic recovery (Seiler and Kjerland, 2006; Stanley et al., 2013), central and muscular fatigue (Noakes et al., 2005; Venhorst et al., 2018), glycogen depletion (Beneke et al., 2011), and gastrointestinal barrier disruption (van Wijck et al., 2012). In view of the difficulties involved in gas exchange analysis, lactate test availability, invasiveness, and accuracy, a search for alternate methods of identifying the limits of low intensity exercise seem worthwhile.
Cardiac interbeat interval variation, commonly referred to as heart rate variability (HRV), has been extensively studied in both resting states (Shaffer and Ginsberg, 2017) as well as during dynamic exercise (Hottenrott and Hoos, 2017; Michael et al., 2017). Certain HRV indexes have been observed to change as exercise intensity rises, potentially providing information regarding an individual’s physiologic status (Tulppo et al., 1996; Casties et al., 2006; Sandercock and Brodie, 2006; Karapetian et al., 2008; Michael et al., 2017; Gronwald et al., 2018, 2019a,b, 2019c). It has also been shown that several of the examined HRV indexes also change during lower intensities (Tulppo et al., 1996; Sandercock and Brodie, 2006; Karapetian et al., 2008; Botek et al., 2010; Michael et al., 2017) making them potentially suitable for zone 1 delineation. However, despite some initial interest, widespread usage for the specific purpose of low intensity training limitation has not occurred. Frequency-domain parameters such as high frequency (HF) power have been noted to be unreliable in a sizable fraction of individuals with up to 20% of subjects not having identifiable breakpoints (Cottin et al., 2007). Time domain measures such as SDNN were found to closely relate with VT1 but little follow-up or verification has been done (Karapetian et al., 2008). The SD1 is another index that has been examined during exercise. It is based on a Poincare plot of each RR interval graphed against the preceding interval and is related to short term trends in RR patterns often assigned to nonlinear indexes (Shaffer and Ginsberg, 2017), although it is mathematically equivalent to another time domain index (Ciccone et al., 2017). While showing potential as a low intensity marker in some earlier studies (Tulppo et al., 1996) other evidence indicates that SD1 was already suppressed in young athletes at the first tested work rate of 60% VO2MAX making it less useful for zone 1 delineation (Blasco-Lafarga et al., 2017).
One nonlinear index, the short-term scaling exponent alpha1 based on Detrended Fluctuation Analysis (DFA a1), has generated interest as both an indicator of autonomic nervous system regulation as well as an overall marker of organismic demands Gronwald and Hoos, 2020). Originally, Peng et al. (1995) developed this method to measure scale-invariant behavior; this involved the evaluation of trends of all sizes in the presence or absence of fractal correlation properties in a heart rate (HR) time series (Yeh et al., 2010). Thus, the DFA method allows for the quantification of the degree of correlation and fractal scale of a HRV signal resulting in dimensionless measures. The short-term scaling exponent DFA a1 is based on the fractal dynamics (self-similarity) of the cardiac beat-to-beat pattern and provides insights into correlation properties of HR time series caused by physiological processes (Peng et al., 1995). DFA a1 values indicate time series correlation properties with approximately 1.5 indicating a strongly correlated pattern and ≤0.5 for anti-correlated pattern with random behavior; approximately 1.0 signifies a mix of uncorrelated and maximally correlated signal components (represents a balance between complete unpredictability (randomness) and predictability (strong correlations), also associated with fractal (self-similar) behavior (Platisa and Gal, 2008). Larger values of DFA a1 represent a smoother time series and smaller values of DFA a1 represent coarser ones (Peng et al., 1995; Goldberger et al., 2002). Within this framework, DFA a1 has been shown to decline as work rate rises, starting from strongly correlated patterns (value of 1.5) at rates well below the first ventilatory threshold (VT1), transitioning (values of 1.0–0.5) through values representing uncorrelated, less complex white noise behavior at moderate to high work rates, then finally showing anti-correlated behavior at the highest intensities (values of <0.5) (Gronwald et al., 2019c; Gronwald and Hoos, 2020). Given this relationship, there may be an opportunity to assist athletes in delineating intensity training zones by observing the change in DFA a1 with increasing exercise intensity (Gronwald et al., 2020; Rogers, 2020).
The purpose of this report is to validate a predefined DFA a1 value of 0.75 with the exercise intensity at VT1 obtained during an incremental treadmill run to exhaustion. This is to be done by a direct comparison of the VT1 intensity based on both absolute VO2 and HR obtained during gas exchange with the same measures derived from analysis of DFA a1 behavior. If it can be shown that a predefined “boundary” value in the DFA a1 index occurs near the VT1, this could establish a basis for further research exploiting a non-linear autonomic nervous system related marker in prospective exercise and training intensity distribution.
Materials and Methods
Seventeen male volunteers aged 19–52, without previous medical history, current medications or physical issues were tested. A background questionnaire regarding medical history was reviewed along with information of the potential testing risks then institutionally approved consent was given. Approval for the study was granted by the University of Derby, United Kingdom (LSREC_1415_02) and conformed to the principles of the Declaration of Helsinki. Participants did not consume caffeine, alcohol or any stimulant for the 24 h before testing. Background data for each subject included, age, body weight, and training volume in hours per week (Table 1). All testing was done in the afternoon and at least 3 h post meal. No exercise was performed the day prior to the test. Two participants with a high degree of cardiac ectopy (frequent atrial premature beats and atrial trigeminy) during testing were excluded from analysis.
Participants performed an incremental VO2MAX test on a motorized treadmill (Woodway, Birmingham, United Kingdom). The treadmill was set for the Bruce protocol with increases in speed and inclination from 2.7 km/h at ten percent grade, increasing by 1.3 km/h and two percent grade every 3 min until volitional exhaustion. A fan was used for cooling.
Gas Exchange Testing and Calculation of the First Ventilatory Threshold
Gas exchange kinetics were recorded continuously using a breath-to-breath metabolic cart (Metalyzer 3B; Cortex Biophysik GmbH Germany). In addition, a Polar H7 (Polar Electro Oy, Kempele, Finland) was wirelessly paired to the Metalyzer cart for the purpose of HR recording concurrent with gas exchange data. VO2, VCO2, PetO2, PetCO2, Ve/VO2, Ve/VCO2, and HR were imported into Microsoft Excel 365 for analysis. The native gas exchange analysis feature of the Metalyzer was not used due to the unreliability of many automated VT1 calculations (Ekkekakis et al., 2008). Graphing of the above parameters were done to derive VT1, VO2MAX, and VO2 vs. time. No averaging was done for either gas exchange parameters or HR. Inspection of the VO2 over time relationship was done to determine any significant plateau of the VO2 curve for estimation of VO2MAX and VO2 linearity. If a significant plateau was found, compensation for calculating both VO2MAX and the VO2 over time equation was done. To reduce the chance of failure to identify the VT1 by gas exchange (VT1 GAS) based on a single method, evaluation was done according to the triple detection method consisting of V slope, Ve/VO2, and excess CO2 from Gaskill et al. (2001) as well as the PetO2 nadir from Binder et al. (2008). Based on the quality and consistency of the plots, the excess CO2 method was chosen to be used for all participants and reviewed independently by two investigators (Figure 1A). VO2 was plotted over the elapsed time of the incremental test to produce a linear regression equation. VO2 at the time of VT1 was based on linear regression from the VO2 over time relationship.
Figure 1. (A) Usage of “Excess CO2” technique to determine VT1 GAS. The intersection of the baseline and first rise in excess CO2 corresponds to the time when VT1 occurs; (B) DFA a1 plotted over time, the area of linear drop of DFA a1 from about 1.0 to 0.5 is used to determine the VO2 (from the VO2 vs. time relation) at HRVT; (C) DFA a1 plotted against HR, the area of linear drop of DFA a1 from 1.0 to 0.5 is used to determine the HR at HRVT. All data taken from Subject #2.
RR Measurements and Calculation of DFA a1 Derived Threshold
A 3-lead ECG (MP36; Biopac Systems Ltd.) with a sampling rate of 1,000 Hz was used to record the subject’s ECG/RR times series. Biopac filter settings were set to 0.05 Hz high-pass filter and 150 Hz low-pass filter. Electrodes were placed in the CM5 distribution after appropriate skin cleansing and shaving if necessary. Sample data from the MP36 was saved as .acq files. ECG files for each subject were imported into Kubios 3.3.2 (Biosignal Analysis and Medical Imaging Group, Department of Physics, University of Kuopio, Kuopio, Finland). Kubios preprocessing settings were at the default values including the RR detrending method which was kept at “Smoothn priors” (Lambda = 500, Tarvainen et al., 2014). For DFA a1 estimation, the root mean square fluctuation of the integrated and detrended data is measured in observation windows of different sizes. The data are then plotted against the size of the window on a log-log scale. The scaling exponent represents the slope of the line, which relates (log) fluctuation to (log) window size (Mendonca et al., 2010). DFA a1 window width was set to 4 ≤ n ≤ 16 beats.
For the detection of a HRV derived threshold, a DFA a1 value of 0.75 was chosen based on this being the midpoint between a fractal behavior of the HR time series of 1.0 (seen with very light exercise) and an uncorrelated value of 0.5 which represents white noise, random behavior (seen with high intensity exercise). A value of 0.75 has also been used as a cut-off value for survival curves and mortality rate assessment during resting conditions (Huikuri et al., 2000).
The following procedure was used to indicate at what level of running intensity (as VO2 or HR) the DFA a1 would cross a value of 0.75: DFA a1 was calculated from the incremental exercise test RR series using 2 min time windows with a recalculation every 5 s throughout the test. Two minute time windowing was chosen based on the reasoning of Chen et al. (2002). The rolling time window measurement was used to better delineate rapid changes in the DFA a1 index over the course of the test. Each DFA a1 value is based on the RR series 1 min pre and 1 min post the designated time stamp. For example, at a time of 10 min into the testing, the DFA a1 is calculated from the 2 min window starting from minute 9 and ending at minute 11 and labeled as the DFA a1 at 10 min. Based on a rolling time recalculation every 5 s, the next data point would occur at 10:05 min (start 9:05 min and end 11:05 min).
Plotting of DFA a1 vs. time was then performed. Inspection of the DFA a1 relationship with time generally showed a reverse sigmoidal curve with a stable area above 1.0 at low work rates, a rapid, near linear drop reaching below 0.5 at higher intensity, then flattening without major change. A linear regression was done on the subset of data consisting of the rapid near linear decline from values near 1.0 (correlated) to approximately 0.5 (uncorrelated). The time of DFA a1 reaching 0.75 was calculated based on the linear regression equation from that straight section (Figure 1B). The time of DFA a1 reaching 0.75 was then converted to VO2 using the VO2 vs. time relation, resulting in the VO2 at which DFA a1 equaled 0.75 (HRVT). A similar analysis was done for the HR reached at a DFA a1 of 0.75. First, ECG data from each 2 min rolling window was used to plot the average HR and DFA a1. The HR at which DFA a1 equaled 0.75 was found using the same technique as above, a linear regression through the rapid change section of DFA a1 values of 1.0 to below 0.5, with a subsequent equation for HR and DFA a1 (Figure 1C). Using a fixed variable of DFA a1 equals 0.75, the resulting HR was obtained. The HR at DFA a1 0.75 (based on ECG data) was then compared to the HR at VT1 GAS obtained from the metabolic cart data (based on the Polar H7).
Visual inspection of the entire test recording was done to determine sample quality, noise, arrhythmia, and missing beat artifact. As mentioned above, two participants with a high degree of atrial ectopy were excluded from analysis. The RR series of the included participants was then corrected by the Kubios “automatic method” and exported as text files for further analysis. Percent artifact reported refers those occurring during the linear regression segment (DFA a1 1.0 to near 0.5).
Statistical analysis was performed for the main variables, VO2 at VT1 derived from gas exchange testing, VO2 at DFA a1 0.75, HR at VT1 obtained from gas exchange testing and average HR at DFA a1 0.75. Standard statistical methods were used for the calculation of means and standard deviations (SD). Normal distribution of data was checked by Shapiro–Wilk’s test. The agreement against the Gold Standard VT1 GAS was assessed using intraclass correlation coefficient (ICC), linear regression, Pearson’s r correlation coefficient, standard error of estimate (SEE), coefficient of determination (R2) and Bland Altman plots with limits of agreement (Bland and Altman, 1999). The size of Pearson’s r correlations evaluated as follows; 0.3 ≤ r < 0.5 low; 0.6 ≤ r < 0.8 moderate and r ≥ 0.8 high (Chan, 2003). The paired t-test was used for comparison of VT1 GAS vs. HRVT for both VO2 and HR parameters. For all tests, the statistical significance was accepted as p ≤ 0.05. Cohen’s d was used to denote effect sizes (small effect = 0.2, medium effect = 0.5, large effect = 0.8; Cohen, 1988). Analysis was performed using Microsoft Excel 365 with Real Statistics Resource Pack software (Release 6.8).
Gas Exchange Testing
Individual gas exchange results are presented in Table 2. Both the VO2MAX as well as the percentage of VO2MAX and HR at VT1 GAS varied considerably among participants. VO2MAX ranged between 41 and 74 ml/kg/min. VT1 GAS was reached between 61 and 86% of the VO2MAX and at HRs between 108 and 183 bpm.
RR Interval Quality
The percentage of artifacts was calculated based on the Kubios automatic correction method for each subject’s test data. Since only a portion of the entire treadmill test was used for the linear interpolation of DFA a1, the artifact percentage listed refers to that section only. Artifact percentage for the linear plotted data series was between 0 and 3%, all consisting of atrial premature complexes (Table 2). There were no missed beats due to noise interference or loss of electrode contact. The two participants originally excluded from analysis had significant ectopy, leading to an uninterpretable DFA a1 pattern.
Comparison of VT1 GAS vs. HRVT
The average VT1 GAS was 39.8 ml/kg/min (±8.9) compared to 40.1 ml/kg/min (±8.6) obtained by HRVT. The average HR at VT1 GAS was 152 bpm (±21) compared to 154 bpm (±20) obtained by HRVT. Strong linear relationships were seen between test modalities, with Pearson’s r values of 0.99 (p < 0.001) and.97 (p < 0.001) for VO2 and HR comparisons respectively (Figure 2). Intraclass correlation between VT1 GAS and HRVT was 0.99 for VO2 and 0.96 for HR. The comparison of VT1 GAS and HRVT showed no differences (VO2: p = 0.347, d = 0.030; HR: p = 0.191, d = 0.091). Bland Altman analysis for VT1 GAS vs. HRVT for VO2 (Figure 3) showed a mean difference of −0.33 ml/kg/min (±1.3) with upper and lower limits of 2.2 and −2.9 ml/kg/min. Bland Altman analysis for VT1 GAS vs. HRVT for HR (Figure 3) showed a mean difference of −1.9 bpm (±5) with upper and lower limits of 8 and −12 bpm.
Figure 2. Regression plots for all subject data. (A) Values of VT1 GAS vs. HRVT for VO2; (B) Values of VT1 GAS vs. HRVT for HR. Bisection lines in light gray. SEE, standard error of estimate; R2, coefficient of determination.
Figure 3. Bland Altman Plot of VT1 GAS vs. HRVT for all participants. (A) Values of VT1 GAS vs. HRVT for VO2; (B) Values of VT1 Gas vs. HRVT for HR. Center line in each plot represents the mean difference between each paired value, the top and bottom lines are 1.96 standard deviations from the mean difference.
This study explored whether values of the nonlinear HRV index, DFA a1, pass through a defined transitional zone at workloads near VT1 during an incremental treadmill test. Since many prior reports have shown DFA a1 to decline during incremental exercise (Hautala et al., 2003; Casties et al., 2006; Platisa et al., 2008; Karavirta et al., 2009; Blasco-Lafarga et al., 2017; Gronwald et al., 2019c), our result showing a similar occurrence is not unanticipated. However, none have attempted to directly examine the possibility that the DFA a1 index has a distinct value at the VT1 work rate. Since many of the prior studies looking at DFA a1 response to incremental exercise intensity have used cycling as the exercise modality, it is also reassuring to see analogous results with treadmill running, adding validity to the behavior of this index during other endurance exercise types. The inclusion of a wide range of subject ages, body weights and fitness abilities, lends strength to the application of our results to the general population and its application in different fields of physical exercise and training.
In a recent perspective review (Gronwald et al., 2020), identification of a low intensity exercise zone based on DFA a1 for the purposes of endurance exercise and training prescription was discussed. The mechanism underlying DFA a1 decline with exercise is felt to be related to autonomic balance and a complex interaction of the two main branches, namely parasympathetic withdrawal, sympathetic intensification as well as other factors (Gronwald and Hoos, 2020). Since VT1 is usually seen at a point of significant parasympathetic withdrawal (Tulppo et al., 1996; Sales et al., 2019), leveraging HRV parameters that reliably reflect this occurrence can be of use during endurance exercise and training. Our methodology to determine HRVT utilized the rapid decline of DFA a1 from 1.0 to below 0.5, seen during progressive exercise intensity. The results presented here appear to indicate that VT1 is reached at a midpoint between a fractal behavior of DFA a1 and a pattern of uncorrelated white noise with random behavior, corresponding to a DFA a1 of approximately 0.75. Bland Altman analysis with limits of agreement showed minimal difference between VT1 GAS and HRVT looking at either VO2 or HR measurements. Correlation coefficients and ICC were high for both VO2 and HR based comparisons.
Although the DFA a1 value of 0.75 was chosen theoretically, a brief review of prior investigation is supportive of this figure. In a study of young men performing a cycling ramp test, an average DFA a1 of 0.49 was associated with a lactate measurement of 2.49, indicating that LT1 had already been exceeded (Gronwald et al., 2019c). Other cycling ramp studies in men of different fitness levels seemed to indicate that DFA a1 crossed the value of 0.75 at about 73–78% of VO2MAX (Hautala et al., 2003; Hottenrott and Hoos, 2017), within the approximate realm of VT1 for many individuals (Gaskill et al., 2001; Pallarés et al., 2016). An examination of the DFA a1 response to incremental cycling exercise in teenage males (Blasco-Lafarga et al., 2017) showed an approximate crossing of the 0.75 value at an average intensity near 65% of maximum, also near published ranges of VT1 (50–65% of VO2MAX) in that age group (Runacres et al., 2019). In the current study, there appeared to be little bias in the VO2 or HR associated with HRVT. Although there was some variability in VT1 GAS vs. HRVT parameters, for the most part, associations were similar to that of other comparisons of threshold approaches such as blood lactate or ventilatory parameters (Pallarés et al., 2016), or assessment of gas exchange techniques for VT1 determination (Gaskill et al., 2001). Several participants had relatively high heart rate at the VT1 (Azevedo et al., 2011) of which we do not have an explanation. A strength of this study is that of RR interval quality. Direct ECG visualization was done and a research grade device with a high sample rate was used. No missed or lost beats were seen, and the only artifact type present was APC aberrancy. In view of reports indicating substantial bias in nonlinear HRV indexes with artifact presence (Giles and Draper, 2018; Rincon Soler et al., 2018), a weakening of DFA a1 derived VT1 accuracy could occur with higher artifact occurrence.
A significant advantage of DFA a1 over other HRV indexes for the determination of a low intensity threshold revolves around the nature of testing. Other HRV metrics proposed to identify VT1 such as SDNN (Karapetian et al., 2008), HF power (Cottin et al., 2007) or SD1 (Tulppo et al., 1996) require testing into high intensity zones since they rely on curve interpretation that displays a demonstrable nadir. With DFA a1, once the VT1 boundary area is reached, little additional increase in exercise intensity should be required. The potential benefit of utilizing a fixed DFA a1 value as the VT1 delineation marker is especially attractive in populations unable or ill advised to enter high intensity regions. In addition, for athletes evaluating low intensity training limits, avoidance of exercise ramps to volitional failure may help avert undue stress in a polarized training model.
Limitations and Future Direction
Given the issues with both availability and accuracy of gas exchange or blood lactate testing in determining VT1 for training zone purposes, an alternate modality that employs relatively simple wearable technology seems attractive. However, while DFA a1 monitoring may be a promising approach, several questions need to be addressed. Although this study was done with a wide range of subject age and fitness characteristics, no female participants were tested. If the DFA a1 index behavior is to be considered as a zone 1 delimiter for the general population, further investigation using female subjects is mandatory. Another area of concern is the transfer of the DFA a1 0.75 breakpoint obtained during incremental testing to that of one found during constant load exercise, including moderate length intervals (5 min). No data is available comparing DFA a1 behavior during an incremental ramp to constant load exercise (Gronwald and Hoos, 2020), making automatic transfer of zone boundaries unclear. Whether the index will remain stable for even longer exercise intervals (>60 min) performed below VT1 intensity is another open question as well as day to day repeatability. Another interesting subject to explore is the impact of athlete overtraining on DFA a1 behavior and VT1 prediction accuracy during exercise. Baumert et al. (2006) did show changes in DFA related scaling behavior after intense training, which may provide both a potential source of HRVT bias and an opportunity to screen for overtraining states. Although it seems that ramp protocol slope has minimal effect on the VT1 gas assessment (Weston et al., 2002; Boone and Bourgois, 2012), the analogous assumption needs to be shown in terms of HRVT thresholds. Another area for investigation is whether DFA a1 cut off values are equivalent between chest belt and research grade ECG recordings. Although in this study, the RR intervals were recorded with a research grade ECG device, it may be possible to reproduce similar results with chest belt ECG recordings. In that regard, two major questions need to be addressed. One is that of exercise associated missed beat artifact with possible faulty interpolation strategies by interpreting software, creating potential bias in the calculated DFA a1 values. As mentioned above, several reports have questioned the degree of bias of nonlinear HRV indexes if artifacts are present in the RR series (Giles and Draper, 2018; Rincon Soler et al., 2018; Stapelberg et al., 2018). Artifacts may be of different types such as missed beats or aberrancy. In the current report, no missed beats were seen, and only relatively rare atrial premature complexes were noted. However, two participants exhibited frequent APC aberrancy, had uninterpretable DFA a1 curves and were excluded from group analysis. Given the relative low numbers of participants, no definitive conclusion can be reached regarding artifact bias, but further investigation into effects of missed beats and aberrancy on the use of DFA a1 to delineate zone 1 transition is needed. Second is the question of DFA a1 value precision obtained by diverse monitoring devices possessing different sample rates and prepossessing strategies. Device sample rates have been shown to variably alter DFA a1 values at rest (Voss et al., 1996; Tapanainen et al., 1999; Singh et al., 2015) but may have more significant effects during exercise. Although no recent caffeine use was noted by history, we have no information on prior long term intake patterns which could affect autonomic balance on abrupt discontinuation (La Monica et al., 2018). Finally, it may be possible to answer many of these questions by “repurposing” prior work already done. For instance, a study by Boullosa et al. (2014) assessed the changes in DFA a1 before and after a typical incremental treadmill ramp to exhaustion. A look back at previously acquired RR recordings during the active ramp portion using the methods discussed here may be a way to rapidly acquire needed information about DFA a1 behavior during dynamic exercise.
DFA a1, an index of fractal dynamics and correlations properties of the heart rate time series, was noted to decline during an incremental treadmill run test to exhaustion. The area of most rapid change of this index occurred near the first ventilatory threshold. The point of DFA a1 reaching a value of 0.75 during the incremental treadmill test was directly associated with the first ventilatory threshold as measured by heart rate and gas exchange VO2. As training intensity below the first ventilatory threshold is felt to have great importance for exercise and training prescription in endurance sport, utilization of DFA a1 activity may provide guidance for a valid low training zone boundary without the need for gas exchange or blood lactate testing. Further study of DFA a1 behavior in female participants, during constant load intervals, index stability over long periods of time and across diverse recording devices is recommended. If investigation into these matters remain consistent with the results presented here, obtaining a low intensity zone boundary by automated analysis of a training session via an appropriate wearable device may be possible.
Data Availability Statement
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.
The studies involving human participants were reviewed and approved by University of Derby, United Kingdom (LSREC_1415_02). The patients/participants provided their written informed consent to participate in this study.
BR and TG conceived the study. DG and ND performed the physiologic testing. BR wrote the first draft of the article. BR and TG performed the data analysis. All authors revised it critically for important intellectual content, final approval of the version to be published, and accountability for all aspects of the work.
Conflict of Interest
DG was employed by company Lattice Training.
The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Azevedo, L. F., Perlingeiro, P. S., Brum, P. C., Braga, A. M., Negrão, C. E., and de Matos, L. D. (2011). Exercise intensity optimization for men with high cardiorespiratory fitness. J. Sports Sci. 29, 555–561. doi: 10.1080/02640414.2010.549613
Binder, R. K., Wonisch, M., Corra, U., Cohen-Solal, A., Vanhees, L., Saner, H., et al. (2008). Methodological approach to the first and second lactate threshold in incremental cardiopulmonary exercise testing. Eur. J. Cardiovasc. Prev. Rehab. 15, 726–734. doi: 10.1097/hjr.0b013e328304fed4
Blasco-Lafarga, C., Camarena, B., and Mateo-March, M. (2017). Cardiovascular and Autonomic Responses to a Maximal Exercise Test in Elite Youngsters. Int. J. Sports Med. 38, 666–674. doi: 10.1055/s-0043-110680
Boullosa, D. A., Barros, E. S., del Rosso, S., Nakamura, F. Y., and Leicht, A. S. (2014). Reliability of heart rate measures during walking before and after running maximal efforts. Int. J. Sports Med. 35, 999–1005. doi: 10.1055/s-0034-1372637
Bourgois, J. G., Bourgois, G., and Boone, J. (2019). Perspectives and Determinants for Training-Intensity Distribution in Elite Endurance Athletes. Int. J. Sports Physiol. Perform. 14, 1151–1156. doi: 10.1123/ijspp.2018-0722
Casado, A., Hanley, B., Santos-Concejero, J., and Ruiz-Pérez, L. M. (2019). World-Class Long-Distance Running Performances Are Best Predicted by Volume of Easy Runs and Deliberate Practice of Short-Interval and Tempo Runs. J. Strength. Cond. Res. doi: 10.1519/JSC.0000000000003176
Casties, J. F., Mottet, D., and Le Gallais, D. (2006). Non-linear analyses of heart rate variability during heavy exercise and recovery in cyclists. Int. J. Sports Med. 27, 780–785. doi: 10.1055/s-2005-872968
Ciccone, A. B., Siedlik, J. A., Wecht, J. M., Deckert, J. A., Nguyen, N. D., and Weir, J. P. (2017). Reminder: RMSSD and SD1 are identical heart rate variability metrics. Muscle Nerve. 56, 674–678. doi: 10.1002/mus.25573
Cottin, F., Médigue, C., Lopes, P., Leprêtre, P. M., Heubert, R., and Billat, V. (2007). Ventilatory thresholds assessment from heart rate variability during an incremental exhaustive running test. Int. J. Sports Med. 28, 287–294. doi: 10.1055/s-2006-924355
Ekkekakis, P., Lind, E., Hall, E. E., and Petruzzello, S. J. (2008). Do regression-based computer algorithms for determining the ventilatory threshold agree? J. Sports Sci. 26, 967–976. doi: 10.1080/02640410801910269
Esteve-Lanao, J., Foster, C., Seiler, S., and Lucia, A. (2007). Impact of training intensity distribution on performance in endurance athletes. J. Strength Cond. Res. 21, 943–949. doi: 10.1519/00124278-200708000-00048
Gaskill, S. E., Ruby, B. C., Walker, A. J., Sanchez, O. A., Serfass, R. C., and Leon, A. S. (2001). Validity and reliability of combining three methods to determine ventilatory threshold. Med. Sci. Sports. Exerc. 33, 1841–1848. doi: 10.1097/00005768-200111000-00007
Goldberger, A. L., Amaral, L. A., Hausdorff, J. M., Ivanov, P. C., Peng, C. K., and Stanley, H. E. (2002). Fractal dynamics in physiology: Alterations with disease and aging. PNAS 99, 2466–2472. doi: 10.1073/pnas.012579499
Gronwald, T., Hoos, O., and Hottenrott, K. (2019a). Effects of Acute Normobaric Hypoxia on Non-linear Dynamics of Cardiac Autonomic Activity During Constant Workload Cycling Exercise. Front. Physiol. 10:999. doi: 10.3389/fphys.2019.00999
Gronwald, T., Hoos, O., and Hottenrott, K. (2019b). Effects of a Short-Term Cycling Interval Session and Active Recovery on Non-Linear Dynamics of Cardiac Autonomic Activity in Endurance Trained Cyclists. J. Clin. Med. 6:8.
Gronwald, T., Hoos, O., Ludyga, S., and Hottenrott, K. (2019c). Non-linear dynamics of heart rate variability during incremental cycling exercise. Res. Sports Med. 27, 88–98. doi: 10.1080/15438627.2018.1502182
Gronwald, T., Ludyga, S., Hoos, O., and Hottenrott, K. (2018). Non-linear dynamics of cardiac autonomic activity during cycling exercise with varied cadence. Hum. Mov. Sci. 60, 225–233. doi: 10.1016/j.humov.2018.06.013
Gronwald, T., Rogers, B., and Hoos, O. (2020). Fractal correlation properties of heart rate variability: A new biomarker for intensity distribution in endurance exercise and training prescription? Front. Physiol. 10:3389. doi: 10.3389/fphys.2020.550572
Hautala, A. J., Mäkikallio, T. H., Seppänen, T., Huikuri, H. V., and Tulppo, M. P. (2003). Short-term correlation properties of R-R interval dynamics at different exercise intensity levels. Clin. Physiol. Funct. Imag. 23, 215–223. doi: 10.1046/j.1475-097x.2003.00499.x
Hottenrott, K., and Hoos, O. (2017). “Heart rate variability analysis in exercise physiology,” in ECG time series analysis: Engineering to medicine, eds H. Jelinek, A. Khandoker, and D. Cornforth (London, UK: CRC Press), 245–257.
Huikuri, H. V., Mäkikallio, T. H., Peng, C. K., Goldberger, A. L., Hintze, U., and Møller, M. (2000). Fractal correlation properties of RR interval dynamics and mortality in patients with depressed left ventricular function after an acute myocardial infarction. Circulation 101, 47–53. doi: 10.1161/01.cir.101.1.47
Jamnick, N. A., Botella, J., Pyne, D. B., and Bishop, D. J. (2018). Manipulating graded exercise test variables affects the validity of the lactate threshold and VO2peak. PLoS One 13:e0199794. doi: 10.1371/journal.pone.0199794
Karavirta, L., Tulppo, M. P., Laaksonen, D. E., Nyman, K., Laukkanen, R. T., Kinnunen, H., et al. (2009). Heart rate dynamics after combined endurance and strength training in older men. Med. Sci. Sports Exerc. 41, 1436–1443. doi: 10.1249/mss.0b013e3181994a91
La Monica, M. B., Fukuda, D. H., Wang, R., Gonzalez, A. M., Wells, A. J., Hoffman, J. R., et al. (2018). Maintenance of Vagal Tone with Time-Release Caffeine, But Vagal Withdrawal During Placebo in Caffeine-Habituated Men. J. Caff. Adenos. Res. 8, 59–64. doi: 10.1089/caff.2017.0039
Mann, T., Lamberts, R. P., and Lambert, M. I. (2013). Methods of prescribing relative exercise intensity: physiological and practical considerations. Sports Med. 43, 613–625. doi: 10.1007/s40279-013-0045-x
Mendonca, G. V., Heffernan, K. S., Rossow, L., Guerra, M., Pereira, F. D., and Fernhall, B. (2010). Sex differences in linear and nonlinear heart rate variability during early recovery from supramaximal exercise. Appl. Physiol. Nutr. Metab. 35, 439–446. doi: 10.1139/h10-028
Meyer, K., Hajric, R., Westbrook, S., Samek, L., Lehmann, M., Schwaibold, M., et al. (1996). Ventilatory and lactate threshold determinations in healthy normals and cardiac patients: methodological problems. Eur. J. Appl. Physiol. Occupat. Physiol. 72, 387–393. doi: 10.1007/bf00242266
Michael, S., Graham, K. S., and Oam Davis, G. M. (2017). Cardiac Autonomic Responses during Exercise and Post-exercise Recovery Using Heart Rate Variability and Systolic Time Intervals-A Review. Front. Physiol. 8:301. doi: 10.3389/fphys.2017.00301
Muñoz, I., Seiler, S., Bautista, J., España, J., Larumbe, E., and Esteve-Lanao, J. (2014). Does polarized training improve performance in recreational runners? Int. J. Sports Physiol. Perform. 9, 265–272. doi: 10.1123/ijspp.2012-0350
Newell, J., Higgins, D., Madden, N., Cruickshank, J., Einbeck, J., McMillan, K., et al. (2007). Software for calculating blood lactate endurance markers. J. Sports Sci. 25, 1403–1409. doi: 10.1080/02640410601128922
Noakes, T. D., Gibson, A. S. C., and Lambert, E. V. (2005). From catastrophe to complexity: a novel model of integrative central neural regulation of effort and fatigue during exercise in humans: summary and conclusions. Br. J. Sports Med. 39, 120–124. doi: 10.1136/bjsm.2003.010330
Pallarés, J. G., Morán-Navarro, R., Ortega, J. F., Fernández-Elías, V. E., and Mora-Rodriguez, R. (2016). Validity and Reliability of Ventilatory and Blood Lactate Thresholds in Well-Trained Cyclists. PLoS One. 11:e0163389. doi: 10.1371/journal.pone.0163389
Peng, C. K., Havlin, S., Stanley, H. E., and Goldberger, A. L. (1995). Quantification of scaling exponents and crossover phenomena in nonstationary heartbeat time series. Chaos 5, 82–87. doi: 10.1063/1.166141
Platisa, M. M., Mazic, S., Nestorovic, Z., and Gal, V. (2008). Complexity of heartbeat interval series in young healthy trained and untrained men. Physiol. Meas. 29, 439–450. doi: 10.1088/0967-3334/29/4/002
Rincon Soler, A. I., Silva, L. E. V., Fazan, R. Jr., and Murta, L. O. Jr. (2018). The impact of artefact correction methods of RR series on heart rate variability parameters. J. Appl. Physiol. 124, 646–652. doi: 10.1152/japplphysiol.00927.2016
Runacres, A., Mackintosh, K. A., and McNarry, M. A. (2019). The effect of constant-intensity endurance training and high-intensity interval training on aerobic and anaerobic parameters in youth. J. Sports Sci. 37, 2492–2498. doi: 10.1080/02640414.2019.1644890
Sales, M. M., Sousa, C. V., da Silva Aguiar, S., Knechtle, B., Nikolaidis, P. T., Alves, P. M., et al. (2019). An integrative perspective of the anaerobic threshold. Physiol. Behav. 205, 29–32. doi: 10.1016/j.physbeh.2017.12.015
Sandercock, G. R., and Brodie, D. A. (2006). The use of heart rate variability measures to assess autonomic control during exercise. Scand. J. Med. Sci. Sports 16, 302–313. doi: 10.1111/j.1600-0838.2006.00556.x
Seiler, K. S., and Kjerland, G. Ø (2006). Quantifying training intensity distribution in elite endurance athletes: is there evidence for an “optimal” distribution? Scand. J. Med. Sci. Sports 16, 49–56. doi: 10.1111/j.1600-0838.2004.00418.x
Singh, M., Singh, M., and Banga, V. (2015). “Optimal ECG Sampling Rate for Non-Linear Heart Rate Variability,” in IJCA Proceedings on International Conference on Advancements in Engineering and Technology ICAET (Korea: ICAET), 22–26.
Stanley, J., Peake, J. M., and Buchheit, M. (2013). Cardiac parasympathetic reactivation following exercise: implications for training prescription. Sports Med. 43, 1259–1277. doi: 10.1007/s40279-013-0083-4
Stapelberg, N. J. C., Neumann, D. L., Shum, D. H. K., McConnell, H., and Hamilton-Craig, I. (2018). The sensitivity of 38 heart rate variability measures to the addition of artefact in human and artificial 24-hr cardiac recordings. Ann. Noninvas. Electrocardiol. 23:e12483.
Stöggl, T., and Sperlich, B. (2014). Polarized training has greater impact on key endurance variables than threshold, high intensity, or high volume training. Front. Physiol. 5:33. doi: 10.3389/fphys.2014.00033
Stöggl, T. L., and Sperlich, B. (2019). Editorial: Training Intensity, Volume and Recovery Distribution Among Elite and Recreational Endurance Athletes. Front. Physiol. 10:592. doi: 10.3389/fphys.2019.00592
Tapanainen, J. M., Seppänen, T., Laukkanen, R., Loimaala, A., and Huikuri, H. V. (1999). Significance of the Accuracy of RR Interval Detection for the Analysis of New Dynamic Measures of Heart Rate Variability. Anna. Noninvas. Electrocardiol. 4, 10–17. doi: 10.1111/j.1542-474x.1999.tb00359.x
van Wijck, K., Lenaerts, K., Grootjans, J., Wijnands, K. A., Poeze, M., Van Loon, L. J., et al. (2012). Physiology and pathophysiology of splanchnic hypoperfusion and intestinal injury during exercise: strategies for evaluation and prevention. Am. J. Physiol. 303, G155–G168.
Venhorst, A., Micklewright, D., and Noakes, T. D. (2018). Towards a three-dimensional framework of centrally regulated and goal-directed exercise behaviour: a narrative review. Br. J. Sports Med. 52, 957–966. doi: 10.1136/bjsports-2016-096907
Voss, A., Wessel, N., Sander, A., Malberg, H., and Dietz, R. (1996). Requirements on sampling rate in Holter systems for analysis of heart rate variability. Clin. Sci. 91(Suppl.), 120–121. doi: 10.1042/cs0910120supp
Weston, S. B., Gray, A. B., Schneider, D. A., and Gass, G. C. (2002). Effect of ramp slope on ventilation thresholds and VO2peak in male cyclists. Int. J. Sports Med. 23, 22–27. doi: 10.1055/s-2002-19267
Yeh, R. G., Chen, G. Y., Shieh, J. S., and Kuo, C. D. (2010). Parameter investigation of detrended fluctuation analysis for short-term human heart rate variability. J. Med. Biol. Eng. 30, 277–282. doi: 10.5405/jmbe.30.5.02
Keywords: detrended fluctuation analysis, ventilatory threshold, aerobic threshold, intensity distribution, intensity zones, endurance exercise, endurance training, polarized training
Citation: Rogers B, Giles D, Draper N, Hoos O and Gronwald T (2021) A New Detection Method Defining the Aerobic Threshold for Endurance Exercise and Training Prescription Based on Fractal Correlation Properties of Heart Rate Variability. Front. Physiol. 11:596567. doi: 10.3389/fphys.2020.596567
Received: 19 August 2020; Accepted: 18 December 2020;
Published: 15 January 2021.
Edited by:Clint Bellenger, University of South Australia, Australia
Reviewed by:Daniel Boullosa, Federal University of Mato Grosso Do Sul, Brazil
Lars Brechtel, MSB Medical School Berlin, Germany
Copyright © 2021 Rogers, Giles, Draper, Hoos and Gronwald. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.