Changes in the Complexity of Heart Rate Variability with Exercise Training Measured by Multiscale Entropy-Based Measurements

Quantifying complexity from heart rate variability (HRV) series is a challenging task, and multiscale entropy (MSE), along with its variants, has been demonstrated to be one of the most robust approaches to achieve this goal. Although physical training is known to be beneficial, there is little information about the long-term complexity changes induced by the physical conditioning. The present study aimed to quantify the changes in physiological complexity elicited by physical training through multiscale entropy-based complexity measurements. Rats were subject to a protocol of medium intensity training (n=13) or a sedentary protocol (n=12). One-hour HRV series were obtained from all conscious rats five days after the experimental protocol. We estimated MSE, multiscale dispersion entropy (MDE) and multiscale SDiffq from HRV series. Multiscale SDiffq is a recent approach that accounts for entropy differences between a given time series and its shuffled dynamics. From SDiffq, three attributes (q-attributes) were derived, namely SDiffqmax, qmax and qzero. MSE, MDE and multiscale q-attributes presented similar profiles, except for SDiffqmax. qmax showed significant differences between trained and sedentary groups on Time Scales 6 to 20. Results suggest that physical training increases the system complexity and that multiscale q-attributes provide valuable information about the physiological complexity.


Introduction
The study of system complexity is very challenging and has attracted much attention in the past few years [1][2][3]. Physiological complexity reflects the interoperability and correct functioning of regulatory processes as a whole, so the higher the complexity, the higher the system ability to adapt to different situations in daily life [4].
Heart rate variability (HRV) series, derived from the electrocardiogram (ECG) or arterial pressure signals, is one of the most important sources of information about system physiological status. Heart rate is actively controlled by the autonomic nervous system and can respond to many situations when the organism is challenged. A number of studies demonstrated that many indices extracted from HRV are powerful risk predictors of morbidity and death, for cardiac and non-cardiac diseases [5][6][7].
One of the most substantial challenges in the quantification of complexity from HRV time series is the difficulty in finding out a single measurement capable of doing this task consistently.
In other words, most of the complexity measurements are capable of extracting some properties that regard complexity itself, but none of them are enough to characterize all the complex traits of a system. Mono-and multi-fractal measurements [8,9], irreversibility estimations [10,11], symbolic methods [12,13], network analysis [14,15], as well as entropy-based approaches [16,17] have been proposed to infer the system complexity.
Multiscale entropy (MSE) is an important example of an approach that has been shown to be quite robust and consistent to characterize the system complexity from HRV time series. Like many other approaches, it has some limitations depending on the situation, and improvements or refinements have been proposed since MSE has emerged [18]. For example, the entropy estimator used in MSE (sample entropy) can be replaced by other estimators, such as permutation entropy [19], fuzzy entropy [20], distribution entropy [21], dispersion entropy [22], Rényi entropy [23] and bubble entropy [24], among others. Some other entropy-based proposals, such as entropy of entropy [25] and multiscale SDiff q (a measure of entropic differences) [26], are markedly different from the MSE original framework, although notably inspired by MSE.
Mild intensity aerobic exercise has been shown to improve several systemic functions and prepare the organism for sudden changes in the body. Experimental models using physical training have demonstrated that gaining physical conditioning, before an induced pathology, can reduce the disturbances caused by the disease [27,28]. In other words, physical conditioning seems to increase the system physiological complexity level. However, controversial findings have been reported about complexity and exercise, and scarce studies applied multiscale complexity approaches to identify how the aerobic training can increase the complexity in healthy subjects [29][30][31][32][33].
In the present study, we applied MSE and two other complexity measurements derived from MSE, namely multiscale dispersion entropy and multiscale SDiff q , to quantify the increase of complexity with physical training in experimental models of healthy rats. Results show that all measurements point to the same direction, but significant findings were obtained only with multiscale SDiff q .

Experimental Protocol
Male Wistar rats (210 g on average) were obtained from the Animal Care Facility at the Campus of Ribeirão Preto of the University of São Paulo. The animals' usage was according to the Ethical Principles in Animal Research adopted by the National Council for the Control of Animal Experimentation, approved by the Local Animal Ethical Committee from the School of Medicine of Ribeirão Preto of the University of São Paulo.
The study divided animals into trained (n = 13) and sedentary groups (n = 12). Since animals could have distinct initial physical conditioning, they were individually tested for maximum velocity (V max ). For the V max test, the animals were placed on a treadmill, with no inclination, and the speed was increased in steps of 3 m/s every 3 min. The stage where the animal fatigued, as well as the time spent on the incomplete stage were noted to calculate the V max of each rat [34].
The trained group underwent a physical training protocol on the treadmill with no inclination for 9 consecutive weeks, 5 days per week. The training protocol consisted of a medium intensity training that initiated at 50% of V max for 20 min and ended, at the ninth week, at 70% of V max for 60 min (Adapted from [35]). At the fifth week, the trained group underwent a new V max test to adjust the training protocol as some animals acquire physical conditioning quicker than others. The sedentary group followed the same protocol, but the treadmill was kept off.

