Heart Rate Variability Alterations During Delayed-Onset Muscle Soreness Inducing Exercise – with Piezo2 Interpretation

1. Introduction

Delayed-onset muscle soreness (DOMS) is an enigmatic pain condition affecting basically everyone. No unequivocal mechanism theory for the cause exists for more than 125 years, when the American Physician Theodore Hough first depicted it [1]. The most circulated theories are the lactic acid, muscle spasm, inflammation, connective tissue damage, muscle damage and enzyme efflux theories [2]. It is known that unaccustomed and/or strenuous eccentric and/or isometric contractions could induce DOMS. The following symptoms are known to be the consequence of DOMS, such as muscle stiffness, swelling, loss of force-generating capacity, reduced joint range of motion, and even more importantly delayed onset of pain sensation and impaired proprioception [3]. The transient pain sensation develops only about 8 hours, tops at about 1 or 2 days, and subsidies within 7 days after the DOMS-inducing exercise bout [4]. So, DOMS is clearly distinguished from pain that is experienced during or right after exercise [5]. Noteworthy that heart rate variability (HRV) is impacted by nociceptive C-fibers and pain [6], but in the case of DOMS pain is not present during the inducing exercise bout and develops only hours later [4].

HRV measures the variation in interbeat intervals from one heartbeat to the next, reflecting changes in the heart rate (HR) over time. Heart rhythm is primarily governed by the synchronized firing of pacemaker cells in the sinoatrial node (SAN) of the heart muscle, which initiates the cardiac cycle. The activity of these pacemaker cells, that generate the so-called orderly “funny” currents by primarily the hyperpolarization-activated cyclic nucleotide-gated (HCN) channels, is regulated by the autonomic nervous system (ANS). The complex actions of the underlying interdependent regulatory mechanisms give rise to the variability of the length of the cardiac cycle over different time scales, supporting the optimal performance of the heart under homeostasis. HRV measurement is one diagnostic tool to detect ANS changes and it is often used in both recreational and competitive sport activities. Moreover, this technique is proven to be effective to follow the proper course of regeneration after strenuous exercise activities. Accordingly, the return of the parasympathetic drive to pre-exercise level, and even higher magnitudes, is of interest on the path to regeneration. Anaerobic exercise delays the return of the parasympathetic tone, or autonomic control, substantially compared to aerobic exercise [7], while strenuous exercise delays this return by two days [8] which is suggested to overlap with the time window of Piezo2 channelopathy [9]. Higher HRV is often associated with better fitness and heart adaptability, while reduced HRV is commonly linked to various pathological conditions, including congestive heart failure, diabetic neuropathy, mental disorders and cancer. As it has been revealed in previous publications about the analysis of HRV data per se, the time- and the frequency-domain, as well as the non-linear HRV parameters, prominently entropy, show the same phenomenon from different viewpoints, hence they are strongly interrelated, but they still carry independent information from each other [10] .

The data on how the ANS exerts its modulation during exercise is accumulating. At low intensities the sympathetic loading increases, while the parasympathetic loading decreases. This parasympathetic tone declines to a point under prolonged sympathetic loading when it is almost entirely withdrawn. Thereafter the sympathetic drive sustains its activity, however in decrements at higher intensities. The current authors hypothesize that a line of demarcation exists when the ANS regulation of exercise activities flips into a disordered state, or a transient point of no-return, reflecting the initiating microdamage during exercise, leading to DOMS. It is important to note that the condition of DOMS is prone to higher injury risk, as the prevailing view of the scientific community [11]. Certain sports and sports modalities, like eccentric/isometric exercise, are highly affected by DOMS. It has been long know that handball is one of those sports that comes with high-injury rate, especially amongst adolescent players [12], including non-contact injuries as DOMS is considered [13]. This is why the current study was conducted on young elite handball players. Therefore, if the aforementioned transient point of no return could be detected by HRV measurements then HRV monitoring devices would gain new diagnostic relevance, serving coaches, athletes and many others. Moreover, the detection of the suspected transient point of no return would help coaches and athletes to change sports modalities from eccentric/isometric exercise to concentric one in order to avoid harsher injury risk.

In support, evidence is accumulating that DOMS starts off even during the inducing unaccustomed and/or strenuous exercise activity [11,14,15,16] and not 8 hours after when the pain arises. HRV alterations is thought to be modulated by nociceptive stimulation, however in DOMS the nociceptive nerve fibers reflect pain only hours after the DOMS-inducing exercise [4]. Indeed, a recent pilot study showed orthostatic imbalance right after DOMS-inducing exercise, detected by an orthostatic stress test that is indicative in reference to ANS dysregulation [15]. This finding not only highlighted that ANS disbalance precedes pain evolvement [15], but increased our suspicion that a transient point of no return should exist during DOMS-inducing exercise activity. However, it is important to note that to our knowledge, no such “primary damage” point has been detected by HRV monitoring device.

