Assessment of Heart Rate Variability Thresholds from Incremental Treadmill Tests in Five Cross-Country Skiing Techniques

This is an electronic reprint of the original article. This reprint may differ from the original in pagination and typographic detail. Author(s): Title: Year: Version: Please cite the original version: All material supplied via JYX is protected by copyright and other intellectual property rights, and duplication or sale of all or part of any of the repository collections is not permitted, except that material may be duplicated by you for your research use or educational purposes in electronic or print form. You must obtain permission for any other use. Electronic or print copies may not be offered, whether for sale or otherwise to anyone who is not an authorised user. Assessment of Heart Rate Variability Thresholds from Incremental Treadmill Tests in Five Cross-Country Skiing Techniques Mendia-Iztueta, Ibai; Monahan, Kristen; Kyröläinen, Heikki; Hynynen, Esa Mendia-Iztueta, I., Monahan, K., Kyröläinen, H., & Hynynen, E. (2016). Assessment of Heart Rate Variability Thresholds from Incremental Treadmill Tests in Five Cross- Country Skiing Techniques. PLoS ONE, 11(1), Article e0145875. https://doi.org/10.1371/journal.pone.0145875 2016 RESEARCH ARTICLE Assessment of Heart Rate Variability Thresholds from Incremental Treadmill Tests in Five Cross-Country Skiing Techniques Ibai Mendia-Iztueta 1 *, Kristen Monahan 1 , Heikki Kyröläinen 1☯ , Esa Hynynen 2☯ 1Department of Biology of Physical Activity, University of Jyväskylä, Jyväskylä, Finland, 2KIHU—Research Institute for Olympic Sports, Jyväskylä, Finland ☯These authors contributed equally to this work. *[email protected] Abstract The assessment of heart rate variability (HRV) thresholds (HRVTs) as an alternative of Ventilatory thresholds (VTs) is a relatively new approach with increasing popularity which has not been conducted in cross-country (XC) skiing yet. The main purpose of the present study was to assess HRVTs in the five main XC skiing-related techniques, double poling (DP), diagonal striding (DS), Nordic walking (NW), V1 skating (V1), and V2 skating (V2).Ten competitive skiers completed these incremental treadmill tests until exhaustion with a minimum of one to two recovery days in between each test. Ventilatory gases, HRV and poling frequencies were measured. The first HRV threshold (HRVT 1 ) was assessed using two time-domain analysis methods, and the second HRV threshold (HRVT 2 ) was assessed using two non-time varying frequency-domain analysis methods. HRVT 1 was assessed by plotting the mean successive difference (MSD) and standard deviation (SD) of normalized R-R intervals to workload. HRVT 1 was assessed by plotting high frequency power (HFP) and the HFP relative to respiratory sinus arrhythmia (HFP RSA ) with workload. HRVTs were named after their methods (HRVT 1-SD ;HRVT 1-MSD ;HRVT 2-HFP ;HRVT 2-HFP-RSA ). The results showed that the only cases where the proposed HRVTs were good assessors of VTs were the HRVT 1-SD of the DS test, the HRVT 1-MSD of the DS and V2 tests, and the HRVT 2-HFP-RSA of the NW test. The lack of a wider success of the assessment of HRVTs was reasoned to be mostly due to the high entrainment between the breathing and poling frequencies. As secondary finding, a novel Cardiolocomotor coupling mode was observed in the NW test. This new Cardiolocoomtor coupling mode corresponded to the whole bilateral poling cycle instead of corresponding to each poling action as it was reported to the date by the existing literature. Introduction Within the last three decades, blood lactate curves and gas exchange measurements from incremental exercise tests have been used for the assessment of endurance capacity and training PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 1/14 OPEN ACCESS Citation: Mendia-Iztueta I, Monahan K, Kyröläinen H, Hynynen E (2016) Assessment of Heart Rate Variability Thresholds from Incremental Treadmill Tests in Five Cross-Country Skiing Techniques. PLoS ONE 11(1): e0145875. doi:10.1371/journal. pone.0145875 Editor: Oyvind Sandbakk, Norwegian University of Science and Technology, NORWAY Received: May 22, 2015 Accepted: December 9, 2015 Published: January 4, 2016 Copyright: © 2016 Mendia-Iztueta et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. Data Availability Statement: The Ethics Committee of the University of Jyväskylä requires that the collected data is given only to authorized researchers as described in the Informed Consent that subjects signed. As a result, the original data cannot be made available. A minimal anonymized data set will be made available from the Research Institute for Olympic Sports for researchers who meet the criteria for access to confidential information. Please contact Esa Hynynen ([email protected]). Funding: This study was partly supported (salary and materials) by the European Regional zones, where two submaximal thresholds have been revealed [1]. However, the assessment of the two submaximal thresholds, being Ventilatory thresholds (VTs) or lactate thresholds (LTs), presents limitations such as laboratory dependence, high costs or an invasive nature [2]. At the beginning of the 21 st century, attempts to overcome these limitations began with the development of a new methodology for the assessment of submaximal thresholds based on Heart Rate Variability (HRV). The assessment of these thresholds, often referred to as HRV thresholds (HRVTs), is based on the links that HRV presents with the autonomous nervous system and respiratory sinus arrhythmia (RSA) [3–6]. The assessment of the first HRVT (HRVT 1 ) has been conducted with complex time varying frequency analysis methods [7–10] but also with simple time analysis methods [2,11]. The later ones are based on the virtual disappearance of the N-N (i.e. normalized R-R) interval variations, as seen in the trends of mean successive difference (MSD) and standard deviation (SD) of N-N, which represents vagal withdrawal. The assessment of the second HRVT (HRVT 2 ) has been conducted with time varying [8– 10,12–14] and non-time varying [13,15] frequency analysis methods. The frequency analysis of HRV decomposes the N-N time-dependent signal into its sinusoidal components, obtaining the so-called power spectral density (PSD). This PSD is then under-divided into low and high frequency (HF) bands. The HF band results from breathing (respiratory sinus arrhythmia or RSA) and mainly vagal cardiovascular control, but it is also related to Cardiolocomotor coupling (LOC) when summit to exercises where the upper body is involved [6,13,15]. The latest studies indicate that the most sensitive methods for HRVT 2 assessment in exercises where the upper body strikes contribute to propulsion are time varying methods that discard the power in the HF band (HFP) relative to the LOC (HFP LOC ) to solely use the HFP relative to the RSA (HFP RSA )[13,15]. To the best of our knowledge, the assessment of HRVTs has been conducted in varied exercise modes, including exercises using upper body movements, such as front crawl swimming [15] and ski mountaineering [13,14], but never in cross-country (XC) skiing. A characteristic of XC skiing is that, like swimming, it consists of different techniques with different relationships between breathing frequency (BF) and poling frequency (PF) [16,17], which would presumably imply different peak patterns in the HF band. The five main techniques used by XC skiers are Diagonal Striding (DS), Nordic Walking (NW), Double Poling (DP), V1 skating (V1) and V2 skating (V2). All techniques seem to share a tendency for an increase in BF-PF coupling with increasing workload intensities. Moreover, it seems that the PF-BF coupling is the weakest in the DS and NW techniques and strongest in the DP technique [13,16–20]. The main aim of the present study was to assess HRVTs in the five main techniques used by cross country skiers, and the secondary aim was to evaluate how the Cardiolocomotor interactions affected these assessments. For these purposes, five discontinuous incremental treadmill tests until exhaustion were conducted. Two time analysis methods were selected for the assessment of HRVT 1 . The HRVT 1 obtained from MSD was named HRVT 1-MSD , and the HRVT 1 obtained from SD was named HRVT 1-SD . For the assessment of HRVT 2 two non-time varying frequency analysis methods were used. The HRVT 2 obtained from HFP was named HRVT 2-HFP , and the HRVT 2 obtained from HFP RSA was named HRVT 2-HFP-RSA . Methods A group of competitive XC skiers performed five incremental tests until exhaustion in a time frame of two to three weeks with at least one day between tests. The tests were performed in randomized order. Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 2/14 Development Fund. The project number was A32156, and the funder's website was https://www. rakennerahastot.fi/. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. Competing Interests: The authors have declared that no competing interests exist. Abbreviations: BF & PF, Breathing and poling frequencies; DP, DS, NW, V1 & V2, Five different cross-country (XC) skiing related techniques: double poling, diagonal striding, Nordic walking, and V1 and V2 skating; fHF &pfHF, Frequencies in HF range corresponding to BF and PF; HF & HFP, High frequency and HF power; HFP LOC & HFP RSA −LOC, and RSA-related components of HFP, respectively; HRVT 1-MSD &HRVT 1-SD , First HRVT determined from MSD and SD trends; HRVT 2-HFP &HRVT 2-HFP-RSA , Second HRVT determined from HFP and HFP RSA trends; LF & LFP, Low frequency and LF Power; LOC, Cardiolocomotor coupling; MSD, Mean successive difference; RSA, Respiratory sinus arrhythmia. Subjects Ten competitive national level XC skiers (5 male and 5 female), with ages ranging from 19 to 30 years participated in the study (Table 1). All athletes were healthy and had just completed their XC ski race season. All subjects gave written informed consent agreeing to the terms and conditions of the study, and all possible risks and benefits were explained. The study was approved by the Ethics Committee of the University of Jyväskylä. The amount of subcutaneous fat was estimated from a four site skinfold measurement (m. biceps, m. triceps, os. scapula and, os. crista iliaca) [21]. Procedure Incremental Tests. NW, DP, DS, V1 skating, and V2 skating were the techniques used for the five incremental tests. All tests were performed on a large motor-driven treadmill (Rodby RL3500E, Rodby Innovations, Vänge, Sweden) located in an indoor laboratory setting. All participants wore a harness that was attached to a rope which hung from a metallic frame in the ceiling above the treadmill for safety. Marwe Classic 800 C and Marwe Skate 610 A roller skis (Hyvinkään Kumi, Hyvinkää, Finland) were used. One Way ski poles (One Way Sport, Vantaa, Finland) were equipped with special rubber tips designed to ensure an optimal grip while roller skiing on the motor-driven treadmill (Biomekanikk AS, Oslo, Norway). Before each incremental test, the subjects performed a 10 to 15 minute warm up with the same technique they were performing in the test on that day. In all protocols, workload increased every three minutes. At the end of every stage, the treadmill was stopped for 10 to 15 seconds for blood lactate samples (not used in this study) from the standing subjects. This small timeframe was included in the three minute stage. The highest 20 second mean values obtained during the tests were defined as peak values. The highest 20 second V O 2 mean value (V O 2 peak) among the five tests was defined as the V O 2 max value for the person if two or more of the following criteria were met: a HR within 10 beats of the age predicted maximal HR (220 –age); a respiratory exchange ratio superior to 1.1; or a plateau in V O 2 with increasing workload [13]. The PF was measured with a chronometer by timing ten poling cycles in the last minute of each stage. In the DP and the V2 protocols inclination was maintained constant at 3% and 5%, respectively, whereas speed increased 2 km/h with every stage from the initial speeds of 8 km/h for men and 7 km/h for women. In the DS and V1 tests, speed was maintained constant at 10 km/ h for men and 9 km/h for women, whereas inclination increased 2% with every stage from the initial 3% inclination. The NW tests followed a protocol that has been used for more than 30 years in Finland (Table 2)[22]. Assessment of First and Second Ventilatory Thresholds. A portable Oxycon Mobile spiroergometer (Viasys Healthcare GmbH, Hoechberg, Germany) was used for the gas measurements. The spiroergometer’s main sampling unit was lightweight and attached to the Table 1. Characteristics of the study subjects (n = 10). Characteristics Mean ±SD Males Females Age (years) 23.4 ±4.5 26.0 ±3.9 Height (cm) 180.2 ±3.1 166.6 ±4.6 Body mass (kg) 72.2 ±2.7 60.2 ±4.1 Fat (%) 10.3 ±1.5 19.4 ±2.9 V O 2 max (ml/kg/min) 73.2 ±1.1 59.9 ±3.5 doi:10.1371/journal.pone.0145875.t001 Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 3/14 upper back of the participants. Gas and volume calibrations were carried out twice before every test. The average breath-by-breath values of the last 90 seconds of every stage became the subject's values for the stage, a timeframe that complied with the interval lengths used for the HRV analysis [15]. VTs were determined by visual analysis of the breakpoints of different variables that were well documented in the literature [23]. For the determination of VT 1 ,V E-V O 2 and V CO 2 -V O 2 were plotted, whereas the determination of VT 2 , was assessed based on the V E- V O 2 and the V E-V CO 2 graphs. In the instances where the referred graphs did not allow determination of VTs, the graphs that plotted Ventilatory Equivalent of Carbon Dioxide (V E/ V CO 2 ) and Ventilatory Equivalent of Oxygen (V E/V O 2 ) with V O 2 were used. Assessment of First and Second Heart Rate Variability Thresholds. A Suunto t6 (Suunto Oy, Vantaa, Finland) HR monitor with a sampling frequency of 1000 Hz for the recording of the R-R intervals was used throughout each test. From the HRV recordings, the R-R interval series relative to the last 90 seconds of each step were extracted, as HF oscillations during these periods were considered stationary [15]. The R-R interval series were then automatically filtered with the software Firstbeat Health 3.1.1.0 (Firstbeat Technologies Oy, Jyväskylä, Finland) for the correction of the eventual artefacts and ectopic beats, which has been proven to be fast, accurate and reliable [24]. The normalized R-R (also named as N-N) time series, were then analysed by the software Kubios HRV 2.1 (Biosignal Analysis and Medical Imaging Group BSAMIG, Kuopio, Finland). HRVT 1 was assessed from two time-domain HRV analyses methods [11]. MSD and the SD were graphically plotted against workload. HRVT 1-SD and HRVT 1-MSD were set by visual interpretation at the point where there was no further decrease in the trends of these two parameters [11]. HRVT 2 was assessed with two time-domain HRV analyses methods [15]. The HFP and the HFP RSA were plotted against workload. One and sometimes two increases were noticed in the trends of both variables. HRVT 2-HFP and HRVT 2-HFP-RSA were set by visual interpretation at the point from which the last increase started [15]. The N-N series were first converted into an equidistantly sampled form, and the HR trend and part of the LF oscillations were then removed [15]. The PSD was estimated by an autoregressive model of order 12 or 18, depending on the type of test. An order 12 was applied to the DS and NW tests. However, an autoregressive model order 18 was considered to be more appropriate for the DP, V1, and V2 tests, because unlike the previously studied disciplines [13–15], these three disciplines use symmetric bilateral upper body strikes. The LF band was set between 0.04 and 0.15 Hz, and the HF band was set between 0.15 and 2.00 Hz. HFP was computed as the spectral power in the referred HF range, whereas HFP RSA was computed as the power from frequencies ranging from 0.04 Hz to a varying cut-off frequency which corresponded to the borderline between the RSA-related and LOC-related HF components [13,15]. The frequencies of the PSD peaks on HFP RSA and HFP LOC were named as fHF and pfHF, respectively. In the occurrence where HFP RSA and HFP LOC could not be divided, the frequency of the density peak was considered to be both fHF and pfHF. To verify that the HFP RSA and HFP LOC components were indeed related to respiration and locomotion, all of the spectrums (n = 321: 45 maximal tests with 5–10 stages each) were visually inspected and fHF and pfHF were compared to the corresponding BF and PF measured by spiroergometer, and chronometer. Table 2. Nordic Walking test protocol, women started from stage 1 and men from stage 2. Stage 1 2 3 4 5 6 78910 Inclination (%) 4 7 9 12 14 17 20 22 23 25 Speed (km/h) 6.0 6.0 6.6 6.6 7.0 7.0 7.0 7.2 7.6 7.8 doi:10.1371/journal.pone.0145875.t002 Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 4/14 The differentiation of the two HF band components in the five different tests occurred as follows. In the HF ranges of spectrums in the V2 tests, two clear PSD peaks emerged, the first peak corresponding to RSA and the second peak corresponding to LOC. In the DS tests, the main peak corresponded to RSA and the last and often non prominent peak's frequency (pfHF) corresponded to the double of PF. In the spectrums of the NW tests, two main peaks emerged. The first peak corresponded to RSA, and the second peak to PF. Lastly, in the vast majority of spectrums obtained from the V1 and DP tests, a single peak emerged in the HF band, preventing to differentiate two HF components. Therefore the decision was made to not split HFP in any of the spectrums of the V1 and DP tests, and the frequencies of their main peaks were decided to be both fHF and pfHF. Statistical Analysis A between methods agreement was used to evaluate whether there was an agreement or bias between the variables (i.e. VTs, BF and PF) determined from the reference methods (i.e. timer for PF and Ventilatory gases for the rest of variables) and the corresponding variables assessed from the alternative HRV-related methods (i.e. HRVTs, fHF and pfHF). The analyses included (a) an evaluation of the relationships between parameters using Pearson’s r correlation coefficients and linear regression lines, (b) an examination of the level of agreement using 95% limits of agreement according to Bland-Altman, and (c) a comparison of mean values using paired t tests. The agreement analysis between VTs and HRVTs was conducted with HR values because of its practical applicability. All the data are reported as the mean ± SD and the statistical significance was set at p 0.05 for all tests. All statistical analyses were performed using IBM SPSS Statistics 20 software (SPSS Inc, Chicago, IL, USA). The normal distribution of the data was verified by One-Sample Kolmogorov-Smirnov Test and the magnitude of the correlations was assessed according to Hopkins' scale [25]. Results In every subject, the test reaching the highest V O 2 peak value met the V O 2 max criteria, meaning that all subjects reached their maximal aerobic capacity in at least one of the five tests. Seven subjects reached the highest V O 2 peak value in the DS test and the remaining 3 subjects reached their highest V O 2 peak value in the NW test. Males' V O 2 max was 73.2 ± 1.1 ml/kg/min, whereas females' value was 59.9 ± 3.5 ml/kg/min. Agreements between Ventilatory and Heart Rate Variability Thresholds All VTs except the VT 1 of one subject's DP test were assessed. The R-R interval data from five maximal tests (2 DP tests, 2 V1 tests and 1 V2 test) were excluded due to excess of artefacts, and hence, none of their HRVTs could be assessed. Moreover, HRVT 1-MSD could not be assessed in one of the DS tests because of the abnormal behaviour of its trend. Additionally, the HRVT 2-HFP-RSA was not assessed in the DP and V1 tests because it was not possible to divide the HFP spectrums into two components. Besides the aforementioned cases where HRVTs could not be assessed, the reminders of HRVTs were assessed, including the HRVT 1 corresponding to the VT 1 that could not be assessed by the Ventilatory gas exchange method. A summary of the agreement between VTs and HRVTs can be seen in Table 3. First Threshold. With regards to the SD method, the NW test was the only test that showed statistical difference between HRVT 1-SD and VT 1 (p = 0.007). These variables where significantly correlated only in the DS (r = 0.77, p = 0.009) and NW (r = 0.68, p = 0.031) tests. The biases show that HRVT 1-SD slightly underestimated VT 1 in the DS and V1 tests (-1 and -2 Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 5/14 bpm, respectively), but VT 1 was overestimated in the other tests. The 95% limits of agreement were narrowest in NW and DS (-4/20 and -17/15 bpm, respectively), and was widest in DP (-32/51 bpm). Based on the regression lines, HRVT 1-SD explained 0.1% of the total variability of VT 1 in DP, 59.5% in DS, 46.1% in NW, 8.1% in V1 and 7.0% in V2. Therefore, it seems like the agreement between HRVT 1-SD and VT 1 was best in the DS test (Fig 1A) since their means were not significantly different, the bias was very small, they were strongly correlated and the regression line explained a good amount of the total variability of VT 1 . With regards to the MSD method, the DP and V1 tests were the only tests that showed statistical differences between HRVT 1-MSD and VT 1 (p = 0.039 and 0.049, respectively). These variables where significantly correlated in the DP (r = 0.88, p = 0.009), V2 (r = 0.77, p = 0.016) and DS (r = 0.75, p = 0.020) tests. All the biases showed a certain underestimation of HRVT 1-MSD over VT 1 ; the biases where smallest in V2 and DS (-3 and -6 bpm, respectively). The limits of agreement were also narrowest in the V2 and DS tests (-18/13 and -21/10 bpm, respectively). Based on the regression lines, HRVT 1-MSD explained 77.6% of the total variability of VT 1 in DP, 56.5% in DS, 5.8% in NW, 3.6% in V1 and 58.5% in V2. Therefore, it seems like the agreement between HRVT 1-MSD and VT 1 was best in the DS (Fig 1B) and V2 tests, since in both tests the variables were not significantly different, were strongly correlated and the agreement intervals were not that large. The equation of the best fitting linear regression line for the V2 test was y = 1.1007x - 19.305. Table 3. Mean (±SD) differences between the first (A & B) and second (C & D) thresholds determined by HRV and conventionally (spiroergometry). A) VT 1 -HRVT 1-SD Agreement (bpm) P R P 95% limits of agreement (bpm) Double poling 9±21 0.281 -0.03 0.955 -32 / 51 Diagonal striding -1 ±8 0.746 0.77 0.009 -17 / 15 Nordic walking 9±8 0.007 0.68 0.031 -4 / 22 V1 skating -2 ±13 0.619 -0.28 0.495 -27 / 23 V2 skating 6±10 0.112 0.26 0.492 -14 / 26 B) VT 1 -HRVT 1-MSD Double poling -9 ±9 0.039 0.88 0.009 -27 / 9 Diagonal striding -6 ±8 0.070 0.75 0.020 -21 / 10 Nordic walking -9 ±16 0.277 0.24 0.502 -37 / 25 V1 skating -9 ±11 0.049 0.19 0.654 -31 / 12 V2 skating -3 ±8 0.364 0.77 0.016 -18 / 13 C) VT 2 -HRVT 2-HFP Double poling -8 ±19 0.250 0.54 0.172 -45 / 28 Diagonal striding 11 ±13 0.027 0.80 0.006 -36 / 15 Nordic walking -4 ±9 0.180 0.44 0.207 -23 / 14 V1 skating -18 ±10 0.001 -0.28 0.510 -37 / 1 V2 skating -8 ±8 0.019 0.81 0.008 -24 / 8 D) VT 2 -HRVT 2-HFP-RSA Diagonal striding -9 ±13 0.055 0.66 0.036 -35 / 17 Nordic walking -1 ±7 0.818 0.82 0.026 -15 / 14 V2 skating -6 ±8 0.064 0.37 0.334 -21 / 10 Mean differences between the first threshold assessed conventionally (VT 1 ) and A) estimated by the Standard Deviation method of HRV (HRVT 1-SD ) and B) estimated by the Mean Successive Difference method of HRV (HRVT 1-MSD ). Mean differences between the second threshold assessed conventionally (VT 2 ) and A) estimated by the High Frequency Power method of HRV (HRVT 2-HFP ) and B) estimated by the Respiratory Sinus Arrhythmia related component of High Frequency Power (HRVT 2-HFP-RSA ). HRV = heart rate variability; r = correlation coefficient; p = statistical significance (p 0.05). doi:10.1371/journal.pone.0145875.t003 Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 6/14 Second Threshold. With regards to the HFP method, the DS, V1 and V2 tests showed statistical differences between HRVT 2-HFP and VT 2 (p = 0.027, 0.001 and 0.019, respectively). These variables where significantly correlated only in the DS (r = 0.80, p = 0.006) and V2 (r = 0.81, p = 0.008) tests. The smallest bias was found in NW (-4 bpm) and the biggest biases were found in V1 and DS (-18 and 11 bpm, respectively). The 95% limits of agreement were narrowest in V2, NW and V1 (-24/8, -23/14 and -37/1 bpm, respectively), and was widest in DP (-45/28 bpm). Based on the regression lines, HRVT 2-HFP explained 28.6% of the total variability of VT 1 in DP, 63.1% in DS, 19.1% in NW, 7.5% in V1 and 65.5% in V2. Therefore it seems that in none of the tests HRVT 2-HFP was a particularly good assessor of VT 2 . The lack of agreement between HRVT 2-HFP and VT 2 of the NW test is well seen in Fig 2A. With regards to the RSA method, none of the three tests showed statistical differences (p 0.05) between HRVT 2-HFP-RSA and VT 2 . These variables were significantly correlated only in the NW (r = 0.82, p = 0.026) and DS (r = 0.66, p = 0.036) tests. The biases of all three tests were negative, representing an underestimation of HRVT 2-HFP-RSA over VT 2 . The smallest bias was found in NW (-1 bpm) and the biggest in DS (-9 bpm). The limits of agreement was widest in DS (-35/17 bpm) and based on the regression lines, HRVT 2-HFP-RSA explained 44.1% of the total variability of VT 2 in DS, 48.3% in NW and 13.3% in V2. Therefore, it seems like the agreement between HRVT 2-HFP-RSA and VT 2 was close to optimal in the NW test (Fig 2B), since apart from being very strongly correlated and presenting no significant difference, the bias was Fig 1. Validity testing of the HRVT 1-SD (A1, A2) and HRVT 1-MSD (B1, B2) for the assessment of VT 1 during diagonal striding test. Bland-Altman (A1 & B1) plots the respective difference in heart rate (HR) between VT 1 -HRVT 1-SD and VT 1 -HRVT 1-MSD for each individual against their respective means. Dashed lines represent the limits of agreement corresponding to ±1.96 SD. The best fitting linear regression lines (A2 & B2) are also displayed together with their equations, enabling a prediction of future HR values. HRVT 1-SD & HRVT 1-MSD = heart rate variability thresholds as determined from standard deviation and from mean successive difference of normalized R-R intervals, respectively; VT 1 = first Ventilatory threshold. doi:10.1371/journal.pone.0145875.g001 Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 7/14 very small, the agreement interval relatively small, and the regression line explained a good level of the total variability of VT 2 . fHF-BF and pfHF-PF Agreements For all tests except for the V2 test, there was no significant difference between the mean of fHF and the mean of the timed BF. These variables were very strongly correlated (r 0.8) in DS, NW and V2, and strongly correlated (r = 0.6–0.8) in DP and V1 (Table 4A). In regards of the agreement between pfHF and the timed PF, the variables were significantly different and uncorrelated (p 0.05) in the DS and V1 tests. Other than this, there were no significant differences in the rest of the tests and the level of correlations were either very strong (NW and V2) or strong (DP) (Table 4B). Relationships between Breathing and Poling Frequencies In the five tests, there was a tendency for an increased coupling between the BF and PF as workload increased and there was always a higher PF than BF during the initial workloads. In the DP test, PF and BF seemed to be quite similar; PF was only significantly higher than BF in the first four stages. However, in the V1 and DS tests the BF trend significantly surpassed the PF trend in a cross-like pattern; significant differences between BF and PF where observed in Fig 2. Validity testing of the HRVT 2-HFP (A1, A2) and HRVT 2-HFP-RSA (B1, B2) for the assessment of VT 2 during Nordic walking test. Bland-Altman (A1 & B1) plots the respective difference in Heart Rate (HR) between VT 2 -HRVT 2-HFP and VT 2 -HRVT 2-HFP-RSA for each individual against their respective means. Dashed lines represent the limits of agreement corresponding to ±1.96 SD. The best fitting Linear Regression lines (A2 & B2) are also displayed together with their equations, enabling a prediction of future HR values. HRVT 2-HFP & HRVT 2-HFP-RSA = heart rate variability thresholds as determined from the high frequency power and from the high frequency power component related to respiratory sinus arrhythmia, respectively; VT 2 = second Ventilatory threshold. doi:10.1371/journal.pone.0145875.g002 Assessment of HRV Thresholds in XC Skiing PLOS ONE | DOI:10.1371/journal.pone.0145875 January 4, 2016 8/14