Comparison of fetal heart rate variability by symbolic dynamics at the third trimester of pregnancy and low-risk parturition.

Fetal heart rate variability (fHRV) is an essential source of information to monitor fetal well-being during pregnancy. This study aimed to apply a nonlinear approach, known as symbolic dynamics (SD), for comparing human fHRV in the third trimester of pregnancy during active fetal state (TT) and active labor at term (P). We performed a longitudinal, prospective, descriptive, and comparative study composed of 42 longitudinal recordings of 5-minutes of fetal heartbeat interval series. Recordings were collected from 21 low-risk, healthy, pregnant women attending the Maternal and Child Research Center (CIMIGen), Mexico City. We calculated relevant linear parameters of fHRV between TT and P stages, such as the percentage of differences between adjacent RR intervals >5 ms (PRR5, related to vagal modulations) and other SD parameters such as the percentage of no variations between three successive symbols (%0V, reflects sympathetic modulations) and the probability of low variability with a threshold of 4 ms (POLVAR4, associated with a low variability). We identified statistical differences for PRR5 between TT and P (37.13% [28.47-47.60%] vs. 28.84% [19.36-36.76%], p = 0.03), respectively. Also, for 0V% (65.66% [59.01-71.80%] vs. 71.14% [65.94-75.87%], p = 0.03) and for POLVAR4 values (0.06 [0.04-0.11] vs. 0.15 [0.09-0.24], p = 0.002), respectively. Our results indicate that during parturition, the short-term fetal fHRV is decreased, showing a decreased vagal modulations and higher adrenergic response of the heart. These autonomic modifications may result from the fetal response to the stressful inflammatory challenge of labor. We thus confirmed that the analysis of the SD applied to fHRV time series could be a potential clinical biomarker to differentiate the fetal autonomic cardiac condition at different stages of pregnancy.

Introduction

Despite that parturition has a short duration in comparison to the span of pregnancy, this period is of high risk for the fetus: 8–13% of critical neonatal neurological troubles occur during this period [1]. Uterine electrical activity and fetal heart rate monitoring, using transabdominal recordings, is widely used to identify fetal distress symptoms [2]. Heart rate fluctuations patterns reflect the autonomic nervous system (ANS) activity. Specifically, the analysis of the fetal heart rate variability (fHRV) has provided valuable information regarding fetal neurodevelopment [3]. The fHRV is regarded as a non-invasive biomarker of the developing cardiac parasympathetic and sympathetic activities [4] because the increment of fHRV along gestation is associated with fetal growth and neural integration [5].

Although fHRV offers crucial parameters to monitor fetal well-being, it is known that the classical linear analysis methods are limited for assessing all complex systems interacting to control heart rate [6]. Thus, nonlinear measures seem more consistent to reflect the prenatal developmental factors influencing cardiovascular regulation [7]. A nonlinear analysis, known as symbolic dynamics (SD), has also been employed to quantify fetal heart rate fluctuations during pregnancy, demonstrating that the use of this technique appears to provide a better and more differentiated understanding of healthy-physiological fetal development [8, 9, 10, 11]. Some studies have confirmed that SD can be a helpful analysis to classify fHRV [6] and to evaluate the individual fetal gestational progress during pregnancy [9]. Additionally, other evidence indicates that, as pregnancy progresses, fetal beat-to-beat heart rate fluctuations show evidence, not only of an increase in complexity, but also of increased regularity [12]. However, limited evidence exists about the fHRV changes occurring at low-risk term parturition, despite the fact that parturition involves a massive inflammatory and stressful processes for the fetus [13, 14].

Some issues still deserve to be evaluated for SD to be considered as a useful tool in the detection of fetal wellbeing during labor. First, it is still unknown whether their different indices convey independent information or contain, by contrast, redundant information despite the formal differences in the corresponding mathematical construction. Second, the analytical value of symbolic parameters in fetuses at term during low-risk labor is still largely unexplored, only few studies with a limited number of tested symbolic parameters have been carried out [8, 9, 10, 11]. Third, the relationship between SD and classic linear indices of fHRV has not yet been explored. Given that the percentage of differences between adjacent RR intervals larger than 5 ms (or PRR5) is a linear parameter that reflects universal fetal developmental characteristics based on heart rate patterns [15], exploring its association with novel symbolic parameters would be relevant in understanding the physiological significance of SD at labor and the third trimester of pregnancy. In this exploratory study, we addressed these points, providing the basis for a critical approach to choose the most informative symbolic index applied to fHRV.