A new neurocentric theory of DOMS hypothesizes that DOMS is a bi-phasic, bi-compartmental non-contact injury mechanism where the primary damage occurs on Type Ia proprioceptive fiber terminals in the muscle spindle under an acute stress response (ASR) time window [17], in the form of an acquired Piezo2 channelopathy [18]. Noteworthy that the bi-phasic injury mechanism in relation to DOMS have been introduced earlier, but in extrafusal muscle territories [19,20]. However, DOMS could be induced without muscle damage and only at high exercise intensities when muscle damage may evolve [21,22]. Moreover, this DOMS theory not only incorporated the importance of the line of demarcation between good stress and bad stress reflected in remodeling and the breach of remodeling correspondingly [18], but also highlighted the role of oxidative stress in DOMS mechanism [17,23]. Furthermore, oxidative stress is shown to be one of the most significant influencing factors of HRV by its time and frequency domains, however not by nonlinear parameters [24]. Interestingly, the conditions where Piezo2 channelopathy is suspected, including DOMS, as the primary damage (the one single initiating cause of aging) are always associated with autonomic disbalance [18] and oxidative stress [25]. Accordingly, the current authors rather followed indirect neurocentric tracing by using HRV measurement tools. Moreover, this onsetting microdamaging event is proposed to display quantum mechanical properties, initiated by proton affinity switch on Piezo2 function that impairs a novel ultrafast proton-based long-distance proprioceptive signaling in the nervous system [18,26]. Indeed, Piezo2 is the principal mechanotransductory ion channel responsible for proprioception, as was shown by the team of Nobel laureate Ardem Patapoutian [27], and in fact DOMS comes with impaired proprioception [3,14]. In support, recent combined, including both transcranial and peripheral, electromagnetic stimulation (also called paired-associative electromagnetic stimulation) showed a positive therapeutic impact after DOMS-inducing exercise on HRV measures [28], and on muscle activation, force and performance measurements [29] as detected by HRV analysis and surface electromyography devices respectively. Additional confirmation of the Piezo2 channelopathy involved neurocentric DOMS theory [13] and the aforementioned therapeutic paired-associative electromagnetic stimulation DOMS studies [28,29], that Piezo2 is found to be the underlying mediator of magnetic stimulation [30] as was theorized earlier [13]. Hence, our deeper understanding of Piezo2 function would not only help us in the treatment of DOMS, but in prevention as well.

Low intensity exercise does not appear to alter low frequency (LF) power of HRV, in contrast to medium to high intensity exercise that decreases LF power to close to zero value as the sympathetic load increases [31]. Furthermore, a recent study revealed that time- and frequency-domain HRV parameters are interrelated by nonlinear indices, notably with sample entropy (SampEn), which is a very sensitive measure of ANS control [10,32]. An earlier paper proposed that the LF power of HRV reflects the activity level of Piezo2 in the baroreceptor of the circulation and the heart [33]. Furthermore, a recent paper suggested that Piezo2 is an ultradian sensor and the only ion channel capable of initiating the theorized novel unaccounted ultrafast proton-based long-distance proprioceptive signaling in the nervous system [18], consequently it is likely manifested in the ultradian rhythm generation [34]. In support, ANS controls ultradian fluctuations through the baroreflex sensitivity and even more importantly it is HR dependent [35]. In support, conditional PIEZO2 and PIEZO1 gene knock-out mice not only endures baroreflex failure, but essentially loses blood pressure and heart rate control [36]. Correspondingly, the current authors assign Piezo2 to ultradian, while Piezo1 to diurnal control of blood pressure and HR regulation. Accordingly, Piezo1 indeed participates in the diurnal regulation of certain homeostatic processes [37]. It has been theorized that the ultradian clock is present at the hippocampal end of the novel unaccounted ultrafast proton-based long-distance proprioceptive signaling arised from rhythmic bacteria induced Piezo2 containing enterochromaffin cells-neural complex as the backbone of the microbiome-gut-brain axis [34]. What’s more, this ultradian rhythm is likely consolidated in learning and memory through the proposed hippocampal clock in synchrony with rapid eye movement (REM) sleep [34] with the contribution of Piezo2 [25]. In favor, the association between the ultradian rhythm and REM sleep has been long observed [38] and this communication is suggested to involve the temporal pattern of the basic rest-activity-cycle (BRC) [39].

After all the goal of our study was to detect alterations in HRV measurements already during DOMS-inducing exercise bout that would reflect upon the suspected initiating neural primary damage, or the transient point of no return in DOMS onsetting mechanism, despite that pain evolves hours after the inducing exercise.

2. Materials and Methods

2.1. Participants

In this study, the sample included 19 adolescent male U16 elite handball players. At the time of the data collection procedure every individual was the player of the Hungarian National Academy of Handball (Balatonboglar, Hungary). Basic anthropometric data of the sample: age: 15.65±0.65 yrs; body mass: 73.68±11.14 kg (Inbody 270, InBody Co., Ltd, South Korean); height:183.53±6.62 cm (SECA 213, Seca GmbH & Co. KG, Germany). Only those players were included in the sample who did not have lower limb pain at the time of the data collection or any kind of lower limb injury prior to the measurements. Moreover, those players were also excluded whose ECG signal was not accurate enough during the dynamometer test or if the electrode fell off from the athlete’s body. Prior to the measurements, the players were informed about the data collection protocol and the possible risks regarding the study. Due to the players young age, every individual’s parent signed a written consent about the the players participation in the study according to the Declaration of Helsinki Ethical Research Principles. The study was approved by the Science Ethics Committee of the Hungarian University of Sports Science (Ethical Approval Number: TE-KEB/18/2022).

2.2. Measurement Procedure

Before and after the data collection the participants filled out a visual analogue scale (VAS) [40,41] questionnaire regarding subjective determination of pain in the lower limb. Moreover, eight and twenty-four hours after the DOMS-inducing protocol the players were asked again about the pain by VAS scale. The VAS scales were evaluated as follows: 0-no pain; 1-3 mild pain; 4-6 medium pain; 7-10 severe pain. The athletes were subsequently divided into two groups: Those of DOMS values of lower than 5 were classified into the ”Low-DOMS”, while those of above 6 into the ”High-DOMS” group. The data of the two groups are also analyzed separately in the paper.