Data Acquisition and Processing
Two to three days after the end of the physical training protocol, rats were anesthetized with a mixture of ketamine and xylazine (50 and 10 mg/kg, ip) and implanted with subcutaneous electrodes for ECG recordings. Two days after surgery, with the animals conscious and under free movement conditions, the electrodes were connected to a bioelectric amplifier (Animal BioAmp FE136, ADInstruments, Bella Vista, Australia), and ECG recordings were acquired (2 kHz) by an IBM/PC coupled to an analog-to-digital interface (ML866 PowerLab 4/30, ADInstruments, Bella Vista, Australia).
ECG was recorded during one hour, so that multiscale measurements could be confidently estimated from HRV series. ECG recordings were processed using computer software (LabChart Pro, ADInstruments, Bella Vista, Australia) that creates HRV series as the sequence of R-R intervals, i.e., the time interval between adjacent R waves. All ECG recordings were carefully inspected, and missing beat detections and artifacts were manually corrected. HRV series are 20,000 beats in length, on average.

Multiscale Sample Entropy
Multiscale sample entropy (MSE) is a widely-known procedure to quantify the irregularity of time series within a time-scale range [36,37]. The MSE algorithm consists of creating multiple scaled versions of the original time series and calculating sample entropy (SampEn) from each scaled time series.
Consider a time series given by u (1) Thus, SampEn is defined as [38]: where: and: The distance function d is given by: In Equations (3) and (5), 1 ≤ j ≤ N − m, j = i. In SampEn equations, m is the pattern length or embedding dimension and r is the tolerance factor assumed for similarity between samples.
To estimate MSE, multiple scaled versions of u are created by a coarse-graining procedure, where each element j in a τ-scaled series is defined by: Next, SampEn is calculated from each scaled time series u τ , resulting in a curve of entropy versus scale. It is worth noting that the higher the time scale (τ), the slower the dynamics that the scaled time series is representing. Importantly, the tolerance factor (r) of SampEn is kept fixed for all time scales (τ) in MSE.
In the present study, we calculated MSE with the most widely-used parameter setting, i.e., m = 2 and r = 15% of the original time series standard deviation. The maximum scale calculated was τ = 20.

Multiscale Dispersion Entropy
Multiscale dispersion entropy (MDE) is similar to MSE and also quantifies the complexity of time series [39]. However, instead of calculating SampEn for each scaled time series, dispersion entropy (DispEn) is used to estimate irregularity.
Consider the same time series given before (u). First, u is filtered by a normal cumulative distribution function (NCDF) with mean µ and standard deviation σ, resulting in a filtered time series u f , which ranges from 0 to 1. This procedure is intended to better treat outliers. Next, u f is mapped into c classes (1 to c), according to z c (j) = round(c * u f (j) + 0.5), a function that linearly maps the range [0, 1] to [1, c].
Each vector y m (i) represents a dispersion pattern. Considering that each value in y m can assume one of the c possible classes, there will be c m potential dispersion patterns.
The probability of occurrence of each dispersion pattern y m (i) in z c can be calculated as the number of times the pattern y m (i) appears on z c , divided by the total number of patterns in z c (i.e., N − m + 1). This procedure will result in a probability distribution for all possible dispersion patterns, p[y m (i)]. Finally, the DispEn is defined as the Shannon entropy of p[y m (i)] [22]: MDE uses the same coarse-graining procedure of MSE. Thus, MDE estimation consists of the creation of scaled versions of the original time series using Equation (7) and the calculation of DispEn from each scaled time series. However, the NCDF function applied to each scaled version is the same as that applied to the first scale, i.e., the original time series. This procedure has a similar effect of keeping r fixed at all time scales in MSE and can be achieved choosing the same µ and σ of the NCDF function at all scales.
Parameters of MDE were set as m = 2, c = 6 and maximum time scale τ = 20. NCDF was generated with µ and σ as the mean and standard deviation of the original time series, respectively.