Therefore, we applied SD to compare fHRV in healthy fetuses during active states of the third trimester and the active phase of labor of low-risk gestations. We explored if the SD of fHRV reflect the fetal autonomic modifications introduced by different factors such as the acute inflammatory and stressful condition of labor.

Materials and methods

Data extraction

We performed a longitudinal, prospective, descriptive, and comparative study with a non-probabilistic sampling. In this research, we used previous data from our database of bioelectrical signals (n = 67) recorded from the maternal abdominal wall at the third trimester of pregnancy (TT) and during the active parturition at term (P) to perform a new analysis based on SD of fHRV [16, 17]. Data were collected from 46 low-risk, healthy, pregnant women attending the Maternal and Childhood Research Center (CIMIGen), Mexico City, Mexico. The Ethics Commission from the Biological and Health Sciences Division (CBS) at Metropolitan Autonomous University (UAM), Campus Iztapalapa (ref. CAEDCBS.01.2017), approved this research protocol. All the data files were obtained from participants under informed consent.

The database includes fetal RR and maternal anthropometric data (maternal weight, waist, and hip circumferences) from participants studied during pregnancy. Recordings satisfied the following inclusion criteria: young women aged 18–32 years old, residents of Mexico City or its metropolitan area, at the third trimester (TT) of gestation (from 35 to 38 weeks), maternal BMI lower than 35 kg/m2 and not presenting any clinical manifestation of the initiation of labor. Our database also included data of some these pregnant women at the active phase of term parturition (P) stage (>39 weeks), showing the presence of three to four uterine contractions in 10 min evaluated by electrical uterine activity, at least 4 cm of cervical dilatation, and 50% of cervical effacement (both assessed by cervical exams). None of them received the administration of epidural anesthesia and presented any clinical manifestation of chorioamnionitis.

The exclusion criteria for both stages involved women with multiple gestations, hypertension disorders in pregnancy, severe allergies, diabetes mellitus, autoimmune disease, renal dysfunction, maternal or fetal infection, and alcohol or drug consumption during pregnancy.

Database acquisition setup

Briefly, abdominal recordings from our database were acquired using a portable maternal-fetal monitor (Monica AN24®, Monica Healthcare, Nottingham, UK) [18]. The transabdominal ECG composite signals were recorded for 5 min using disposable electrodes (Ambu® BlueSensor VL) in a bipolar configuration while women maintained a semi-Fowler's position either at TT or P (Figure 1). The electrodes were positioned after cleaning the abdominal area with an alcohol swab and after carefully abrading the skin with sandpaper tape to reduce skin impedance.

Figure 1

Transabdominal electrodes attachment. Electrode (1) is the ground reference, it was positioned towards the back, while electrodes (3), (4), and (5) were positioned on the maternal abdomen approximating an arc that resembles the arc of the participant's uterus fundus. Electrode (2) is placed at a location to approximate the symphysis pubis of the subject. Electrode (2) be placed between 2 cm and 5 cm above (rostral) of the symphysis pubis. The figure is modified from Escalante-Gaytán et al., 2018 [56].

Fetal ECG recordings of both pregnancy stages were extracted and visually inspected (Figure 2) by using the Monica DK software (Monica Healthcare, Nottingham, UK), which also displays values of fetal heart rate and uterine activity. We applied, as the selection criterion of fHRV time series at TT, the identification of fetuses in active state by a visual inspection showing a good variability (at least equal to 5 bpm) and at least one fetal movement reflected in the manifestation of a heart rate acceleration [19]. We also considered the manifestation of regular uterine activity (four uterine contractions in 10 min) in correspondence with fHRV time series extraction for P (Figure 2). For both stages of pregnancy, we only selected those fHRV beat-to-beat time series having a duration of no less than five continuous minutes according to the standards of short-term analysis of HRV [20]. We discarded incomplete longitudinal data sets with only one successful time series recovered at either the TT or P stages and including misdetections of consecutive five minutes of fetal RR intervals. The fHRV time series were reconditioned by a filtering approach to exclude ectopic beats and artifacts [21].