Prior to the dynamometer test measurements, a universal warmup protocol was conducted for every participant including a 5 minutes of bicycle ergometer task with minimal load and 5 minutes of basic gymnastic exercises performed without any additional weight. After the warmup, the measurements were executed with a Humac NORM isokinetic dynamometer (CSMi, Computer Sports Medicine, Inc, Stoughton, MA) controlling the movement of the knee in the dominant leg of the participants. All the participating athletes were familiar with the Humac NORM isokinetic dynamometer device therefore no prior training exercise was necessary. 10 participants executed a concentric contraction (CON) protocol and nine an eccentric (ECC). The measurement procedure in the isokinetic dynamometer was carried out as follows: 15 repetitions in 6 sets of CON or ECC while the motor moved the limb with an angular velocity of 60 deg/sec between a knee angle of 10-100 degs. Between the sets a 1 minute rest interval was included. The participants were instructed to generate maximal intensity of resistance while performing the ECC protocol and similarly maximum intensity of push while the CON protocol was carried out.

2.3. Data Acquisition and Preprocessing

Electrocardiograms (ECG) were recorded by MOX3 wearable ECG-Actigraph devices (Maastricht Instruments, Maastricht, The Netherlands), which non-invasively detect single-channel ECG (sampled at 400 Hz frequency), and the participant’s motor activity. The dataset included ECG signals resampled to 2000 Hz using cubic spline interpolation to enhance R-peak detection accuracy. [Accelerometry.eu. MOX3 Activity Monitor + ECG Logger. Available online: https://www.accelerometry.eu/products/wearable sensors/mox3/ (accessed on 21 May 2024).]

2.4. R-Peak Detection and RR Interval Calculation

R-peaks were identified using a custom detector, from which RR intervals and heart rate (HR) time series were derived. Artifactual beats and ectopic intervals were removed through a sequential segmentation procedure combined with a 30-points moving-median-based outlier rejection MATLAB routine “isoutlier” (MathWorks, 2020).

2.5. Spectral and Complexity Analysis

Interpolated RR intervals, RR-differences, and HR sequences were generated at 2 Hz using linear interpolation. Time-frequency analysis was performed using short-time Fourier transform (STFT) with overlapping windows (200-sample length, 10-sample step). Total power (TP), low-frequency (LF) and high-frequency (HF) spectral components were extracted for each window. Additionally, non-linear dynamic complexity was quantified using SampEn. Entropy, in general, characterizes the unpredictability of the RR-interval time series (Higher entropy = more complex signal), and among its various measures, SampEn is the most frequently used one. More formally, SampEn(m, r, N) = ln(A/B), where A=N(d[RR_m+1(i), RR_m+1(j)] < r) and B=N(d[RR_m(i),RR_m(j)] < r). N is the number of vector pairs, while in our case the embedding dimension m=2, and tolerance r=0.2∙SD (standard deviation) [42]. Root mean square of successive differences (RMSSD) and standard deviation (SDNN) were also computed in each window.

2.6. Significance Analysis

In order to quantify the group-level differences between the calculated HRV parameters of Low-DOMS and High-DOMS groups, a 2D Kolmogorov-Smirnov test was applied. Subsequently, we calculated the mean HRV values of the two groups as a function of HR, then the difference between the respective HRV(HR) curves representing the two groups were subtracted to visualize differences.

3. Results

3.1. The HR-Dependence Approach

In our study, we considered typical time-domain HRV parameters (RMSSD and SDNN), frequency-domain parameters (LF and HF power), and the non-linear parameter, Sample Entropy (SampEn), with special respect to the latter, as a function of HR. The primary reason for this is that all the HRV-parameters are dependent on the actual HR they measured at, and this dependence is generally pivotal, namely, it usually masks the effects of all the other factors on HRV [32,43,44]. On the example of RMSSD, Figure 1 illustrates our concept of depicting the HRV parameters of a selected typical individual, as a function of HR (see Figure 1b), which considerably simplifies the representation of the otherwise complex data set of HR(t) and HRV(t) functions (such as in Figure 1a), as well as the interpretation of the observed phenomena according to the ANS effects.

3.2. The Effect of Exercise on the HRV Parameters of the Athletes

In first approximation, we analyzed all the participants’ data to highlight differences in the various HRV parameters between time domains before and during the targeted exercise sections. Figure 2. shows the collected HRV data of all the athletes, as a function of HR, distinguishing the data got from the preparation/warm-up (Period A), and the targeted physical exercise (Period B) periods (blue and red dots, respectively). During period A, SampEn decreased gradually as a function of HR (blue dots in the figure), up to a point around 120 bpm (turning point, HR_tp), from where the entropy was to increase again. This type of entropy-HR relationship was observed for all athletes, but the HR range during warm-up was individual for each athlete, so it shows different parts of the entire curve, for each person. It can also be seen that the SampEn value covers a wide range at a given HR. We also plotted the SampEn(HR) function during exercise (Period B), where the load also includes the rest between the 6 subsequent sets. During exercise, it can be observed that the entropy was to decrease compared to the state during period A. In other words, increased exercise loading characteristically resulted in lower entropy at a given HR. This is true for the entire HR-range examined, including the range below and above the turning point. Another interesting observation is that the entropy at a given HR was not only lower, but its value remained within a much narrower range than in the case of no load (Period A) (Figure 2a).

The RMSSD, SDNN, LF and HF parameters show a different pattern than SampEn, but very similar to each other, so we describe them together. The parameters measured during the warm-up (Period A) decreased linearly in a semi-logarithmic representation as a function of HR, which confirms their exponential dependence [43]. This was also true for each athlete individually. However, under intensive exercise load (Period B), the values of each of these HRV parameters increased compared to the reference set at Period A. The increase in the time- and frequency-domain parameters compared to the reference curve (Period A) due to exercise (Period B), however, suggests that ANS activity was increased under exercise over the entire examined HR range (the sympathetic - parasympathetic balance is upshifted compared to the normal state) [44].

3.3. DOMS-Related Differences in Typical HRV Parameters

As a next step, we performed an analysis to reveal whether the degree of DOMS had an interrelationship with HRV-parameters determined before or during exercise sessions. First, the HRV(HR) functions were determined for the preparation and warm-up period (Period A) from the RR time series data of both the Low-DOMS and High-DOMS groups, separately, then the natural logarithms of the respective group-averaged HRV(HR) curves were subtracted from each other (always the Low-DOMS from the High-DOMS), after a proper normalization to the number of individuals in the two groups (Figure S1). Although the difference curves were noisy, a 2D Kolmogorov-Smirnov test could make a clear distinction between the data sets of the two groups. In case of all parameter, H=1, meaning that the two distributions are different (SampEn: p=10^-19, SDNN: p=10^-26, RMSSD: p=10^-24, LF: p=10^-37, HF p=40^-23). A difference between the respective HRV(HR) curves representing the two groups could be visualized for all the ΔHRV parameters as a function of HR, including an apparent sign change around the turning-point HR value (HR_tp).

Namely, during Period A, at HR values lower than HR_tp, in average, there appear greater values of all HRV parameters for the High-DOMS group than for the Low-DOMS one, while above the turning point, the relation reverses. Performing the same type of evaluation for the data obtained during exercise (Period B) (Figure S2, S3, and Figure 3), a similar picture arose for all HRV parameters but the entropy, where, interestingly, the difference curve describing the ΔSampEn(HR) function was practically mirrored to the abscissa, inasmuch as at HR values lower than HR_tp, the averaged entropy values during exercise are somewhat higher for the Low-DOMS group, while above HR_tp, the entropy values of the High-DOMS group clearly dominate, increasingly with growing HR (Figure 3).

4. Discussion

The current study analyzed HRV measures during DOMS-inducing exercise. Our objective was to find alterations in these HRV parameters even during the DOMS-inducing exercise despite pain evolves only post-exercise. We used an isokinetic dynamometer to induce DOMS in our research on 19 young male elite handball players who were continuously measured for HRV parameters during the employed exercise protocol. In order to establish DOMS-related differences in the HRV parameters of athletes, we represented characteristic HRV-parameters as a function of HR. While the HR-dependence of time- and frequency-domain parameters could be described by an exponential-like function, entropy (SampEn) showed a V-shaped function, with a minimum (“turning point”, HR_tp), separating descending and ascending intervals (before and after HR_tp, respectively).

In the HR-dependent representation, it was apparent that intense exercise loading upshifted the time- and frequency-domain parameters of all the athletes in the entire HR range, contrary to SampEn values that were systematically down-shifted, indicative of an upregulated sympathetic tone during the applied intensive exercise protocol. The standard deviation of the distribution of the group-level SampEn values were, in turn, found to be decreased, as well.

All the HRV-parameters of the high-DOMS group showed significant group-level differences compared to those of the low-DOMS group, which, however, change sign around HR_tp. The averaged SampEn(HR) function showed a non-linear character under the inducing exercise, with SampEn values lower for the high-DOMS than for the low-DOMS group below HR_tp, and vice versa above it. The HR-functions of the corresponding time-and frequency-domain parameters showed an opposite tendency, i.e., turning from positive to negative around HR_tp. Less prominent, but still significant, differences between the high- and low-DOMS groups could also be observed during the warm-up period, as well.

Time-domain parameters, such as RMSSD and SDNN, report on the size (amplitude) of HR-fluctuations, non-linear parameters, such as entropy or DFA-exponent, inform us about the stochastic structure of HRV [10]. Frequency-domain parameters, such as LF and HF power, on the other hand, form a link between the other two types of representation [10]. Whichever parameter one chooses, nevertheless, to characterize the momentary HRV, it will strongly depend on the HR value as it was demonstrated before [32,45,46,47,48]. There are a number of publications (e.g., [45,46,47,48]) discussing this problem, but lacking the exact mechanistic explanation that describes the relationship between HRV parameters and HR. Moreover, there has been no practical way to define a single, ”heart-rate-corrected” HRV parameter that could properly account for the HR-dependence. Moreover, HRV-entropy reflects the orderliness of the “funny” current modulation, however it is not fully understood what are the underlying physiological processes that it depends on, and its HR-dependence has not been clarified either. Factors that can influence HRV-entropy can come from the ANS function, cellular processes (e.g. mitochondrial function, etc.), or inflammatory factors during exercise, etc. [32]. As for the HR-dependence of SampEn, some people associate the value of the turning point with the intrinsic frequency [48], which is also supported by the age-dependence of the entropy, though, there is no unequivocal evidence for this [32]. A decrease in entropy at a given HR under exercise may indicate altered ANS function, and especially altered cellular processes. The current authors would like to suggest as mechanistic explanation for the above mentioned unexplained HR dependence that the ANS controls ultradian fluctuations through baroreflex sensitivity in a HR dependent manner [35] and Piezo2 (and Piezo1 as diurnal sensor) exerts fine regulation over HR [36], as ultradian sensor. Not to mention that DOMS may come with a transient proprioceptive neural switch initiated on Piezo2 as the primary damage, including alterations in neurocellular metabolic and energy generation [26].

In an earlier study signals were recorded for several-hours-long period under natural, daily-routine circumstances by a wearable ECG-device (i.e., Holter-monitoring, where othervise no special measuring conditions are required), and data were used to derive a reference curve, called ”Master Curve” (MC) [43]. The MC, with the help of non-linear method based on modified Poincaré plots, was able to describe the HR-dependence and appeared to be rather conservative for each individual from the daily to monthly scales [43]. It was suggested in this earlier paper that the momentary deviations detectable on the minutes scale from MC, could potentially be usable as a measure of changes in the dynamics of the heart rhythm, accompanying a disbalance of the ANS. Later on, this concept was applied to reveal stress- and relaxation-related physiological changes [44]. On the other hand, long-term dependence of the MC and other HRV parameters was described and interpreted as an age-related remodeling of the cardiac system [32], in line with the acquired Piezo2 channeloapthy induced quad-phasic non-contact injury model [9].

Earlier it was proposed that LF power is not a reflection of sympathetic tone, as widely viewed even today, but rather the modulation of cardiac autonomic outflows by the baroreflex [49]. Moreover, in line with this proposition it has been also suggested that LF power reflects mostly upon the excitatory activity level of Piezo2 in the baroreceptors and the heart [33]. Indeed, research is emerging to detect the role of Piezo2 in the heart as well [50,51]. Since a crosstalk is suspected between the proprioceptive system and the ANS through Piezo2 crosstalk [18,33], therefore the novel unaccounted ultrafast proton-based long-distance proprioceptive signaling in the nervous system [25] could be reflected in the cardiac autonomic outflows by the baroreflex in the form of proton release by Piezo2.

The LF power of HRV has been observed to be dampened with age as an indication of decreased ANS regulation on HR [49] and that is also in support of the aforementioned finding, leading to the age-related remodeling of the cardiac system [32]. Moreover, it is known that under medium to high intensity exercise loading the LF power and SampEn takes an inverse course when LF power decreases almost to 0 value then SampEn is the greatest. This point is suggested to be the Piezo2 inactivation moment in order to prevent pathological hyperexcitation to maintain homeostasis [52]. Interesting to note that Type Ia proprioceptive glutamatergic neurons also GABAergic [53] and GABA likely promotes this inactivation on Type Ia proprioceptors, leading to the hyperexcitation of ASIC3 containing Type II proprioceptive fibers [18]. Under physiologic conditions this increased sympathetic loading induced Piezo2 inactivation reflected in the near 0 value of LF power leads to a triggered automaticity in the SAN pacemaker cells likely regulated by Ca_v1.3 channels on a different pathway, and not by Piezo2 anymore [52]. Furthermore, the Piezo2 inactivation at higher sympathetic loading intensities, induced by GABAergic inhibition [53], likely reflects the desynchronization of the parasympathetic vagal and sympathetic regulation, hence the cardiac sympathetic activation is left unopposed [52]. Moreover, the gradual inactivation of Piezo2 under this proposed ASR may trade the fine motor movements for the “fight-or-flight” response [17,33]. In addition, the current authors suggest that this trade for the “fight-or-flight” response also sacrifices the fine modulation of blood pressure and HR, or the abovementioned ultradian sensory regulation, through the baroreceptors.

In support, a recent traumatic brain injury (TBI) study revealed the essential contribution of Piezo2 in the defensive arousal response (DAR) [54]. Important to note that mild TBI is suggested to be an analogous bi-phasic non-contact injury mechanism, like DOMS is, where the primary damage likely involves an acquired proprioceptive neuron terminal Piezo2 channelopathy, too [18]. This DAR mechanism is pivotal for survival and it is instigated by perceived threat and evoked by visual and auditory cues in the presence of motor abilities [54] and could be similar to the ASR, implicated in DOMS [17]. A recent unpublished preprint paper emphasizes the ultrafast matching of the Piezo2 initiated eye-brain, auditory/vestibular-brain and proprioceptive muscle-brain axes within the hippocampal hub [14]. In addition, the proprioceptive impairment, as a result of DOMS, is reflected in a tendency of mimicking a positive Romberg test [14]. Correspondingly, the aforementioned TBI study showed with genetic manipulation that reintroducing PIEZO2 promoted the reduction in escape latency and an increase in escape speed during DAR [54]. In fact, neural Piezo2 was shown to be the one activating the DAR in association with enhanced motor abilities [54]. The current authors suggest that this TBI study may be supportive of Piezo2’s ultradian sensory and ultradian rhythm generation function, and Piezo2 channelopathy is why, DOMS alters the response to postural perturbations [55], and significantly increases the medium latency response of the stretch reflex [56]. Moreover, the Piezo2 initiated heart-brain oscillatory axis has been also theorized [33,57]. The lack of the integration of the principal ultradian backbone of this oscillatory heart-brain axis into the hippocampal hub is indicative in the degenerative heart failure condition with hippocampal damage that is associated with short-term memory loss [58], reinforcing the role of hippocampus in learning and memory. Here, we should add that Piezo2-induced repeated firing, in turn, causes local membrane tension changes and cytoskeletal rearrangements (actin dynamics), not to mention the osmotic shifts and swelling/shrinking of the neuron cell body. These mechanical and structural changes are sensed again by mechanosensitive ion channels like Piezo2, in association with pressure pulse detection capability of Piezo2 [59], completing the Piezo2-mediated feedback loop. Such positive-feedback loops serve as bases of rhythm generation [60] and memory imprinting processes [61], as well.

It is worthy of considering that the aforementioned proton-release capability of Piezo2 is symmetry breaking, causing the collapse of the disordered symmetric state in order to accomplish an ordered, however not symmetric, state as acute intensive exercise loading increases. Moreover, protons could initiate the novel unaccounted ultrafast proton-based cross-frequency-coupled long-distance proprioceptive signaling in the nervous system between proprioceptive Piezo2 and hippocampal Piezo2 at a distance, constructing the proprioceptive backbone of the muscle brain-axis [34]. As a distant analogy, reentrant superconductivity may arise in superconductor and ferromagnetic layers at a higher (also lower) temperature phase, leading to pairing oscillate in the space where the superconducting state presents a higher order compare to the normal state [62]. Similarly, this is why Piezo2 could be coined as a principle cross-frequency-coupler, or entrainer, under stress. However, the microdamaging event leading to acquired Piezo2 channelopathy breaks this function as the primary damage, resulting in proton reversal [18]. We suspect that this primary damage point may be reflected in HRV measures because the ultradian rhythm is likely transduced to the heart by the cross-frequency coupling of the muscle-brain, heart-brain axes and the ANS through Piezo2-Piezo2 crosstalk in a HR dependent manner. Carrying on with solid-state-physics analogies, the Piezo2 on proprioceptive terminals was coined to have a low-frequency Schottky semiconductor barrier diode like function [25]. Indeed, super-Schottky diode exist within superconducting tunnel junctions, not to mention in superconductor-semiconductor junction [63]. Moreover, it is theorized that two-channel point-contact tunneling may exist in superconductors [64], giving further hints for the possible schemes of the novel unaccounted ultrafast proton-based long-distance proprioceptive signaling in the nervous system [26].

Thereupon, the question arises what SampEn suggests in reference to this long-distance proton-based signaling mechanism. SampEn is an approximation of the complexity of time-series dynamic physiological signals [42]. The principality of Piezo2 in proprioception [27] is suggested to come from its enigmatic feature that only Piezo2 could initiate the novel unaccounted ultrafast proton-based long-distance proprioceptive signaling within the nervous system [25]. The current authors view this novel ultrafast proprioceptive signaling could be reflected in the cardiac autonomic outflows by the baroreflex, or LF power, in the form of proton release by Piezo2. Hence, SampEn may represent the Piezo2 modulating capability of excess entropy, and Piezo2 channelopathy or proton reversal impairs this modulating capability. By accounting so, the current authors suggest that from low- to moderate load range, due to the increasing Piezo2 activity, the orderliness of events increases monotonically by the extent of load, implying that entropy (measured as SampEn in HRV terms), in turn, decreases up to the point where Piezo2 channels get inactivated upon the high strain of membrane due to increasing external load. From here on, the ordering action of Piezo2 channels cease when the load further increases, and entropy starts to increase again, due to other physiological side effects (e.g., local heating or increased reactive oxygen species (ROS) production). We tentatively assign this turning point of the V-shaped SampEn(HR) curve to the inactivation of the Piezo2 channels. Note that below HR_tp, LF-power and SampEn correlate, while above it, one finds an anti-correlation (see Figure 3 and Figure 5). Noteworthy also, that the involvement of heat shock protein 70 (Hsp70) activation through the Hsp70/TLR4/Interleukin-6/TNF-α pathway implicated in DOMS [65] and theorized on the route to Piezo2 channelopathy [16]. In addition, both the spontaneous activity dependent interactions and the potentially life-threatening activity contribution, initiated by DRA or ASR, of ultradian rhythmicity increases the temperature of brown adipose tissue by approximately 1^OC [39]. Moreover, this elevation of temperature is associated with approximately 0.8^OC increase of brain and 0.8^OC increase of body temperature [39], explaining the “heat of battle” response during DOMS-inducing exercise [17]. Accordingly, the low activity level of Piezo2, as the equivalent of high excitatory energy, keeping the SampEn low on a reverse course, reflecting a peak in excess entropy. However, under acute intensive exercise loading the excitatory activity level of Piezo2, as the equivalent of low excitatory energy, leaves SampEn on high course, reflecting lower excess entropy.

As for the connection between HRV parameters and DOMS, the first conclusion one may deduce is that below the turning point (HR_tp), all the time- and frequency-domain parameters are greater for the High-DOMS group during the whole training, implying an increased ANS activity at low HR (Figure 3 and Figure 4). Above HR_tp, however, an increasing level of autonomic attenuation develops as HR increases. This appears to correlate with the hypothesized increasing inactivation of Piezo2 above HR_tp. However, the entropy under intensive exercise load behaves the opposite way around: at low HR values, the SampEn values of the low-DOMS group are somewhat higher, but above HR_tp, the the SampEn level of the High-DOMS group start to dominate steeply, as HR rises. We suspect it to reflect an increased level of Piezo2 channelopathy among the athletes of the High-DOMS group.

The delayed onset of pain mystery and the delayed involvement of nociceptive C-fibers in DOMS is an intriguing phenomenon when it comes to the link between pain and HRV alterations, therefore the understanding of the underlying mechanism is important. Earlier it was suggested that this delayed pain sensation of DOMS and its associated delayed movement limitation in the form of reduced joint range of motion and loss of force-generating capacity, initiated by DAR/ASR, had provided 8 hours of pain-free limitless escape from danger in the wild evolutionarily [17]. The primary damage phase of the neurocentric DOMS mechanism theory may be initiated by a Piezo2 channelopathy on Type Ia proprioceptive terminals within the muscle spindle, leading to the secondary damage phase with harsher tissue damage in the extrafusal space with the involvement of the wide dynamic range (WDR) neurons on the spinal dorsal horn and other ion channels [18]. It is also important to note that acquired Piezo2 channelopathy is likely the initiating peripheral input source that drives central sensitization on spinal nociceptive neurons, even in the case of DOMS as well [18,66,67]. Thereupon, acquired Piezo2 channelopathy is proposed to be the autonomous pain generator [18]. Indeed, the evolvement of pain and sensitization are lost as a consequence of loss-of-function mutations on PIEZO2 [68]. Furthermore, the aforementioned new neurocentric DOMS hypothesis heavily relied on the so called gate control theory of pain [17], conceived by Melzack and Wall [69], and the activation of the aforementioned WDR neurons on the spinal dorsal horn [18]. Moreover, a research finding showed a gate control mechanism of pain on the dorsal root ganglion (DRG) by proprioceptive neurons [53]. Interestingly, this study also theorized the presence of quantum tunneling in reference to the presented pain gating mechanism [53]. However, the principal gate control of pain mechanism may reside further upstream on Type Ia proprioceptive terminals in the case of DOMS where the acquired channelopathy of Piezo2 may lead to neural switch and pain evolvement [9]. However, the pain evolves only as a result of the secondary damage phase of DOMS with nociceptive C-fiber contribution [18], when the aforementioned gate control mechanism of pain becomes critical downstream on the DRG, leading to WDR activation in the spinal cord [53]. This means two things in support of the bi-phasic neucentric DOMS theory. First, in the presence of intact intrafusal Type Ia proprioceptive terminal open gate cannot prevail for pain evolvement, hence the primary damage or Piezo2 channelopathy is a must. Second, extrafusal secondary damage is also essential for the activation of C-fibers in order for pain evolvement in DOMS [23].

An interesting translational pilot study investigated the sympathetic regulation of the DOMS effect [70]. One of the key findings was the evolved neurogenic inflammation with the involvement of the ANS, particularly the sympathetic nervous system (SNS) [70]. Correspondingly, the more the SNS was activated, the more the pain and inflammation was augmented [70]. Even more importantly, the applied stellate ganglion block (SGB) was capable of getting in the way of full DOMS evolvement [70]. The current authors translate these results, based on a recently published paper [18], that the acquired Piezo2 channelopathy on Type Ia proprioceptive terminal instigated the so-called inflammatory and even the gateway reflex during DOMS-inducing exercise. However, with the applied SGB the pathway to the secondary damage phase was completely interrupted, equivalent of irreversible Piezo2 channelopathy with no pain and sensitization [18,52], hence the so-called gateway reflex could not fully evolve. The data of our current research is in support of increased SNS activity as a result of DOMS-inducing exercise bout and the correlation between SNS activity and pain (Figure 3., Figure S1). We suggest that the increased SNS activity is likely driven by unopposed Ca_v1.3 ion channel due to the impaired Piezo2-Piezo2 crosstalk between the proprioceptive muscle-brain axis, the heart-brain axis and the ANS.

We touched upon the relevance of DAR and ASR as important underlying stress mechanisms. Essentially, personal differences of the stress response may exist, especially between trained and untrained individuals [71,72]. Since we followed a neurocentric tracing in regards to DOMS, therefore allostatic stress of neurons should be considered as well under prolonged sympathetic loading of eccentric (forced lengthening) contractions [18]. More specifically, two states of Piezo2 may prevail under allostatic stress, namely inactivated intact Piezo2 and the acquired channelopathy of Piezo2 due to the prolonged eccentric nature of DOMS-inducing exercise [18]. The fundamentals of this two distinct states under allostatic stress was laid down by János (Hans) Selye by naming good stress as eustress and bad stress as distress [73]. Correspondingly, underlying Piezo2 channelopathy under allostatic stress is analogous to Selye’s bad stress, or the gateway to pathophysiology, leading to impaired Piezo cross-talk and the inducement of the so-called gateway reflex [18]. On the contrary, the inactivated Piezo2 under allostatic stress is analogous to Selye’s eusrtress that induces the so-called inflammatory reflex and leads to remodeling and adaptation within homeostasis [18]. The adaptive mechanism is served by ASIC3 ion channels when it comes to proprioception, as secondary proprioceptive ion channels [9,74]. Noteworthy, that research showed this, secondary [74], protective role of ASIC3 in DOMS [75]. If we consider Piezo2 as a principle cross-frequency-coupler under stress [18] then good stress could be translated as frequency-coupled, or fine-tuned stress buffering (intact underlying Piezo function and crosstalk, modulation primarily taken over by secondary ASIC3), in contrast to bad stress that is decoupled or disbalanced (impaired underlying Piezo2-Piezo1 and Piezo2-Piezo2 crosstalk even when the modulation is taken over by secondary ASIC3). Since Piezo ion channels are evolutionarily conserved, therefore it is important to note that Piezo is the one buffering the mechanical stress via modulation of intracellular calcium handling in Drosophila heart, while the functional mutation of PIEZO fails to buffer mechanical stress, leading to pathological remodeling [76].

This demarcation line between good and bad stress is likely reflected in oxidative stress as well. Even the new neurocentric DOMS mechanism theory implicated the mitochondrial electron transport chain generated free radical production as an important underlying factor in the evolvement of the primary damage [17] that was later coined as Piezo2 channelopathy [23]. In fact, lost function of Piezo2 impairs nitric oxide synthases and instigates remodeling [77] and this is suggested to be analogous to the acquired Piezo2 channelopathy, or the primary damage, of DOMS. Another important consideration that ROS production in mitochondria could induce high-frequency oscillations [78]. Thereby this, ROS-dependent mitochondrial oscillatory signaling transduction pathway should be accounted in cardiomyocytes and myocytes under oxidative stress, like in case as a consequence of DOMS. Even more importantly the impaired low-frequency Schottky semiconductor barrier diode like function of Piezo2, as a result of Piezo2 channelopthy, may fail to modulate these ROS-dependent mitochondrial high frequency oscillations. In support it is evident that DOMS increases ROS production [79]. Moreover, the dual role of ROS in hippocampal learning and memory, more specifically in long-term potentiation or the aforementioned memory imprinting, and hippocampal neurotoxicity or even neurodegeneration has been long known as well [80,81]. Not to mention this dual role of ROS in the mitochondria of the heart is also present [82]. Noteworthy that the neurocentric DOMS theory emphasized that the primary damage could evolve at nerve terminals, like the Type Ia proprioceptive one, where the mitochondria content is high [17], leading to microscopically undetectable functional microdamage of Piezo2 [26]. Critical pathway of this hypothetical microdamage is the electron leakage serving ROS production essentially through the electron transport chain [17,25] and proton motive force [26], like in the case of the heart [82]. Since Piezo2 channelopathy is theorized to entail a proton affinity switch [26], therefore not only could explain the dual role of ROS, but the loci of Piezo2 in neuron terminals packed with mitochondria indicative of the abovementioned ultrafast axes. Acquired Piezo2 channelopathy was coined as the principal gateway to pathophysiology or the primary damage [9], as the common root cause of aging initiation [83], therefore increased ROS production may contributes to the inducement of the inflammatory reflex within homeostasis in support of remodeling, but the proton affinity switch may induce the gateway reflex as a breach of remodeling [18]. Correspondingly, Piezo2 channelopathy may not only the result of proton affinity switch and resultant neural switch or miswiring, but may also disrupt quantum tunneling of protons and electrons from mitochondria [26], leading to increased ROS production, and resultant increased entropy and aging acceleration. This later mentioned pathophysiology is analogous to the one observed in aging of brain mitochondria [84]. In addition, proton affinity switch may be the reason why reverse electron transport may also prevail. Indeed, reverse electron transport increases mitochondrial ROS production [85].

Compelling finding shows that paired associative transcranial (TES) and peripheral electromagnetic stimulation (PES) diminish DOMS and facilitate the reinstitution of force generation [86]. In the contrary, PET alone without paired associative TES was not capable to provide any DOMS related remedy [87]. Hence, the combined TES and PES not only in support of the neurocentric DOMS theory, but also implies that the muscle-brain axis is also critical in DOMS mechanism [74]. It has been presumed that sensing electromagnetic field induced oscillating energy is conveyed by a common receptor and in return this receptor could activate frequency dependent biological pathways with the involvement of interfacial water [88]. Piezo2 was theorized to be this common receptor and the one to initiate the ultrafast signalling in the nervous system [9]. In support of this theory, a recent study indeed shows that Piezo2 is the underlying mediator of magnetic stimulation [30]. This proton-based ultrafast signalling was also proposed to have a role in the force modulation of motoneurons [26]. Overall the paired associative TES and PES seems to reduce the secondary damage effect of DOMS, since the beneficial effect was detected in 72 hours after DOMS inducement [86], in which the impairment of the Piezo2-Piezo2 crosstalk is implicated not only along the proprioceptive terminal-hippocampal axis [25], but between the intrafusal and extrafusal compartment as well [23]. As for the concrete biophysical mechanism of this non-synaptic, bi-compartmental Piezo2 cross-communication, we suggest the rapid proton-conducting Grotthuss-mechanism to be a proper candidate [89]. It practically refers to a concerted action of a large number of quasi-simultaneous proton hopping events along a H-bonded network of water molecules. The muscle spindles, with proprioceptive terminal Piezo2 content, are enclosed with a connective tissue capsule, which contains hygroscopic glycosaminoglycan (GAG) residues attached to collagen bundles [90]. These GAG residues, such as chondroitin sulfate, heparan sulfate or hyaluronic acid side chains, are all of kosmotropic nature according to the Hofmeister-terminology [91], forming a strongly H-bonded network of water molecules, giving rise to conduction of protons released by the Piezo2 channels. On the other hand, the collagen bundles are themselves piezoelectric [92], i.e., they possess a large electric dipole moment along their axis upon external load, providing the asymmetry that will be able to define the directionality of the proton conduction. In support it is noteworthy that muscle spindles should be considered as a continuum with extrafusal space, and not as an entirely isolated structure [23]. Indeed, intrafusal muscle fibers for example reach beyond the muscle spindle capsule into extrafusal territory and tethered to the extrafusal connective tissue [93]. Since, extracellular matrix damage and impairment of the selective barrier of the muscle spindle capsule are likely part of the secondary damage phase of DOMS [18,23], therefore combined TES and PES treatment seems to reinvigorate Piezo2 from its proton reversal induced channelopathy by giving rise to proper proton conduction not only along the muscle-brain axis, but along the damaged bi-compartmental communication structures of Piezo2 ion channels. Finally, but not surprisingly, paired associative TES and PES even improves HRV parameters substantially [28], likely restituting the aforementioned Piezo2-Piezo2 crosstalk with the ANS.

5. Limitations

The aim of the study was to use the CON group as a control group initially. However, the 15 repetitions in 6 sets of maximum intensity CON protocol induced DOMS as well. Unfortunately, it turned out that 6 sets were too excessive for the CON group and the outcome is in agreement with earlier observation that at high exercise intensities there is no difference in the change in muscle soreness from pre-exercise to post-exercise intervals [94]. Since, the sample size of the ECC would have been too low, therefore the samples of the CON group was also included in the study.

Moreover, the research group must have collaborated with the coaches of athletes and their athletic program, therefore some of the athletes were exposed to the study exercise protocol after days of training, while some of the athletes were exposed right after summer recess. This difference may have impacted the VAS score levels, due to repeated bout effect of DOMS with diminished symptoms [95], however it did not impact the mere fact of DOMS inducement.

6. Conclusions

DOMS is a perplexing pain condition, still DOMS-specific alteration in HRV has not been detected in the scientific literature so far during DOMS-inducing exercise activity, when pain is not present yet, in contrast to after exercise activity when DOMS fully evolves. We suggest, and then prove, that HRV parameters considered as a function of HR can be used as sensitive measures of DOMS-related phenomena, already before the full development of pain. Consequently, the findings of the current study provide further indirect evidence for an initiating neural microdamage, also called primary damage, that prevails already under DOMS-inducing exercise and the diagnostic detection of this point empowers control over avoiding further injury risk in sport and exercise activities.