Wavelet Transform Analysis of Heart Rate to Assess Recovery Time for Long Distance Runners

: The diagnostics of the condition of athletes has become a ﬁeld of special scientiﬁc interest and activity. The aim of this study was to verify the effect of a long (100 km) run on a group of runners, as well as to assess the recovery time that is required for them to return to the pre-run state. The heart rate (HR) data presented were collected the day before the extreme physical effort, on the same day as, but after, the physical effort, as well as 24 and 48 h after. The Wavelet Transform (WT) and the Wavelet-based Fractal Analysis (WBFA) were implemented in the analysis. A tool was constructed that, based on quantitative data, enables one to conﬁrm the completion of the recovery process that is related to the extreme physical effort. Indirectly, a tool was constructed that enables one to conﬁrm the completion of the recovery process. The obtained information proves that the return to the resting state of the body after a signiﬁcant physical effort can be observed after two days entirely through the analysis of the HR. Certain practical measures were used to differentiate between two substantially different states of the human body, i.e., pre- and post-effort states were constructed. The obtained results allow for us to state that WBFA appears to be a useful and robust tool in the determination of hidden features of stochastic signals, such as HR time signals. The proposed method allows one to differentiate between particular days of measurements with a mean probability of 92.2%.


Introduction
The process of the differentiation of patients' conditions as a function of their physical effort and body positions is an important step in medical diagnostics regarding both healthy and unhealthy individuals [1,2]. Recently, the diagnostics of the condition of athletes has become a field of special scientific interest and activity. In the context of medical diagnostics, this is a group of healthy individuals that were subjected regularly to substantial physical effort. Therefore, it is important to evaluate any changes in their body functions based on pre-and post-effort medical assessment. This can be done by regular monitoring of body functions that are represented next as time series,

Study Description
The main purpose of the research is to investigate the changes in the circulatory system of long distance runners before and after an extreme physical effort. A group of volunteers from the Sport Centre of the Gdansk University of Physical Education and Sport have been selected as the test subjects. They were obligated to take part in regular measurements of their HR before the physical effort, just after the effort, as well as 24 and 48 h after the effort. The volunteers took part in a long distance run.
The run covered 100 km on a 3375 m loop with the total altitude difference not being greater than 2 m. The weather conditions during the run were as follows: air temperature of 4-5 degrees C and relative air humidity 86% with no wind. The run started at 7:20 in the morning and it lasted for approximately 10 h. The shortest recorded finishing time was 9 h and 11 min. and the longest was 12 h and 8 min., not including 3-min. breaks for collecting research data. Twelve hours before the race the runners ate dinner, and thenm after a night rest, at 6:30 in the morning, they ate a light breakfast. The runners were dressed according to their individual needs in sports shoes, cotton tracksuits, as well as gloves and hats made of natural materials.
During the run, each athlete wore on their chest a telemetric HR meter (Polar Electro, Kempele, Finland), which recorded their heart rate every 5 s. The time intervals of running and breaks between sections were measured by an electronic Timex clock (Zug, Switzerland, 2009) with the current reading visible on a board that was placed next to the starting line. During the run, they consumed individually prepared sets of beverages and nutritional products whenever they reported the need. The administration of food took place at designated stations. A standard meal consisted of low-mineralized water, energy drinks, sandwiches with cheese or ham, energy bars, and bananas. Runners consumed products with an average energy value of 1151 kcal per person, whereas their average energy demand was at the level of 5610 kcal, which is equivalent to oxygen consumption of 1.6 l/min. Blood samples were collected immediately before the run and after 25, 50, 75, and 100 km of running, as well as after 12 and 24 h of rest. The parameters of acid-base balance (arterial blood gas-ABG) were measured in arterialized blood samples that were drawn from the fingertip. During the analysis of blood parameters, the individual and mean values were taken into account, as well as the reference values (standards): pH (7.35 ÷ 7.45), base excess in extracellular fluids (BE ecf.; −2.3 to +2.3 mEq/L), bicarbonate level (HCO 3 ; 21 ÷ 27 mmol/L), oxygen partial pressure (pO 2 ; 75 to 100 mmHg), carbon dioxide partial pressure (pCO 2 ; 32 to 45 mmHg), oxygen saturation of haemoglobin (O 2 sat.; 95% ÷ 98%), and lactate level (0.5 ÷ 2.22 mmol/L). The arterial blood gas levels were determined with the use of the analyser type ABL 835 FLEX, which was produced by Radiometer Medical ApS (Brønshøj, Denmark). For the determination of the lactate concentration in blood samples, the enzymatic method with Randox reagent kit was used and reading of the extinction had been done at 37 degrees C while using a spectrophotometer type EPOLL 200, produced by Serw-med s.c. (Warsaw, Poland, 2006).
The HR of the runners was measured continuously for 12 min. and thanks to that 450-sample-long time series were collected. All of the HR signals were recorded based on the time intervals between successive R-waves of electrocardiographic signals (ECG). Their measurements were performed with NiccomoTM device (Medis, Ilmenau, Germany). Eight spot electrodes and a cuff manometer were used. The haemodynamic parameters were measured continuously and their values were recorded beat-to-beat. The quality of all the measurements was high and no relevant errors were noted. The HR of all volunteers measured before the run was around 64-72/min., while just after the run was around 141-158/min. For further analysis, only HR measurements from 13 male amateur runners were selected. It should be mentioned that the runners represent different training levels that cover between 25-35 years of running. Apart from running, they are also active in other sports, such as: swimming, cycling, or strength sports.
The study was approved by Independent Bioethics Commission for Research at Medical University of Gdansk. Written informed participation consent was obtained from all of the volunteers taking part in the programme. All of the experiments presented were carried out in accordance with the approved guidelines.