Figure 2

Representative fetal heart rate variability (fHRV) traces and uterine activity expressed in arbitrary units during an active stage of the third trimester of pregnancy (a) and active parturition (b) from the same participant (37 and 39 gestational weeks, respectively). Vertical marks provided by an accelerometer are depicted to indicate the presence of gross maternal movements.

fHRV linear analysis

The following time-domain indices were calculated for the third trimester of pregnancy and parturition using the Kubios HRV software (University of Eastern Finland, Kuopio, Finland): the mean value of RR intervals () and heart rate (), the root mean square of successive differences (RMSSD) and the standard deviation of RR intervals (SDRR). Complementary linear indices were also included such as the corrected SDRR (cSDRR) [22] and the percentage of differences between adjacent RR intervals >5 ms (PRR5) related to the manifestation of fetal vagal modulations [11].

Symbolic dynamics: construction of symbolic sequences

Previous studies have shown that SD is an efficient approach to analyze the dynamics of HRV data [23, 24]. Some detailed data is lost in this process, but the dynamic behavior can still be analyzed [25]. SD analysis of fHRV has shown to be a reliable marker of the fetal autonomic condition [8, 9, 10, 11]. The different transformations to symbolize the fHRV time series produce symbolic series that contain specific dynamical data. However, it remains to be clarified whether the dynamical data retained by each transformation is different [9].

The first step of SD is the construction of the symbolic sequences; the fHRV time series were transformed into a symbolic series using two different approaches:

First symbolization: σ-method

The fHRV time series x, x, x,…,x is transformed into the symbolic sequence s,s,s,…,s, s ε A (Eq. (1)) on the basis of the alphabet A = {0,1,2,3}. The transformation to symbols refers to four given levels where μ denotes the mean RR interval.

The parameter “a” describes the threshold above and below the mean concerning the mean. In this study, we set a = 0.05, i.e., the thresholds were 5% above or below the mean [26].

Second symbolization: binary Δ-coding-method

We applied a binary coding method using the symbols ‘‘0″ and ‘‘1’’ (Eq. (2)) to indicate differences of RR intervals below or above a threshold [26], i.e.:where Δms indicates the time interval. In this study, the thresholds were set to Δms = 3, 4, and 5 ms. These thresholds are lower in comparison with studies dealing with RR data from adults because the mean RR interval of the fetus is also considerably lower than the mean RR interval of adults [27].

Analysis of symbolic dynamics

Quantification of symbolic series σ

In this study, we analyze the frequency distribution of length three words, i.e., substrings that consist of three symbols from the alphabet A={0,1,2,3}. The following SD indices were determined using the probability distribution of the 43 = 64 possible word types [24, 28, 29] from the first symbolization or σ-method:

WSDVAR: standard deviation of the word sequence.

FORBWORD: forbidden words, number of seldom (p < 0.001), or never occurring word types. Lower FORBWORD values are associated with higher complexity.

WPSUM02: the relative portion of words consisting only of the symbols ‘0’ and ‘2’; a measure of a decreased HRV.

WPSUM13: the relative portion of words consisting only of the symbols ‘1’ and ‘3’; a measure of an increased HRV.

FWSHANNON: Shannon entropy of the word distribution.

FWRENYI025: Renyi entropy of the word distribution with weight coefficient α = 0.25 to quantify smaller probabilities.

FWRENYI4: Renyi entropy of the word distribution with weight coefficient α = 4 to quantify larger probabilities.

The sequences of length k = 3 were categorized according to their amount of variations between successive symbols as follows [9]:

0V%: Percentage of no variations between three successive symbols, i.e., all three symbols are equal (considered to reflect the sympathetic modulation only).

1V%: Percentage of one variation between three successive symbols, i.e., two symbols are equal (considered to reflect the sympathetic and parasympathetic modulation).