Multiscale SDiff q
An alternative proposal for multiscale complexity measurement is the multiscale SDiff q analysis [26]. Although still inspired by MSE in the sense of multiscale analysis, multiscale SDiff q do not use the entropy values over scales directly to characterize complexity. Instead, differences of entropy between the time series and its uncorrelated version, i.e., surrogate data, are used to represent the complexity. The difference of entropy is evaluated for a range of q-values, which is a parameter derived from nonadditive mechanical statistics [40,41]. The so-called nonadditive q-entropy has three regimes, namely classic additive when q = 1, sub-additive when q > 1 and super-additive when q < 1.
SDiff q accounts for the difference between the SampEn q of a given time series and the mean SampEn q of a set of surrogate series. SampEn q is a generalization of SampEn inspired by nonadditive statistics, which introduces the nonadditive parameter q to SampEn. Its equation is given by [42]: where log q is defined as [43]: and [Z] + = max{Z, 0}. The definitions of U m (r) and U m+1 (r) are the same as presented in Equations (2) and (4) for SampEn.
To calculate SDiff q , one has to follow the steps: • From a given time series u, S surrogate series are generated from u. The surrogate series is obtained by simply shuffling u [44]; • Next, values A = U m (r) and B = U m+1 (r) are calculated from u; • Values of U m (r) and U m+1 (r) are also calculated from each surrogate instance, obtaining their mean values C = U m (r) and D = U m+1 (r); • Finally, SDiff q is defined by Equation (11) below: Both SampEn q and SDiff q are parametrized in q so that they represent a curve of entropy, or entropy difference, as a function of q. From SDiff q curves, three attributes (q-attributes) are obtained to characterize the time series dynamics, namely SDiff q max , q max and q zero . The SDiff q max represents the maximum value for SDiff q in the range of q. The q max and q zero represent the q-value where SDiff q finds its maximum and zero values, respectively. q max is the q parameter that gives the largest entropic separation between the actual time series and its surrogate versions, whereas q zero is the q parameter where original and shuffled dynamics have the same entropy. For more details on the calculation of q-attributes, please refer to [26,42].
The extension of SDiff q to a multiscale measurement is straightforward. Scaled versions of the original time series are created using the same coarse-graining procedure of MSE, given by Equation (7). Then, for each scaled time series, the SDiff q curve is calculated and q-attributes are obtained, so that it ends up with multiscale q-attributes.
Multiscale SDiff q parameters were set with the same values chosen for MSE, i.e., m = 2, r = 0.15 and maximum time scale τ = 20. The number of surrogate instances generated for each time scale was S = 20, and the nonadditive q parameter ranged from −2 to 2 to estimate the q-attributes.
It is worth emphasizing the fact that q-attributes represent the SDiff q behavior. Furthermore, the q parameter comes with the power law equation proposed for nonadditive entropy (q-entropy) [40,43]. Therefore, one can say that q max and q zero indicate where this power law results in maximum entropy differences regarding surrogates and where this difference is null (zero-crossing), respectively.

Statistical Analysis
We assessed mean MSE, MDE and multiscale q-attributes values in two range segments: short (1 to 5) and long (6 to 20) time scales. Those variables were checked for normality by the Shapiro-Wilk test. Differences between trained and sedentary groups were verified by Student's t-test or the Mann-Whitney rank sum test when required. Significance was assumed when p < 0.05.

Results
The curve profiles of MSE and MDE were very similar for both trained and sedentary rats ( Figure 1A,B). Likewise, no difference was found between the groups in the mean values of MSE and MDE grouped by short (1 to 5) and long (6 to 20) time scales ( Figure 1C,D), although for higher scales, there was a tendency of increasing differences among groups ( Figure 1A,B).
The curve profiles of q max and q zero were very similar to each other ( Figure 2B,C), which in turn were also very similar to MSE and MDE ( Figure 1A,B), regardless of the experimental group. On the other hand, those curves are entirely different from the profile of SDiff q max (Figure 2A). For q max and q zero , the curve values decrease for the first two or three scales; after that, they start to increase ( Figure 2B,C). However, in the case of SDiff q max , values increase for, approximately, the first six scales, and then, the values are virtually stable (Figure 2A). A significant difference was found between trained and sedentary rats in the mean q max at long time scales (6-20) ( Figure 2E). No difference was observed among groups in the mean SDiff q max ( Figure 2D) or mean q zero ( Figure 2F).  Bar graphs show mean q-attributes values, obtained from SDiff q max (D), q max (E) and q zero (F), grouped by short (1 to 5) and long (6 to 20) time scales. SDiff q max : maximal SDiff q ; q max : q value where SDiff q is maximal; q zero : q value where SDiff q is zero; HRV: heart rate variability. Bars represent the mean ± standard error. * p < 0.05 when compared to the trained group.