Wavelet Analysis
The Wavelet Transform (WT) employed in the current paper is a process of source signal decomposition and its following representation by a linear combination of some base functions, which are known in the literature as wavelets. As mathematical functions, wavelets are characterised by zero mean value, finite signal power, rapid decay, as well as a finite time length. These features of wavelets make them well localised bases in both time (or space) and frequency. For that reason, wavelets are especially useful in the analysis of biomedical signals, in which singular points are omnipresent and where the problem of decomposition and reconstruction of original non-periodic and non-stationary signals is extremely important.
Based on the literature available [34][35][36] related to the analysis of discontinuous, non-stationary signals of high variability, it was decided by the authors to use the Haar wavelet (synonymous with Daubechies 1 [37]) in the present work. These wavelets were used in wavelet decomposition processes that were subsequently used for HR assessment.
The first decomposition process was based on the use of the Continuous Wavelet Transform (CWT), which generated certain coefficients describing the resemblance between a chosen wavelet and the signal under investigation. The schematic diagram that is presented in Figure 1 illustrates the methodology employed by the authors based on the application of CWT, which led them to the calculation of the sums of wavelet coefficients.
WBFA was chosen as the second method. Its choice was supported by the fact that the analysed HR signals can be understood as time series (stochastic processes). On the other hand, stochastic processes can be well characterised by a physical quantity, known as a spectral power density (PSD).
Because the PSD represents the signal power that is associated with particular frequencies ω, it becomes possible to investigate the frequency content of stochastic processes as well as identify any periodic behaviour that is associated with them [38]. It should be also underlined, at this point, that HR signals are characterised by the PSD distributions that correspond to 1/ f processes (pink noise), which means, in practice, that their amplitudes are inversely proportional to frequency f [39,40].
The WBFA allows for an investigation of the HR correlation in a similar manner to that presented by Tan et al. [41]. In their work, they analysed, by fractal methods, the influence of various drugs on the heart rate variability (HRV) in the case of 10 healthy individuals. They pointed out the fact that increasing the time signal lengths during the assessment, even up to 1.5 h, had no significant influence on the quality of the results. Moreover, they suggested that the most appropriate signal lengths should stay around 15-20 min. This kind of information is vital for the material that is presented in this work, as all analysed time signals of HR include 450 samples, which is approximately equal to 12 min.
The self-similarity estimation of the signals under investigation was made based on the spectral exponent γ, calculated by WBFA, by: where d m,n denotes the wavelet coefficients that were obtained by the application of the Discrete Wavelet Transform (DWT) [42]. The wavelet coefficients are dot products of the time signal and a series of base functions. The schematic diagram that is presented in Figure 2 illustrates the methodology employed by the authors, but, based on the application of DWT using Daubechies base wavelet (db3) and eight decomposition levels, which led them to the calculation of the spectral exponent coefficients.
In order to decompose a given source time signal s(t) by DWT, the signal is first separated into a detail D1 and its approximation A1. In the following step, the number of time samples is reduced by half and next detail D2 and its approximation A2 are obtained. This process can be continued further; however, as the number of time samples is halved, such a discrete process can be carried out on a finite number of decomposition levels. This constraint is not present in the case of CWT that was described above earlier.
Differences between the scalograms representing two days of observation, thus a different condition of the patients, can indicate the capture of certain subtle, but significant changes. They can also confirm the robustness of the applied methodology and the numerical tool developed. In that context, the essence of the application of CWT or DWT comes from the curves that represent the sums of the wavelet coefficients as a function of the scale. In contrast, the essence of the application of WBFA comes from the spectral coefficients, which carry all necessary resulting information. On every stage of numerical simulations standard procedures for calculations of wavelet coefficients were used (in the case of the continuous as well as discrete wavelet transforms), which are integral parts of MATLAB software. For these calculations we employed MATLAB cwt function. Based on the documentation provided by the MATLAB software developer it can be stated that this function uses the norm L2.