2LV%: Percentage of two like variations between consecutive symbols, i.e., two successive decrements or increments of symbols (e.g., ‘123’ or ‘321’); reflects sympathetic and parasympathetic modulation with a vagal predominance.

2UV%: Percentage of two unlike variations between successive symbols, i.e., one decrement followed by an increment or vice versa (e.g., ‘131’ or ‘213’); reflects, exclusively, the vagal modulation.

Quantification of binary symbolic series

We created sequences of length k = 6, i.e., the 26 = 64 different binary patterns. A measure based on this kind of symbolic sequences is the probability of a low variability (POLVARΔ) that it is equal to the probability of occurrence of the subsequence of six zeros ‘‘000000″ in the symbolic string [26] and it reflects decreased HRV. In our comparison of features, we explored different , i.e., POLVAR3, POLVAR4, and POLVAR5, respectively. For instance, POLVAR4 is the probability of occurrence for subsequence “000000” with equal to 4ms in Eq. (2).

Fractal methods

The detrended fluctuation analysis (DFA) provided the scaling exponent α1 as detailed elsewhere [30]; this method has proven useful in revealing the occurrence of long-range correlations in time series. The first step is the integration of the RR time series. Subsequently, the integrated series are divided into segments having an equal number of n segments. The local trends Yn are obtained for all segments by a least-squared linear fit and subtracted from Y(k) to reduce the non-stationary artifacts. The average root-mean-square fluctuations, F(n), are then calculated:

The relationships on a double-log graph between F(n) and time scales n are approximated by a linear model F(n) ~ nα, so providing the scaling exponent α1 as the slope of the plot covering the short-term range of n from 4 to 11 intervals.

Additionally, we applied a magnitude and sign analysis (MSA) [31]. The original RR sequences are processed to obtain series of increments by taking the differences between adjacent intervals (RRi+1−RRi). These series (ΔRR) are decomposed into magnitude |ΔRR| and sign series ±(ΔRR). After subtracting their respective means, the magnitude and sign series are integrated, and DFA is again applied as described above. The slope of F(n)/n covering the range from 4 to 11 intervals then provides magnitude and sign scaling exponents (α1(MAG) and α1(SIGN), respectively). The α1(SIGN) exponent includes information about the temporal directionality of the original series or the way series’ increments (decrements) alternate, indicating if a further positive or negative increment (decrement) is more likely to occur given a current increment (decrement). Results suggest that the interaction between the sympathetic and the parasympathetic systems is reflected by the sign of the heartbeat increments [32].

All these calculations were obtained using Matlab® software (the MathWorks, Inc. Natick, Massachusetts, USA).

Statistical analysis

The statistical analysis was carried out using the GraphPad Prism version 7.00 for Windows, (GraphPad Software, La Jolla California USA). Continuous variables are expressed as median (25th–75th percentile), and categorical variables were expressed as numbers and percentages. Normal distributions were tested by Shapiro-Wilk tests. Paired t-test or Wilcoxon signed-rank were used to compare fHRV parameters for TT and P stages. Clinical characteristics of both stages were compared using a Wilcoxon signed-rank test. Finally, Spearman's nonparametric correlation coefficient (rsp) and linear regression were estimated for analyzing the associations between the differences of PRR5 and the differences of SD parameters at both stages. Significance was considered by p < 0.05 for all analyses.

Results

A total of 67 recordings at the TT and P stages were inspected from our database. However, after applying the above selection criteria, we analyzed data collected from 21 participants that included both stages (i.e., 42 beat-to-beat fHRV time series). The median gestational age was 36.5 [35.3, 37.7] weeks at the third trimester and 39.6 [39.0, 40.3] weeks during active parturition; no complications occurred in newborns as confirmed by median birth weight 3.3 [3.0, 3.6] kg and median Apgar scores at 1 min of 8.0 [8.0, 8.1] and at 5 min of 9.0, [9.0, 9.1] points. Given the absence of clinical symptoms such as maternal tachycardia, fetal tachycardia, and purulent amniotic fluid, no clinical signs of chorioamnionitis were noted at the time of parturition. No malformations were observed in any of the newborns. Other fetal and maternal clinical characteristics are shown in Table 1.