Discussion
The characterization of system physiological complexity from a univariate variable, such as HRV, is a hard task. Previous studies have reported on MSE as a powerful tool to assess the complexity of HRV [37,[45][46][47][48]. Many studies have proposed and evaluated modifications in MSE, given its success in characterizing complex dynamics. Some of them are based on the replacement of sample entropy by another entropy measurement, such as MDE, attempting to improve the accuracy of MSE in specific situations. In the present study, we applied MSE and MDE to account for the complexity changes due to physical training in rats. However, neither MSE nor MDE were able to detect any difference between HRV complexity from trained and sedentary rats.
On the other hand, multiscale SDiff q is a recent proposal of complexity measurements (q-attributes), inspired by MSE, but with a different theoretical background. It relies on nonadditive statistics and uses the difference of q-entropy between the actual and surrogate HRV time series to characterize the complexity. Interestingly, from all the multiscale measurements studied, only q max was able to distinguish the complexity of HRV between trained and sedentary animals. Moreover, the difference was found only at long time scales (6 to 20). Recent studies have pointed out that short time scales of MSE are more associated with the vagal control of HRV, whereas long time scales seem to be more related (although not exclusively) to the sympathetic control of HRV [46,49], reinforcing the existence of long-term memory in the components of the autonomic nervous system. Extending this interpretation to SDiff q , one could say that the difference between sedentary and trained HRV is more related to differences in the sympathetic control. This seems a reasonable assumption, given that (1) q-attributes use the same coarse-graining procedure of MSE to create the scaled time series and (2) physical training promotes, among other benefits, a lower sympathetic activity and modulation [50,51].
Even though there is a significant difference in q max between trained and sedentary groups, the difference is not huge. An interesting question to ask is: how much is changed in the physiological complexity with physical training? Another question would be: how do the interactions between physiological systems change in a physically trained animal? One has to bear in mind that all those multiscale measurements represent a general view of the system function. In other words, those complexity measurements extract the overall complexity of the system, which is the result of several mechanisms contributing to the homeostasis. Considering that the sedentary animals are healthy, a tremendous increase would not be expected in the complexity after physical training, given that most of the regulatory mechanisms are supposed to be already working at a high complexity level. Therefore, results suggest that systemic changes induced by physical training increase the system complexity to a slightly higher level.
The ability of those multiscale measurements to quantify the overall system complexity of HRV is a distinguishing feature. Many classical HRV indices seek to extract information related to the sympathetic or vagal autonomic modulation, not to mention that they are all linear models. Those indices are usually very sensitive to the environment and behavioral conditions and cannot represent the physiological complexity [4]. For example, during one hour of ECG recording, the rat may explore, sleep, groom, dig and other typical rat behaviors. All those situations will change the autonomic balance, and it is difficult to say what is the real sympathetic and vagal modulation of the rat during the whole one-hour period. On the other hand, applying the multiscale complexity measurements during the whole period, it was possible to identify that the dynamics of HRV has higher complexity in the trained rat compared to the sedentary one, even though the rat can change its physiological state several times during the recording. It is worth noting that all multiscale approaches were also applied to differential HRV series, but no difference was found between trained and sedentary animals, for any measurement [52].
The classical concept of entropy, e.g., SampEn and DispEn, relies on the quantification of the irregularity of a given series. The more irregular (unpredictable) the series, the higher the entropy. Thus, the entropy of any series is supposed to be lower than the entropy of its shuffled version (surrogate), even though the correlation properties of the dynamics were broken when samples are shuffled. However, with q-entropy, it is possible to achieve the same entropy values for both situations (q zero ). Therefore, if we consider the classical entropy (q = 1), surrogate data will always be assigned to a higher entropy value, but if we consider q near 0.5 (super-additive), the two dynamics will be assigned the same q-sample entropy. More interestingly, there are some values of q where the actual dynamics is assigned higher entropy regarding its surrogate (also for super-additive q). Hence, q max can be interpreted as the nonadditive parameter that maximizes the complex properties present in the actual dynamics.
In summary, results with multiscale SDiff q confirmed previous findings that q max and q zero provide similar, although not equivalent information, which is quite different from SDiff q max . Furthermore, MSE, MDE, q max and q zero presented very similar curve profiles, despite their different theoretical definitions, and q max was the only measurement that detected differences in the physiological complexity after physical training. There is no doubt that MSE represents a relevant tool for complexity analysis. This study reinforces that multiscale SDiff q is an alternative tool for characterizing the complexity of HRV time series, which can add information in some situations where MSE is not accurate enough. Multiscale SDiff q could also be used to help to characterize the complexity of HRV time series in different pathophysiological conditions, as well as in situations where the signal source is other than HRV.