Statistical Analysis
The curves and coefficients that wereobtained thanks to the analysis performed in the steps described above must be qualitatively assessed. In order to do that, the authors employed statistical tools that are commonly used by many researchers and they are a gold standard in analysing signal processing test results.
The fundamental analysis of the distributions of obtained values was performed by the use of the box plots, which allowed for the detection of important statistical features, including information regarding the mean and extreme values as well as the limitations that are related to 50% of test trials. This kind of approach is especially useful in the demonstration of results that were obtained by WBFA [29,42].
Apart from that, the normalised determination coefficient f it was used, which is a measure of the absolute fit between two curves c 1 (x) and c 2 (x). The normalised determination coefficient f it is very often used in technical problems [43,44]. It can be noted that, in the case of a perfect fit between curves c 1 (x) and c 2 (x), the value of this coefficient is equal to 100% according to: where c 1 (x) denotes the mean of c 1 (x) and denotes a vector norm.
A test of the hypothesis regarding the normal probability distribution of the set of n samples x 1 , x 2 , . . . , x n was performed by the application of the Shapiro-Wilk (SW) test [45], which is based on the following relation: where w i is a constant of the SW test, while x m is the sample mean.
Here, it should be added that the value of the SW test should be compared with the critical test value that is determined by the size of sample set.

Data
The volunteers were obligated to take part in regular measurements of their heart rate before the run, just after the run, as well as 24 and 48 h after the run. During the run, all of the volunteers had run with an average speed of 2.65 m/s and their heart rate reached an average rate of 143.25 bpm (Supplementary Table S1).
Each of HR time signals obtained during the measurement procedure [46] was initially pre-processed in such a manner that the resulting time signals had the same time lengths for signal comparison purposes. As a result, a set of 13 time signals was generated for the period of the four days measurement programme. Figure 3 presents exemplary HR signals. Typical heart rate signal recorded by electrocardiography (ECG). Particular HR signals were pre-processed in order to reduce their time length to 12 min. Different colours indicate particular a measurement day: 1st day, one day preceding the extreme physical effort; 2nd day, the day of the effort; 3rd day, one day following the effort, 4th day, two days following the effort. During the exercise time, the lactic acid level was measured. The samples were taken before the run, every 25, 50, 75, and 100 km of the run, as well as after subsequent 12 and 24 h of rest after the run. Figure 4 presents the pattern of lactic acid levels.  Table 1 presents the results regarding the physiological status of the subjects after 100 km run. The correlation between HR and lactic acid (LA) was 0.69.