Table 1

Maternal and newborn clinical characteristics of the study cases (n = 21).

Description
Maternal age (y)	25.0, [19.5, 27.5]
Maternal heart rate-third trimester (bpm)Maternal heart rate – parturition (bpm)	89, [80.0,95.0]93, [85.0,99.0]
Gestational age - third trimester (weeks)	∗36.5, [35.3, 37.7]
Gestational age – parturition (weeks)	39.6, [39.0, 40.3]
Uterine activity at labor (contractions in 10 min)	4.0, [4.0, 5.0]
Cervical dilatation at labor (cm)	4.0, [4.0, 7.0]
Cervical effacement at labor (%)	70.0, [60, 80]
Maternal weight (kg)	66.5, [63.2, 78.4]
BMI (kg/m2) at the third trimester	27.0, [25.2, 31.1]
Waist circumference at the third trimester (cm)	102.0, [93.5, 109.5]
Hip circumference at the third trimester (cm)	107.0, [96.5, 113.5]
Caesarean section (n, %)	7 (33)
Fetal gestational age-parturition (weeks)	39.6, [39.0, 40.3]
Apgar
1 Minute	8.0, [8.0, 8.1]
5 Minutes	9.0, [9.0, 9.1]
Birth Weight (kg)	3.3, [3.0, 3.6]
Head Circumference (cm)	33.5, [32.7, 34.5]
Length (cm)	50.0, [49.0, 51.0]
Sex (n, % male)Cases of chorioamnionitis	10 (47)0/21

Values expressed as median (interquartile range) unless otherwise indicated.

∗p < 0.0001 between the third trimester of pregnancy and parturition (Wilcoxon signed-rank test).

In this study, a total of 24 indices were obtained by applying both linear and SD analysis on filtered fHRV time series. The estimated fHRV features are presented in Table 2 using the median and interquartile range. By performing the statistical analysis, we found that the PRR5 was the only linear parameter that exhibited a significant decrement at P in comparison with TT (p = 0.03, Table 2). In contrast, SD parameters from the first symbolization presented significant differences between both stages with low p-values (p < 0.03, Table 2). For instance, WPSUM02 revealed decreased fHRV in P in comparison with TT (p = 0.02, Table 2). Other parameters from the σ-method, such as 1V% and 0V%, also showed significant modifications introduced by parturition (p<0.03, Table 2).

Table 2

Comparison of longitudinal fetal heart rate variability (fHRV) parameters between the third trimester of pregnancy and parturition (N = 21).

Parameter	Third trimester (TT)	Parturition (P)	Significance (p-value)
RR¯ (ms)	425.70 [406.50–448.15]	431.70 [411.60–452.00]	0.80
HR¯ (bpm)	140.95 [133.88–147.60]	139.00 [132.75–145.78]	0.85
RMSSD (ms)	8.70 [6.32–10.96]	6.52 [5.50–9.03]	0.11
SDRR (ms)	12.30 [9.05–13.18]	8.87 [7.15–12.48]	0.13
cSDRR (ms)	123.68 [107.18–147.85]	98.65 [75.24–134.49]	0.17
PRR5 (%)	∗37.13 [28.47–47.60]	28.84 [19.36–36.76]	0.03
WSDVAR	∗1.51 [1.36–1.96]	1.09 [0.86–1.85]	0.03
FORBWORD	43.00 [41.0–44.00]	43.00 [42.00–46.00]	0.09
WPSUM02	∗0.54 [0.37–0.63]	0.77 [0.46–0.84]	0.02
WPSUM13	0.20 [0.15–0.35]	0.09 [0.04–0.31]	0.06
0V%	+65.66 [59.01–71.80]	71.14 [65.94–75.87]	0.02
1V%	+27.14 [22.00–32.41]	21.58 [19.22–24.03]	0.02
2LV%	0.43 [0.00–0.63]	0.00 [0.00–0.43]	0.08
2UV%	6.84 [5.20–9.28]	6.89 [5.09–9.07]	0.12
FWSHANNON	2.37 [2.17–2.54]	2.17 [1.84–2.39]	0.07
FWRENYI025	3.12 [2.96–3.22]	3.02 [2.86–3.19]	0.25
FWRENYI4	1.56 [1.39–1.70]	1.39 [1.11–1.73]	0.10
POLVAR3	∗0.03 [0.01–0.06]	0.06 [0.04–0.10]	0.01
POLVAR4	∗0.06 [0.04–0.11]	0.15 [0.09–0.24]	0.002
POLVAR5	∗0.13 [0.08–0.21]	0.26 [0.14–0.38]	0.006
α_1	1.27 [1.03–1.39]	1.25 [1.10–1.39]	0.61
α_2	1.28 [1.08–1.37]	1.16 [1.08–1.31]	0.46
α_1(MAG)	0.75 [0.63–0.81]	0.68 [0.65–0.80]	0.82
α_1(SIGN)	0.48 [0.36–0.55]	0.46 [0.39–0.54]	0.90