Wavelet Decomposition of HR Data
In the first step, a HR time signal was pre-processed by the determination of the distribution of HR values in the assumed range. Figure 5 shows the result of signal processing in the case of one athlete. Figure 5a-d show the HR time signals for four consecutive days of the study before initial signal processing as raw HR signals, while, in Figure 5e-h, the corresponding distribution of these signals in the inverse domain, i.e., the number of samples vs HR. Such a procedure allowed for the authors to present the non-stationary and discontinuous raw HR signals, independent of time by the extraction of their global features. At the same time, the signals obtained are free of undesired edge phenomenon, which greatly simplifies the subsequent wavelet analysis by the selected wavelet as the beginnings and ends of the analysed signals are smooth and zero valued ( Figure 6). a b c d e f g h  In the second step of the current HR signal assessment procedure, CWT was applied in order to determine the necessary wavelet coefficients. This was carried out for 500 scales. The results of the wavelet transformation in the form of scalograms were presented in a graphical form, which represent the changes in wavelet coefficients as a function of the scale. Figure 7 shows the typical scalograms, based on the Haar wavelet, for four different days of measurements. Figure 7. Scalograms of heart rate (HR) signals that were obtained for Haar wavelet related to the all days of the measurement programme. Scalograms show the result of the Continuous Wavelet Transform, based on HR signals for a selected volunteer. Each plot indicates a particular measurement day: 1st day, one day preceding the extreme physical effort; 2nd day, the day of the effort; 3rd day, one day following the effort, 4th day, two days following the effort.
In the following step, the calculated wavelet coefficients were summed up related to the succeeding scales in order to reveal any quantitative differences between the particular HR signals measured. Figure 8 presents the calculated curves representing the state of a selected volunteer that was related to four measurement days.
Based on the results that are shown in Figure 8, it was decided to compare curves that represent all days of the measurement programme, for the two specific groups of volunteers. For each group, the mean distribution of HR values representing particular days was determined. The first group included four athletes, for which the results from all four days of the experiment were available. Figure S1, while the results that were obtained on their basis are shown in Supplementary Figure S2. The second group included all athletes, regardless of the completeness of HR data. Supplementary Figures S3 and S4 present the distributions of HR values and the corresponding sums of wavelet coefficients for the whole group of athletes, respectively. Figure 8. Curves representing sums of wavelet coefficients. Each curve was obtained as a result of summing up wavelet coefficients after the application of the Continuous Wavelet Transform of HR signals for a selected volunteer during four days of measurement programme. Different colours indicate a particular measurement day: 1st day, one day preceding the extreme physical effort; 2nd day, the day of the effort; 3rd day, one day following the effort, 4th day, two days following the effort.

Statistical Validation
The significance of the results of the research programme was also statistically verified. For that purpose, for each combination of the two curves, representing a particular measurement day, a set of coefficients f it was calculated, according to Equation (2), which contains information regarding their cross-similarity. Next, in order to check the hypothesis about the normal probability distribution of data in each set, the Shapiro-Wilk test was used. Figure 9 shows the probability density functions for each day of the measurement programme.