Note: values are expressed as median (interquartile range).

between the third trimester of pregnancy and parturition (Paired T-test).

between the third trimester of pregnancy and parturition (Wilcoxon signed-rank test).

Furthermore, all POLVARΔ measures calculated from the second symbolization or binary Δ-coding-method confirmed lower fHRV during P in comparison with TT of pregnancy (Table 2). Among all POLVARΔ parameters studied, the POLVAR4 presented the lowest p-value (p = 0.002, Table 2).

The binary series S(x) reflect the acceleration and deceleration of fHRV. For instance, Figure 3a exhibits a representative example of the absolute values of differences of fetal RR adjacent intervals (|RRi+1−RRi|) or (|ΔRR|) at the third trimester of pregnancy, whereas Figure 3b corresponds to data from the same participant at active parturition (for visualization only the first 50 beats are shown). The dashed line indicates a threshold Δ = 4 ms, related to the probability of low variability (POLVAR4). The lower diagrams show the corresponding binary series; we can observe a higher number of consecutive “000000” patterns in the binary series of Figure 3b in comparison to Figure 3a, indicating lower fetal variability at parturition.

Figure 3

Examples of the binary Δ-coding-method transformation: (a) Binary transformation using absolute values of the differences of fetal RR adjacent intervals (|RRi+1−RRi|) or (|ΔRR|) derived from a fetus at the third trimester of pregnancy; b) |ΔRR| from the same case at active parturition at term. The dashed line indicates a threshold Δ = 4 ms, related to the probability of low variability (POLVAR4). The lower diagrams show the corresponding binary series.

We found significant negative correlations of the differences of PRR5 (ΔPRR5) and the differences of relevant symbolic dynamics parameters between TT and P (Table 3). Specifically, the ΔPRR5 was positively correlated with Δ1V% (rsp = 0.77, p < 0.0001) and negatively correlated with Δ0V% (rsp = -0.83, p < 0.0001). Interestingly, ΔPRR5 was negatively correlated with all POLVARΔ measures (p < 0.0001).

Table 3

Spearman's correlation and linear regression analysis between the differences of the parameter PRR5 (ΔPRR5) and differences of relevant symbolic dynamics parameters between the third trimester of pregnancy and parturition stages (N = 21).

Parameter	Spearman rsp	p-value	R2	F value
ΔWSDVAR	0.24	0.2862	0.09	1.96
ΔWPSUM02	-0.28	0.2115	0.11	2.33
Δ0V%	-0.83	<0.0001	0.61	30.87
Δ1V%	0.77	<0.0001	0.56	24.85
ΔPOLVAR3	-0.78	<0.0001	0.58	27.29
ΔPOLVAR4	-0.73	<0.0001	0.64	34.81
ΔPOLVAR5	-0.90	<0.0001	0.86	123.5

Discussion

This longitudinally study compared linear and SD indices of fHRV between the third trimester of pregnancy (TT) and active parturition at term (P). The results indicate that the fHRV in active labor shows a significantly lower variability. Besides, our results from SD first symbolization show that the percentage of 0V% was significantly higher at parturition as compared with the third trimester of pregnancy, suggesting that the fetal adrenergic response of the heart is prominent during labor at term, because an increased sympathetic activity is associated to higher values of 0V% [28, 33]. This activity was consistently accompanied by diminished effects of the cardiac parasympathetic modulation at parturition, as indicated by lower PRR5 values as well as decreased fHRV indicated by the WSDVAR and WPSUM02 values. Interestingly, there were no significant changes in the between both stages despite the ongoing uterine activity at active labor.