Wavelet-Based Fractal Analysis
The self-similarity of HR signals was evaluated by the calculation of the spectral exponent γ of 1/ f processes. It should be pointed out that, in the case of 1/ f processes, the observed changes in the spectral exponent γ can correspond to subtle changes in signal features. For that reason, and based on the results that were obtained thanks to the application of CWT, it was assumed that the spectral exponents γ related to different days of observation should also be different. Their values were established based on the calculation of wavelet coefficients by the use of DWT for 450-sample-long time signals. According to Equation (1), the spectral exponent γ denotes an angle between a trend line for the points that were obtained by the wavelet transformation, which is illustrated in Figure 10a.
The calculated values of the spectral exponent γ corresponding to the whole four day measurement programme and related to four selected volunteers are, as follows: 0.67, 0.81, 0.74, and 0.66 for day 1, 0.92, 1.03, 0.5, and 0.67 for day 2, 0.44, 0.79, 0.67, and 0.55 for day 3, 0.75, 0.82, 0.59, and 0.61 for day 4 (all corresponding to athletes 1-4). It can be seen that, in the case of three out of four volunteers, significant changes in the values of the spectral exponent γ can be noticed, which are associated with day 3, i.e., one day following the extreme physical effort. In order to analyse these changes in more detail, box plots were prepared, illustrating similarities and differences between the four volunteers, as shown in Figure 10b.  Additional chemical analysis was made in order to confirm the runners' changes of physical state. Figure 4 and Table 1 present the acid-base balance parameters that were gathered during the study. Changes in concentration of lactic acid in blood of the runners confirm the above presented claims that signs of an extreme physical effort can still be differentiated two days after the effort.

Discussion
The observed changes in the scalograms indicate the existence of differences between the physical condition of the volunteers taking part in the programme, as shown in Figure 7. However, this approach did not allow any quantitative evaluation of the features of measured signals, which made the clear evaluation of the physical condition of a particular volunteer impossible.
Sums of all wavelet coefficients, as presented in Figure 8, indicated changes between particular obtained curves. It also indicated that the differentiation of the physical condition of a particular volunteer by the use of appropriately constructed analytical criteria was possible. Moreover, the distance between the curves associated with post-effort-related states carry information about any recovery processes taking place.
The investigation of all the results obtained and that were related to particular measurement days allows for us to note that the sums of wavelet coefficients that were calculated and that corresponded to a particular effort state were concentrated in certain common domains. This is clear evidence of not only the existence of representative signal features for a set of signals measured on the same day, but also the successful identification of these features.
In most cases, the greatest changes were observed between the curves related to day 1 and day 2, as shown in Supplementary Figures S2 and S4. They differentiate these two effort-related states since the sets of curves associated with them are clearly separated. Similar trends can be observed in the case of day 2, day 3 and day 4. They indicate similarities between the physical conditions of the pre-effort and full-recovery states related to the sets of the curves representing day 1 and day 4. The curves representing day 3 and 4 are close to one another; however, they are fully separated from the other curves. Based on that, it can be stated that a full recovery of the human body takes place not earlier than after two full days following an extreme physical effort.
The statistical tests that were carried out confirmed the validity of the wavelet transform-based approach used in the paper. It can be stated that an appropriate selection of the base wavelet, as well as the scaling factor that is responsible for the number of scales, allowed for the determination of certain common features of selected sets of data that represent the physical condition of a particular volunteer. The first important indicators here are the means of the probability distributions that are related to different effort-related states: one day preceding the effortx 1 = 1.74 × 10 5 , the day of the effortx 2 = 2.40 × 10 5 , as well as two days following the effortx 4 = 1.83 × 10 5 . By knowing the valuex 1 , it is possible to differentiate day 1 in relation to day 2, with a 92% probability and 2.36 × 10 4 standard deviation. On the other hand, by knowing the valuex 2 , it is possible to differentiate day 2 from day 1 with a 93% probability and 2.23 × 10 4 standard deviation. Finally, based on the valuex 4 and its relation tox 2 , day 1 and day 4 can be differentiated and the probability of full recovery to the pre-effort state can be calculated as 91.6%. The remaining statistical data carried no useful information due to significant values of errors.
The obtained results allow for us to state that WBFA appears as a useful and robust tool in the determination of hidden features of stochastic signals, such as HR time signals. The statistical tests of the spectral exponent γ, as presented in Figure 9b and the above listed results, reinforce the statement. Based on these tests, the application of the wavelet transformation, together with fractal analysis, can be fully justified. The results obtained for the scaling factor related to the pre-effort state concentrate in a narrow interval of values.
Similarly to CWT, the key elements that influence the effectiveness of WBFA investigation of the HR signals is their quality and an appropriate choice of wavelet transform parameters. These HR time signals that contain too few samples limit signal decomposition levels, thus preventing the extraction of the signal details determining the value of the spectral exponent and/or scaling coefficient. The appropriate selection of the base wavelets in that process is equally important, as an inappropriate base wavelet may lead to great discrepancy and inaccuracy of the obtained results. Based on the simulation programme carried out, it can be stated that the best results were obtained for the 3rd-order Daubechies wavelet at the maximal decomposition level, which, in the current case, was equal to 8.
The HR as such has been extensively investigated for many years and has been found as a significant risk factor, especially in the case of cardiovascular events. Its predictive ability has been proven at rest and during and after various physical efforts. Both HR and HRV may reflect autonomic imbalance associated with different physiological and pathological states. However, the HRV provides information on quantities, i.e., HR and its variability. The question that arises here is what part of that information remains hidden, because of the non-linear relationship between the R-R intervals and the dynamic HR frequency in pre-and post-effort states that are related to extreme physical effort. It seems interesting to employ the wavelet transform-based analysis in the detection process of additional modulation of the HR. Thanks to this approach, it has been found that dynamic HR fluctuations have been present in the HR time signals measured until the second day following the extreme physical effort.
The results of the study, as revealed by the wavelet analysis, stand in opposition to common medical knowledge. A previous study, which concerned the influence of a 100 km run on the ANS parameters, revealed that, during the next day after the run, all of the measured parameters came back to their initial values [16]. Therefore, there was a clear indication of the lack of noticeable changes during further periods of recuperation. This observation is backed by several other studies [5][6][7][8], in which various cardiovascular indices measured ceased to change during the first day after a long run. In turn, our study revealed that the lack of change in the case of dynamic HR fluctuations, and likely many other cases, is associated with the selection of the analysis method, rather than the true response of a human body.