These results were confirmed by the SD second symbolization using the POLVARΔ measures because a decreased fHRV at parturition was found in comparison with the third trimester of pregnancy. Hence, POLVARΔ indices, related to the binary Δ-coding-method, could provide more detailed information about cardiac autonomic nervous control during parturition in comparison to the linear and the σ-method. Other studies have already established that binary SD representations of fHRV reflect individual gestational development [9, 10]. Noteworthy, changes at P in comparison with TT were not possible to be identified in none of the scaling exponents (α1, α2, α1(MAG) or α1(SIGN)), which have been used in previous studies to assess maternal HRV during pregnancy [34].

Notwithstanding that for healthy fetuses, various linear and SD fHRV measures increase as pregnancy progresses [9, 35] and that we found significant differences between the stages of the third trimester and active labor (p < 0.0001, Table 1). We cannot attribute the decreased fHRV or increased 1V% values to differences in the gestational age. Instead, among several factors involved, our results seem related to the fact that uterine contractions represent by themselves stressful conditions for the fetus during active parturition [36]. Yet, as this study analyzed only healthy outcomes and uncompromising fetuses, the decrease in the SD parameters of fHRV during labor are considered to represent the manifestation of such labor's stressful environment.

Other prenatal studies support the use of beat-to-beat fHRV analysis as a non-invasive, continuous, sensitive, and specific approach to characterize the fetal inflammatory response [37]. This is in line with previous results indicating that pregnancy at term and labor entail different inflammatory scenarios [16]. A healthy pregnancy involves an inflammatory process characterized by a fetal immune tolerance, while parturition is considered as an exacerbated maternal inflammatory sterile state that aids in expelling the fetus naturally [38]. Recently, we have demonstrated that a lower complexity and variability evaluated by SD has been associated with inflammatory conditions such as the endotoxemia induced by lipopolysaccharide in rats [39]. Therefore, the autonomic changes in fHRV reported here may result from a systemic inflammation that is also occurring in the fetus and not only for the mother. Further studies are needed to validate this potential explanation. Interestingly, the exploration of fHRV has provided a deeper understanding of the cholinergic anti-inflammatory pathway (CAP) [38].

To our best knowledge, this is the first longitudinal study to report fetal autonomic cardiac modifications in low-risk active parturition at term in comparison with the third trimester of pregnancy assessed by symbolic dynamics and using transabdominal recordings. Only few studies have investigated changes of linear or nonlinear fHRV indices between pre-labor and labor stages. Lim et al. reported that that the complexity of fHRV in the second stage of labor was significantly lower than that manifested previous labor, suggesting that parturition is a highly stressful situation for the fetus [40]. Other findings support the hypothesis of a strong role of the autonomic nervous system in the final minutes of parturition, which has been associated with a reduced central nervous system activity [41].

Contrary to our results, a recent study found increased parasympathetic fetal activity in labor in comparison to a pre-labor stage assessed by spectral analysis of fHRV [42]. This discrepancy is most likely caused by the differences in acquisition techniques, i.e., Doppler ultrasound used in such study vs. transabdominal ECG recordings used here. According to some authors, the comparison between results of fHRV analysis obtained by different methods should be interpreted carefully because fHRV signal acquisition and processing techniques vary widely [43]. It is known that Doppler ultrasound does not offer beat-to-beat fHRV time series, and fHRV time series are typically averaged over three beats and sampled at 4 Hz. Van Laar et al. [44] considered that fHRV high-frequency spectral analysis is only reliable if fHRV is acquired on a beat-to-beat basis. Other differences between studies may be explained by additional methodological differences such as the duration of segments analyzed (3 min of averaged fHRV vs. 5 minutes of beat-to-beat fHRV), and the evaluation of regular uterine activity in correspondence with the fHRV signal (tocodynamometer transducer vs. electrohysterographic estimation).