Conclusions
The technical aspect of the research programme that was carried out by the authors allowed for the construction of certain practical measures used to differentiate two substantially different states of the human body, i.e., pre-and post-effort states. In order to assist in the analysis, the Wavelet Transform and the Wavelet-based Fractal Analysis were implemented. Indirectly, the authors also managed to construct a tool, entirely based on quantitative data, which enables for one to confirm the completion of the recovery process related to extreme physical effort that is entirely based on the quantitative data. The information that was obtained by the tool presented proves that the return to the resting state of the body after a significant physical effort can be observed after two days entirely through the analysis of the heart rate.
The main conclusion of the study was the slow return of a runner's body to its original (resting) state after a long, extreme physical effort on the next day after the run. Previous studies concerning the influence of 100 km run on the ANS parameters revealed that, on the next day after the run, all of the measured indices come back to their initial values. The studies suggested that the 100 km run seems to be safer and, hence, more advisable to amateur athletes than, e.g., a marathon, which has become a popular mass-event. In our present study, a slow return of a runner's body to its original state is demonstrated, which reveals the usefulness of our novel analytic techniques and shows that the full recovery of measured parameters was observed within 48 h after the 100 km run.
Supplementary Materials: The following are available at http://www.mdpi.com/2079-9292/9/12/2189/s1, Figure S1: HR value distributions for the group of four volunteers, Figure S2: Comparison of curves representing sums of wavelet coefficients for the group of four volunteers, Figure S3: HR value distributions for the entire group of volunteers, Figure S4: Comparison of curves representing sums of wavelet coefficients for the entire group of volunteers, Table S1: Heart rate and running speed of amateur long distance runners (mean ± standard deviation) measured during subsequent mileages on a distance of 100 km.