In this study, we introduced and compared the PRR5 with SD parameters because it reflects universal fetal developmental characteristics based on heart rate patterns [15]. Accordingly, we found significant correlations between the variations or differences of the linear indice ΔPRR5 (related to vagal modulations) and several Δ of SD indices between the TT and P. For example, a negative correlation of ΔPRR5 with Δ0V% (related to sympathetic activity) and a positive correlation with Δ1V% (related to parasympathetic and sympathetic activity). Additionally, the negative correlations with all the ΔPOLVAR measures show that a decreased vagal activity is associated with reduced fetal variability, being ΔPOLVAR5, the differences with the highest negative correlation with ΔPRR5. Interestingly, several studies have reported the correlation of SD indices with sympathetic and parasympathetic activities at the heart and vessels, especially in humans [33, 45]. In fact, Voss et al. found that POLVAR5 is a marker of reduced HRV and loss of complexity in patients with ischemic heart failure [46].

Finally, the measurement of fHRV by various methods remains a fascinating research topic, but it is not a recurrent clinical tool yet [47]. The findings of this study may contribute to the inclusion of the use of SD of fHRV as a potential biomarker to evaluate the maturation of the ANS in term fetuses at active parturition. Specifically, the SD of HRV has been evaluated by other authors as a biomarker of sleep apnea [48], stress [49], and inflammation [39].

Limitations

Our study has the following limitations. It is recognized that our final group of participants, providing fHRV measures at both third trimester and active labor, was small (n = 21) due to the inherent difficulty of successfully recording longitudinal fetal electrocardiograms. For that reason, we were not able to discriminate between male or female fetuses, which is relevant because the fetal sex seems to influence fHRV [50]. Additionally, in the present report, we did not acquire simultaneous ultrasound images to visualize fetal body or eye movements helping to discriminate among fetal behavioral states [51, 52]. The fetal behavioral states were likely different at TT compared to P, which could introduce significant effects on fHRV. Given the small sample size in the present study, our findings and interpretation should be confirmed and cautiously taken into consideration. Yet, a post hoc power analysis indicated that our study had adequate statistical power (at least 80%) to detect a minimum difference between both stages using some parameters.

We did not apply a nonlinearity test to evaluate the dynamics of fHRV. However, the presence of nonlinear dynamics in fetal RR fluctuations has been proved [53], thereby supporting the use of SD analysis as applied here. Further work will involve the reduction of redundancies for the binary patterns, as well as the application of pattern sets, pattern classes, and measuring of fetal breathing movements as is described by Cysarz et al. [10, 54, 55].

Conclusions

Our results suggest that during active labor, the fetal short-term cardiac variability is decreased (as indicated by WSDVAR, WPSUM02, and POLVARΔ); and it is possibly accompanied by decreased fetal vagal modulation (as indicated by PRR5) and higher fetal adrenergic response of the heart (as indicated by 0V%). These fetal autonomic modifications could become manifested in healthy fetuses at term to contend the stressful inflammatory environment of parturition. We then confirmed that the analysis of the SD as applied to beat-to-beat fHRV time series could be a potential source of clinical biomarkers to differentiate the fetal autonomic cardiac condition at different stages of pregnancy.

Declarations

Author contribution statement

C. Montalvo-Jaramillo and J. Reyes-Lagos: Conceived and designed the experiments; Performed the experiments; Wrote the paper.

A. Pliego-Carrillo: Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data.

M. Peña-Castillo and J. Echeverría: Conceived and designed the experiments; Analyzed and interpreted the data.

E. Becerril-Villanueva, L. Pavón, R. Ayala-Yañez and R. González-Camarena: Contributed reagents, materials, analysis tools or data.

K. Berg and N. Wessel: Analyzed and interpreted the data.

G. Pacheco-López: Conceived and designed the experiments; Contributed reagents, materials, analysis tools or data.

Funding statement

J. Reyes-Lagos was supported by the Mexican Public Education Ministry (Secretaría de Educación Pública: SEP) (511-6/18-8991). G. Pacheco-López was supported by the Mexican Science and Technology Council (CONACYT) (Grant No. 267761 DAAD-PROALMEX).

Competing interest statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.