Distance to Healthy Metabolic and Cardiovascular Dynamics From Fetal Heart Rate Scale-Dependent Features in Pregnant Sheep Model of Human Labor Predicts the Evolution of Acidemia and Cardiovascular Decompensation

The overarching goal of the present work is to contribute to the understanding of the relations between fetal heart rate (FHR) temporal dynamics and the well-being of the fetus, notably in terms of predicting the evolution of lactate, pH and cardiovascular decompensation (CVD). It makes uses of an established animal model of human labor, where 14 near-term ovine fetuses subjected to umbilical cord occlusions (UCO) were instrumented to permit regular intermittent measurements of metabolites lactate and base excess, pH, and continuous recording of electrocardiogram (ECG) and systemic arterial blood pressure (to identify CVD) during UCO. ECG-derived FHR was digitized at the sampling rate of 1,000 Hz and resampled to 4 Hz, as used in clinical routine. We focused on four FHR variability features which are tunable to temporal scales of FHR dynamics, robustly computable from FHR sampled at 4 Hz and within short-time sliding windows, hence permitting a time-dependent, or local, analysis of FHR which helps dealing with signal noise. Results show the sensitivity of the proposed features for early detection of CVD, correlation to metabolites and pH, useful for early acidosis detection and the importance of coarse time scales (2.5–8 s) which are not disturbed by the low FHR sampling rate. Further, we introduce the performance of an individualized self-referencing metric of the distance to healthy state, based on a combination of the four features. We demonstrate that this novel metric, applied to clinically available FHR temporal dynamics alone, accurately predicts the time occurrence of CVD which heralds a clinically significant degradation of the fetal health reserve to tolerate the trial of labor.


INTRODUCTION
Monitoring fetal heart rate (FHR) during labor is a common clinical routine worldwide, aiming to asses fetal well-being and ensure safe delivery. The main objective is to decide on timely operative delivery or uterine relaxation to prevent brain injury and adverse outcomes (1). In clinical practice, fetal well-being is assessed by obstetricians principally by visual inspection of cardiotocograms (CTG, bivariate time series of beat-per-minute FHR and uterine activity). The interpretation is guided by a set of rules combining a collection of features, aiming to probe various aspects of the CTG, such as baseline FHR, FHR variability and deceleration shape and timing as well as the relation of the various FHR features to the patterns of uterine activity. One such set of features and rules was defined by the International Federation of Gynecology and Obstetrics (2,3). Applying such procedure has however been documented as yielding significant inter-and intra-observer variability (4), one of many causes of the failure of the present FHR monitoring to predict fetal brain injury (5)(6)(7).
These short-comings in FHR monitoring during labor triggered significant efforts to develop computerized and automated assessment of FHR patterns intrapartum. Beyond the direct computation of the FIGO features themselves [cf., e.g., (8)(9)(10)], from digitized CTG usually sampled at 4 Hz in clinical practice, a large variety of features stemming from advanced signal processing and information theory tools has been computed for FHR assessment. These advanced features, however, have not reached performance benchmarks to lead to a consensus in the research and medical communities. These observations leave open a significant number of issues ranging from the choice of relevant FHR features and the construction of decision rules for such features to the assessment of the relationships between FHR time series and fetal well-being. Interested readers are referred to Georgieva et al. (11) (and references therein) for a recent (lack of) consensus overview and the interdisciplinary discussions.
Besides the need for large labeled databases to make machine learning on FHR data effective (12), a recurrent issue is associated with the ground truth being based on pH from the immediate post-birth umbilical cord pH measurements. However, it has been documented that fetal brain injury poorly correlates with measures of acidemia at birth such as pH (7,11,13). First, pH is only available after delivery hence when FHR is no longer available. Second, brain compromise due to hypoxia-ischemia can ensue when the fetal cerebral blood flow is persistently reduced, e.g., due to precipitous drop in cerebral perfusion pressure resulting from cardiovascular decompensation (CVD) (7,14,15).
In a recent series of experiments, to better assess the relations between FHR, systemic arterial blood pressure (ABP) and fetal health state (including the impact of chronic hypoxia) sheep fetuses were surgically instrumented and subjected to an umbilical cord occlusion (UCO) protocol in Frasch et al. (16), in a well-established animal model of human labor. CVD onset was observed at individually variable times, regardless the presence of chronic hypoxia (7,(16)(17)(18)(19). This animal experimental model generated the dataset used in the present study.
Consequently, based on this dataset, the goal of the present work is to assess whether FHR monitoring permits detection of the individual onset of CVD accounting for the presence of chronic hypoxia prior to the onset of UCOs in some fetuses. More particularly, we aimed to assess the sensitivity of FHR temporal dynamics, probed by four scale-dependent features, to CVD, metabolites and pH measurements.
We propose four features which all have the temporal scale of the signal as a parameter. We compute these four quantities from the whole FHR signal to probe its dynamics along the complete experiment. The first quantity measures the average variation of the FHR over the prescribed time scale. The second one measures the FHR variability over the time scale as the standard deviation. The third one is the ratio of the first two and provides a normalized version of the average variation. The fourth one is very similar to Approximate Entropy or Sample Entropy and provides a measure of the information content of the FHR signal at the given time scale. These four quantities can be computed with any signal and give robust results even with the clinically relevant low sampling rate of 4 Hz. This feature choice is also designed to allow computation within short-time time windows, thus permitting to achieve a sliding-window, time-dependent analysis of FHR, which may eventually be exploited to perform real-time FHR monitoring on noisy data.
We show the importance of coarse time scales (2.5-8 s) and construct an individual self-referencing "distance to healthy state" metric based on combination of the four features. We then demonstrate the use of the novel composite distance metric to predict individual CVD from FHR time series alone.

MATERIALS: SHEEP ANIMAL MODEL AND UMBILICAL CORD OCCLUSIONS
Fetal sheep model of labor and surgical preparation. The anesthetic and surgical procedures, postoperative care of the animals and the UCO model of labor have been previously described Frasch et al. (16). Briefly, 14 near-term ovine fetuses (123 ± 2 days gestational age (GA), term = 145 days) of the mixed breed were surgically instrumented. Animal care followed the guidelines of the Canadian Council on Animal Care and was approved by the University of Western Ontario Council on Animal Care.
Polyvinyl catheters were placed in the right and left brachiocephalic arteries, the cephalic vein, and the amniotic cavity. The fetal arterial lines were used for measuring ABP, sampling arterial blood gases, metabolites and cytokines. The fetal venous line was used for administration of fluids and postoperative antibiotics. Stainless steel electrodes were sewn onto the fetal chest to monitor ECG. A polyvinyl catheter was also placed in the maternal femoral vein. Stainless steel electrodes were additionally implanted biparietally on the dura for the recording of electrocorticogram, ECOG, as a measure of summated brain electrical activity [results reported elsewhere (16)]. An inflatable silicon rubber cuff (in vivo Metric, Healdsburg, CA) for UCO induction was placed around the proximal portion of the umbilical cord and secured to the abdominal skin. Once the fetus was returned to the uterus, a catheter was placed in the amniotic fluid cavity. Antibiotics were administered intravenously to the mother (0.2 g of trimethoprim and 1.2 g sulfadoxine, Schering Canada Inc., Pointe-Claire, Canada) and fetus and into the amniotic cavity (1 million IU penicillin G sodium, Pharmaceutical Partners of Canada, Richmond Hill, Canada). Amniotic fluid lost during surgery was replaced with warm saline. The uterus and abdominal wall incisions were sutured in layers and the catheters exteriorized through the maternal flank and secured to the back of the ewe in a plastic pouch. Postoperatively, animals were allowed 4 days to recover prior to experimentation and daily antibiotic administration was continued intravenously to the mother (0.2 g trimethoprim and 1.2 g sulfadoxine), into the fetal vein and the amniotic cavity (1 million IU penicillin G sodium, respectively). Arterial blood was sampled for evaluation of the fetal condition and catheters were flushed with heparinized saline to maintain patency. Animals were 130 ± 1 day GA on the first day of the experimental study.
Umbilical cord occlusion protocol. The experimental protocol has been reported (16,20,21). Briefly, all animals were studied over a ∼6 h period. Fetal chronic hypoxia was defined as arterial O2Sat <55% as measured on postoperative days 1-3 and at baseline prior to beginning the UCOs. The first group comprised five fetuses that were also spontaneously hypoxic (n = 5, H/UCO). The second group of fetuses was normoxic (O 2 Sat more than 55% before UCOs) (n = 9, N/UCO). As reported, after a 1-2 h baseline control period, the animals underwent mild, moderate, and severe series of repetitive UCOs by graduated inflation of the occluder cuff with a saline solution (16).
During the first hour following the baseline period, mild variable FHR decelerations were performed with a partial UCO for 1 min duration every 2.5 min, with the goal of decreasing FHR by ∼30 bpm, corresponding to a ∼50% reduction in umbilical blood flow (22,23). During the second hour, moderate variable FHR decelerations were performed with increased partial UCO for 1 min duration every 2.5 min with the goal of decreasing FHR by ∼60 bpm, corresponding to a ∼75% reduction in umbilical blood flow. Animals underwent severe variable FHR decelerations with complete UCO, i.e., ∼100% reduction of umbilical blood flow, for 1 min duration every 2.5 min until the targeted fetal arterial pH of <7.00 was detected, at which point the repetitive UCO were terminated. A summary of timings is reported in Table 1. These animals were then allowed to recover for 48 h following the last UCO. Fetal arterial blood samples were drawn at baseline, at the end of the first UCO of each series (mild, moderate, severe), and at 20 min intervals (between UCO) throughout each of the UCO series, as well as at 1, 24, and 48 h of recovery. When pH <7.00 was measured, the UCO were stopped and this time point noted as the end of the occlusions. We then obtained the precise pH = 7.00 time point by linear interpolation from this last measured pH value. All blood samples were analyzed for blood gas values, pH, lactate and base excess (BE) with an ABL-725 blood gas analyzer (Radiometer Medical, Copenhagen, Denmark) with temperature corrected to 39.0 • C. Plasma from the 4 ml blood samples was frozen and stored for cytokine analysis, reported elsewhere (24). After the 48 h recovery blood sample, the ewe and the fetus were killed by an overdose of barbiturate (30 mg sodium pentobarbital IV, MTC Pharmaceuticals, Cambridge, Canada). A post mortem was carried out during which fetal sex and weight were determined and the location and function of the umbilical occluder were confirmed. The fetal brain was perfusion-fixed and subsequently dissected and processed for later immunohistochemical study (25).
Data acquisition and pre-processing. A computerized data acquisition system was used to record fetal systemic arterial and amniotic pressures and the ECG signal (26). All signals were monitored continuously throughout the experiment. Arterial and amniotic pressures were measured using Statham pressure transducers (P23 ID; Gould Inc., Oxnard, CA). Fetal systemic ABP was determined as the difference between instantaneous values of arterial and amniotic pressures. A PowerLab system was used for data acquisition and analysis (Chart 5 For Windows, ADInstruments Pty Ltd, Castile Hill, Australia). Pressures, ECOG and ECG were recorded and digitized at 1,000 Hz for further study. For ECG, a 60 Hz notch filter was applied. R peaks of ECG were used to derive the heart rate variability (HRV) times series (26). The time series of R-R peak intervals were then uniformly resampled at 4 Hz (26). A representative FHR signal is shown in Figure 1A: a visual assessment of the whole FHR trace reveals that FHR variability increases when UCO strength is increased.
Metabolites data (pH, lactate and BE) is obtained by blood sampling performed at specific times during the experiment (vertical dashed lines in Figure 1B). In order to have metabolites data at any time, we assume a linear drift between two successive measurements and consequently perform a linear interpolation between two measurements times. We thus obtain a piece-wise linear time series sampled at 4 Hz, depicted in Figure 1B as a black curve. Using this interpolated data, as noted above, the time t pH when pH = 7.00 is computed in each fetus as indicated in Table 1.
Fetal cardiovascular decompensation (CVD). CVD has been reported in detail in Frasch et al. (15,16) and Gold et al. (27). The visual representation of CVD can be found in these publications, e.g., in the Figure 2 in Frasch et al. (16). The reader can readily observe the pronounced pathological hypotensive responses to the UCO-triggered FHR decelerations during the CVD. This behavior is in stark contrast to the normally observed ABP increases during the occlusions which compensate the hypotension caused by the FHR decelerations. As we reported, once this pattern conversion from hypertensive to hypotensive responses occurs, it persists until the UCOs are stopped. Its effects are also seen directly in the brain electrical activity (15,16). It is hence easy to reliably visually identify the timing of the onset of CVD in each recording. Consequently, during UCOs, by expert visual inspection, we noted the individual time point t CVD at which three successive hypotensive ABP responses to UCO-triggered FHR decelerations occurred. Quantitatively, with hypotensive ABP response we refer to the failure of ABP to rise during UCO-triggered FHR deceleration above the preceding baseline value when compared to the average ABP rise during the  UCO series prior to the CVD. We refer to this animal-specific time point as the ABP sentinel corresponding to the timing of CVD. As an illustration, t CVD is reported in Figure 1B as a vertical black line.

Sliding Window Analysis
Analysis of FHR and metabolites data are performed in sliding time-windows of size T = 20 min. The time-windows are shifted by dT = 5 min, thus implying a T − dT = 15 min (75%) overlap. The k-th time-window thus corresponds to time ranging in [kdT, kdT + T]. This sliding window analysis permits the assessment of the temporal evolution of cardiovascular responses to changes in UCO strength.

Scale Dependent Features
Using FHR, x t , four quantities, whose definitions rely on the choice of a time scale τ , are computed, for each time-window k: increment mean m k (τ ), increment standard deviation σ k (τ ), the corresponding Student ratio R k (τ ), and the entropy rate h k (τ ). The first one, m k (τ ), measures the average change of the FHR signal over the time lag τ . It can be pictured as the derivative of the signal on the time-scale τ , averaged in the timewindow k. The second one, σ k (τ ), is the standard deviation of the signal estimated on chunks of the signal of duration τ , and then averaged in the time-window k. The third one, R k (τ ), is a normalized version of the first quantity m k (τ ): the average variation is now expressed in standard deviation units, before being averaged along the time-window k. The last one, h k (τ ), in a measure of the information or complexity of the signal at scale τ , similar to Approximate Entropy or Sample Entropy, estimated at scale τ over the time-window k. As described below, the first three quantities are averages over the time-window k of dynamical quantities m t (τ ), σ t (τ ), R t (τ ) defined at scale τ ; these quantities, along with FHR, are shown for illustration purposes in Figure 2 for arbitrarily chosen 20-min window, time scale τ and animal.

Mean Variation (or Local Trend) at Time-Scale τ
For all t in window k, the average increment over a time scale τ is computed as: These m t (τ ) are then averaged across window k, for all non- where ⌊T/τ ⌋, the floor of the fraction T/τ , indicates the number of time-intervals of size τ available in the time-window of size T. An illustration of the methodology is given in Figure 2B for the window k = 0: values of m t (τ ) are depicted in black, and the single value m k (τ ) is represented in red. The quantity m k (τ ) measures the average variation -either an increase or a decrease -of FHR on the time scale τ . The average is indeed a double average: first over all time scales smaller than τ , according to Equation 1, and second over all available intervals available in the kth time-window of size T=20 min, according to Equation (2). m k (τ ) can also been interpreted as the averaged derivative of the signal after a low-pass filtering using a finite impulse response with cut-off frequency 1/τ .

Standard Deviation at Time-Scale τ
This indeed is nothing but the variance of x t , computed over the set of values in the time-interval [t − τ ; t]: We then average its square root over the ⌊T/τ ⌋ non-overlapping time intervals of size τ available in the kth time-window of size T: This quantity measures the average -in the kth time-window of size T -amplitude of the fluctuations of x t over τ consecutive points. The methodology is illustrated in Figure 2C.

Normalized Local Trend at Time-Scale τ
The Student ratio, or normalized local trend, at time-scale τ , is defined for each time interval [t − τ ; t], as: It is averaged across all available non-overlapping intervals in the kth time-window: This quantity, up to a factor √ τ , would correspond to a random variable drawn from the distribution of the t-value if the data x t were independently drawn from a Gaussian distribution. It can be interpreted as the average variation over a time step τ , normalized by the local standard deviation; as such, it provides a normalized measure of the trend of the signal x t to depart from its expected value when observed across a duration τ .

Entropy Rate at Time-Scale τ
One commonly used feature in heart rate analysis, both for adults and fetuses, is sample entropy (SampEn) (28)(29)(30), an elaboration on approximate entropy (ApEn) (31,32). It was shown recently that the entropy rate provides a related tool to probe FHR with better performance than ApEn or SampEn to detect acidosis (33)(34)(35).
The entropy rate of order 1 in the kth time-window at timescale τ is defined as: where denotes the Shannon entropy (36) of either a vector x = (x t , x t−τ ) or a scalar x = x t . h k (τ ) is computed using all the pairs of points (x t , x t−τ ) available in the k-th time-window, and following Theiler's prescription (34) to avoid spurious correlation. h k (τ ) measures the extra information conveyed by the vector (x t , x t−τ ) when (x t ) is known, or in other words, the extra information given by the knowledge of the signal at an earlier time t − τ . The entropy rate probes the dynamics of the signal, and to better focus on this dynamical aspect, we compute it on the normalized signal (x t − x t )/ (x t − x t ) 2 , where . stands for the time average on the window of size T.

RESULTS AND DISCUSSION: FEATURES, TIME-SCALES AND DISTANCE TO HEALTHY STATE
In section 4.1, the four features -computed in overlapping time-windows -evolution in time are firstly presented and studied with respect to their relations to UCO strength. Because these features are computed at a given time-scale τ , they offer a description of the FHR dynamics at this time-scale. We thus explore the correlation between the features at a given time-scale τ and the measured values of the metabolites -including the pH. This global analysis, presented in section 4.2, is performed using all available time-windows and all available animals. We then reduce the dimensionality of the analysis by averaging results over the long-term time-scales, as defined and presented in section 4.3. This allows us to examine more clearly how the features evolve jointly with the UCO strength for the entire cohort, while quantifying the variability between animals. We then examine quantitatively in section 4.4 how these long-term features correlate with metabolites. We then combine them in an appropriately normalized vector; we are then able to describe the large variability across the subjects in the population as the variability of this vector in the early stages of the experiments. This allows us to define a measure of the degradation of the health state of an animal as the distance from healthy state. Finally, we propose in section 4.5 to use this "individual" distance as a novel indicator -or sentinel -to alert for the degradation of the health status due to CVD. We also show that this indicator/sentinel matches very well with pH measurements.

Features and UCO Strength
We first examine on a single animal how the four FHR features evolve throughout an experiment, depending on the time-scale τ . The values obtained in the k-th time-window [kdT; kdT + T] are assigned to the date t k = kdT + T/2 at the center of the timewindow. The dynamical evolutions of m k (τ ), σ k (τ ), R k (τ ) and h k (t, τ ) are depicted in Figure 3 for a large band of time scales τ . Such a time-scale representation reveals qualitatively that when the UCO strength is increased, m k (τ ), R k (τ ) and h k (τ ) decrease along time, while σ k (τ ) increases along time. This agrees with the previous studies where the decrease of the entropy rate h k (τ ) was associated with fetal acidosis (33)(34)(35).
Qualitatively, although the four features barely evolve in time for smaller values of τ (below 2 s, bottom of the images in Figure 3), a noticeable time evolution can be observed for large values of τ and especially in the severe UCO regime. To better observe the dependence of the four features on the scale τ , we plot in Figure 4 their evolution with τ for the time points when blood sampling was performed. Figure 4 therefore presents the evolution of the four features along the vertical color lines indicated in the images of Figure 3.
We observe in Figure 4 that the evolution of m k (τ ) is rather linear in τ , but the slope depends on the time, and hence on the UCO level. We observe almost no evolution of R(τ ) with τ , but the value of R(τ ) depends on time, so on the UCO level. On the contrary, both σ (τ ) and the entropy rate h(τ ) present a distinct change of their evolution with τ below and above τ = 2.5 s, which emphasizes the distinction between short (< 2.5 s) and large (> 2.5 s) time scales, in accordance with previous literature (26,37,38). We use this information on the time scales as follows.

FHR Features, Arterial Metabolites and pH
We now examine, for a fixed time-scale τ , how the features relate to the health state of the animal, as described by the metabolites and pH. To do so, we use all time-windows of size T on one side, and interpolated metabolites data on the other side. We compute the correlation between any of the four features (for a fixed τ ) and any of the biochemical measurements, by averaging over all time-windows (average over k) and over , depending on the scale τ for animal 473726. The time in abscissa is kdT + T/2, the location of the kth time-window of size T where the quantity is computed, and the ordinate represents the scale τ . Vertical color lines indicates the times at which blood sampling was performed (same color code as in Figure 1: green in the mild UCO regime, magenta in the moderate UCO regime, and red in the severe UCO regime). In the severe UCO regime and for larger time scales τ , stronger variations are observed.  Figure 3: each curve corresponds to a time-window of size T for which a fetal arterial blood sample was taken. The color of the curve represents the corresponding UCO level, with the same color code as in Figures 1, 3: blue is the baseline prior to any UCO, green in the mild UCO regime, magenta in the moderate UCO regime, red in the severe UCO regime, and then black in the recovery regime (after UCO). Vertical black dashed lines indicate the time-scales 2.5 and 8 s.
all animals. Results are plotted in Figure 5 as a function of the scale τ .
As suggested by Figure 4 and confirmed by Figure 5, we can isolate two bands of time scales: shorter scales τ < 2.5s (i.e., high frequency band, above 0.4 Hz, short term time scales, labeled ST) and larger ones τ > 2.5s (low frequency band, below 0.4 Hz, long term time scales, labeled LT).
For any of the four features and any of the three biochemical measurements, the correlation in the range [2.5 − 8] s is not only the largest-in absolute value-but also the most stable: it fluctuates less and does not depend much onτ . Above 8 s, all correlations decrease in absolute value, which may be attributed in part to poorer statistics: the number ⌊T/τ ⌋ of available time-intervals of size τ in a time-window of size T decreases, which impacts the averages, see, e.g., Equation 2).
As a consequence, we choose in the following to restrict the long term (LT) range to τ ≤ 8 s in order to have enough statistical power.

Long-Term Scales Averaged FHR Features
For the sake of simplicity, we now eliminate the dependencies of our features on τ and focus on the LT range. To do so, we compute the area under the curve (AUC) of our four FHR features in the range 2.5 < τ < 8 s. For a given time-window indexed by k, we compute: and we define accordingly σ LT k , R LT k and h LT k . These features depend only on time, via the index k of the time-window in which they are computed.
Time evolutions of these four LT features are depicted in Figure 6 for the complete set of 14 animals. For some animals, there may be missing data due to experimental conditions, and hence there may be less consecutive time-windows of size T available than expected in a given UCO region; in that situation, we have then chosen to assign the dark blue color (arbitrary) for the quantity-see, e.g., the second line (a hypoxic animal), where no data is available in the mild UCO region, and only 4 windows are available in the severe UCO region.
Using the first column (on the left of the vertical green line) of each subfigure as a reference, we observe that every quantity evolves as the UCO strength is increased. Although very few changes are observed in the mild UCO region, much larger variations are observed in the severe UCO region. After the stopping of UCOs (on the right of the vertical black line), we observe that the four features seem to regain their original value, which we interpret as indicating the recovery of the animal, typically after 1 window of size T, so typically 20 min after the end of UCOs.

Distance to Healthy State, Metabolites and pH
We now explore how our four FHR features relate to the metabolites' levels, and especially to the pH value, which is a widely used indicator of fetal well-being. We report in Table 2 the correlation coefficient between each of the four features m LT , σ LT , R LT and h LT on one hand, and the three biochemical measurements pH, BE and lactate on the other hand. To increase the statistical power, we use all available time-windows of size T and so all linearly interpolated values of the three biochemical measurements.
We observe that the four FHR features correlate well with the pH and BE, while the correlation with the lactate is smaller. All features but σ LT -the LT amplitude of fluctuationshave a correlation coefficient with pH that is at least 0.50, and a correlation coefficient with BE that is at least 0.43. This is interesting, as R LT appears strongly correlated with m LT while relatively uncorrelated with σ LT .
We believe that each of the four FHR features contributes a particular piece of information about FHR and we therefore aggregate them as follows. For a single animal and a single time-window indexed by k, we consider the vector where each component is normalized by its standard deviation computed over all animals and over all available time-windows of size T. The four values (m RMS LT , σ RMS LT , R RMS LT , h RMS LT ) used for this normalization are hence the same for all animals and all time-windows; they are reproduced in the third line of Table 3.
For a given animal and for a given time-window indexed by k, we use the L 2 norm in R 4 to project any vector u k into a positive real number u k as follows. For each animal, we assume it is in a healthy condition when the experiment is started (so the FHR is fluctuating around the baseline) and we use the first time window of size T as a reference. We thus define the distance between u k which describes the state in the k-th time-window and u 0 which describes the state in the first time-window [0; T]: We interpret this distance D k for a single animal as a measure of the deviation from the animal's "healthy" state during the experiment. We report in Table 3 global statistics -obtained by considering all animals -of the four FHR features used as the four components of the vector u t .
The third line of Table 3 reports the values m LT RMS , σ LT RMS , R LT RMS and h LT RMS used to normalize the vector u k . Their amplitude is notably different, and the normalization is necessary to ensure that each component of u k contributes equally to its norm || u k ||. Whereas, this normalization uses all available data (using all times and all animals at once), it is important to stress that we have accounted for the large variability from one animal to another by defining D k with a reference relative to the very animal  Table 1); chronically hypoxic ones are above and normoxic ones are below. Each column represent a time-window of size T and time is increasing from left to right. The region of mild UCO starts at the vertical green line and lasts for four windows, up to the vertical magenta line, followed by 4 windows in the moderate UCO region, and then up to six windows in the severe UCO region, and up to three windows in the recovery region.  under consideration. The variability of the reference point can be seen in the fourth line of Table 3: it accounts for a large part of the RMS values used in the normalization. Comparing the first two lines of Table 3 brings an additional observation leading to the same conclusion: the position of the healthy state u 0 is on average over the animals (second line of the table) sensibly different from the position of u k averaged over all animals and all times (first line of table). Using D k instead of || u k || removes a large part of the inter-animal variability and definitely improves the relevance of the distance, as measured by the correlation with the metabolites, see Table 2.  We present in Figure 7 the 14 trajectories of the vector u(k) in its phase space, for the complete cohort. u k has 4 coordinates so there are 6 different projections in a plane defined by two variables. Each subplot in Figure 7 corresponds to one of these possible projections. Along each trajectory, the color changes to indicate the interpolated pH value. Although the trajectories wander in a large region of the phase space, their color-coding seem to only depend on the distance from the origin: blue (larger pH) close to the origin, and orange or red (lower pH) outside of the circle defined by D = 2. During an experiment, the UCO's strength increases and, as a consequence, the pH decreases. We observe that the distance D k appears to increase concomitantly, and more precisely we observe its correlation with the pH value. The correlation coefficients between the distance D and the biochemical measurements, computed over all animals, are reported in Table 2 (gray-colored cells). We observe that among all FHR features we have computed, the distance D is the one that is the most correlated with pH, as well as with the other metabolites. This confirms that using all four FHR features simultaneously-by considering the vector u-not only mitigates the various evolutions of single features with the metabolite value but also aggregates their correlations.

Distance to Healthy State as a New Sentinel for CVD
To better illustrate the relation between the dynamical featuresespecially the distance D k -on one hand, and the health status as assessed by the metabolites-especially the pH-and blood pressure responses to UCOs on the other hand, we examine in detail in Figures 8, 9 how these are co-evolving for each individual animal. From now on, we discard any indication of the UCO level. Figure 8 presents jointly the time-trace of the FHR signal (Figures 8a,e) and the evolution of the distance D (Figures 8b,f), which we color-code using the pH value, as in Figure 7 for two typical normoxic animals. We also present two projections of the trajectory of u k in its phase space in order to illustrate the evolution of each of the four quantities.
Because the vector u k has been carefully normalized, and appropriately centered to define the distance D k , this last quantity has no dimension and can be compared to absolute values. The particular value D = 1 (gray circle) defines the standard deviation range in a healthy situation and the value D = 2 (black circle) corresponds to variations with an amplitude of 2 standard deviations. Looking at the trajectory projections in Figures 8, 9, we see that during the early stages of the experiments the trajectory remains close to the origin, hence the distance D remains small, albeit fluctuating, and the pH value remains close to its normal value (bluish color, indicating a pH close to 7.4). More interestingly, we see that when the pH decreases down to 7.2 (greenish color), the trajectory usually reaches the black circle, hence the distance increases up to 2. Finally, we observe that when the trajectory is outside the black circle, hence D > 2, the pH has low values but more importantly, values of pH≤7.00 (orange to red color) are only observed on the trajectory much later after the trajectory wandered outside the black circle.
The very same observations can be made for hypoxic animals, see Figure 9 for two examples.
We now examine quantitatively for each animal the time evolution of the distance from its own reference healthy state. In every time-window [kdT, kdT + T] of size T indexed by k, we have a value D k which we assign to time t = kdT + T/2; we plot the time-evolution of the distance in Figures 8b,f, 9b,f with a color that indicates the pH value. This allows for the following interesting observations. First, before the occurrence of UCO, the distance fluctuates with a typical standard deviation of 1. This confirms that the normalization step is valid, albeit it uses values averaged over all animals and all available timewindows. Second, we see that the distance D is substantially larger when UCOs are performed, and more precisely, we see that D increases as the UCO strength is increased or UCOs are being applied to the animal with the same strength but for a longer period of time. As such, the distance D seems to be a good representation of the health condition of the animal. More interestingly, we showed ( Table 2) that D is highly correlated with the pH value throughout the complete experiment, but we now observe that a large value of D, above 2, proves to be a very good indicator of a low pH value signaling an acidemia with pH < 7.2.
To further test the ability of the distance D to alert on the fetal condition, we now try to relate the large values of the distance D to the onset of the fetal CVD, i.e., failure of the fetus to mount a hypertensive arterial blood pressure response to UCOs and the UCO-induced FHR decelerations, a prerequisite to maintaining an adequate cerebral perfusion pressure. To do so, we use t CVD as the reference time when CVD occurs (vertical red line in Figures 8b,f, 9b,f), which offers a valuable benchmark for an early detection of hypotensive blood pressure response.
t CVD appears on group average for a pH of 7.20 and 60 min prior to pH nadir of <7.00, but shows a considerable interindividual spread. To quantify the hypotensive behavior at t CVD , we report the individual pressure differential ABP at t CVD in Table 4 following the same approach as reported in Table  1 of Frasch et al. (16). Here, ABP=ABP max − ABP is the difference between ABP max , the maximal ABP during a UCO, and ABP , the mean ABP between UCOs. The average of ABP is 4±6 mmHg for hypoxic fetuses and 4±8 mmHg for normoxic fetuses. Overall, we see no difference in ABP between hypoxic and normoxic groups (p = 0.97). The corresponding drop in ABP during the CVD period compared to the preceding UCO period is −19 [−24; −1] mmHg, i.e., during CVD, fetuses failed to mount hypertensive response to UCO-triggered FHR decelerations with a median drop of 19 mmHg compared to ABP preceding the CVD. These values clearly indicate the pathological hypotensive responses of the sheep fetuses during the UCOs at t CVD and onward until the end of the UCOs. The noted inter-individual variability in ABP values is subject of ongoing research.
Here we see on phase-space projections in Figures 8c,d,g,h, 9c,d,g,h that the criterion D ≥ 2 offers a similar early alert on the deterioration of the animal condition with regard to CVD timing. Looking at either the phase space representation or the time traces of the distance D, we see that this quantity evolves continuously in time, on typical time-scales larger than 20 min, the duration we have chosen to compute our quantity. The distance increases over the duration of the experiment and one can easily measure the time t D at which D crosses the value D = 2 (red circle, or horizontal red line in Figures 8, 9). Unfortunately, the distance D is very sensitive and it can be seen on the examples that it is possible for D to reach values larger than 2 early in the experiment. To overcome these events-and, hence, to make our new sentinel less sensitive-we arbitrarily adjust our criteria and require D k > 2.5 for at least 3 consecutive time-windows, so for a long enough duration of about 40 min. Table 1 presents the various timings corresponding to the various UCO regimes for each animal, together with an estimate of the pH nadir time, while Table 4 presents a summary of our findings, together with the CVD time (ABP sentinel), the two of them appearing before pH≤ 7.00. Each trajectory corresponds to an animal and is colored to indicate the pH value at the time k: in this way, we observe the joint temporal evolution of u k and of the pH throughout the experiment. Trajectories have been centered by subtracting u 0 , according to Equation (11) to account for the variability between animals: the thick black dot at the origin thus represents the starting point of all trajectories. The gray circle corresponds to D = 1 and the black circle to D = 2.
The agreement between the CVD time and the distance time is very satisfying: the difference between t CVD and t D is not only always smaller than the difference between t pH and t CVD , but also smaller than 20 min, the size of the time-windows we have used.
However, for one animal (number 473360, last line in Table 4), a large discrepancy is observed. A closer examination of both the data and our distance for this animal is given in Figure 10 and allows us to discuss the sensitivity of our measure. We have used the 4 Hz FHR dataset which was also studied in earlier literature. This dataset is obtained from the R-R intervals data at 4 Hz, which is interpolated from the raw ECG-derived R-R intervals data recorded at 1,000 Hz. As can be seen in Figure 10, the genuine 1,000 Hz dataset (in red) is missing some values during short intervals and the resampling process, which uses splines interpolation, creates arbitrary values for the 4 Hz FHR dataset (in black) within such intervals. This results in additional values which exhibit large and fast fluctuations which are nonphysiological. Whereas, most of these do not impact the value of the distance D (see Figures 10e,f,g,h), there is a time interval (at about t = 172 s, see Figures 10b,d where D is unexpectedly large, reaching a value around 4. This is concomitant with a sharp drop in FHR, as can be seen in Figures 10a,c. This sharp drop is exacerbated on the 4 Hz signal compared to the 1,000 Hz signal, and is very localized in time, which leads to a later decrease of D, contrary to the pathological situation reported in Figures 10g,h where D remains at a large value. As a consequence, we obtain a false positive sentinel time t D which corresponds to this event and is hence much earlier than t CVD , although in agreement with previously reported results using the same 4 Hz FHR dataset (19). We conclude that splines interpolation should be avoided, and we suggest instead not to add or create artificial data points when genuine data is not available. Additionally, each of the quantities we propose, and hence the distance D, can still be computed, as they are all robust with respect to missing data, as seen for example in Figures 8e,f.
The robustness with respect to missing data is 2-fold. First, all quantities we compute do not require equi-sampled data: this is in contrast to a power spectrum for example, where missing points prevent the estimation and jeopardize the estimated value if using interpolated values. For our quantities, missing data only impacts the number of points used to compute averages, as can be seen in Figure 2. Second, having missing data only reduce the number of points over which statistics are computed: a reduced number of points increases the bias and the variance of the estimators. As can be seen in Figures 10b,d,f,h, the distance D evolves smoothly which suggest the standard deviation is not strongly impacted. However, one may wonder if an increased bias impacts the reported values, especially when a lot of data is missing. We report in Table 5 the average fraction of missing data points in a time-window of size T=20 min: increasing the UCO strength is typically associated with an increase of missing data. Let's focus on the entropy rate h k (τ ), which is algorithmicaly the most complex quantity: it has been reported that the bias of h k (τ ) not only behaves as 1/ √ N (39,40), similar to the bias of a sliding average over N points like m k (τ ), but also that this bias is small. A time-window of 20 min should contain N = 20 × 60 × 4 = 4, 800 points, and even a reduction of 50% of available data points should leave more that 2,000 points so a bias smaller than 1% (40). We are thus confident that the reported results are not an indirect measure of the number of missing data points.
A deeper examination of each experiment, using animal's systemic arterial blood pressure data, should clarify the relationship between the increases of the distance from healthy condition and the incipient arterial hypotension. This work is out of scope of the present article which focuses on the dynamics of the clinically relevant 4 Hz-sampled FHR signal.
As such, we propose that the distance D has the potential to serve as an individual biomarker of the incipient CVD, i.e., an early sentinel of the fetal brain injury.
In the present work, four FHR features, whose definitions depend on the timescale, are computed on the whole FHR dataset derived from an animal model of human labor to quantify the evolution of FHR temporal dynamics. That means that in our approach we do not rely on considering UCO periods only, but are able to process the entire FHR signal as it would be available in real-time in clinical setting. These quantities are local statistical averages that probe the variation, the amplitude of fluctuations and the information content at a given time scale. They are purely statistical quantities that can be computed even when some data points are missing. Firstly, we qualitatively related the variations of such timescale-dependent quantities to the UCO strength; secondly, we quantitatively computed their correlation to metabolite and pH measurements.
As to the etiology of CVD, we propose a role for the Bezold-Jarisch reflex, a vagal cardiac depressor reflex, as part of a complex dynamic interplay, based on the observations of acidemia-triggered inflammation in fetal sheep (16,17,25), and studies in adult species linking rising systemic acidemia and inflammation with worsening cardiac contractility, impaired beta-adrenergic and potentiated bradycardic responses (25,(73)(74)(75)(76)(77). We suggest that the integrated ability of the four FHR features introduced in this study to track the individual evolution of acidemia and cardiovascular responses stems from capturing the individual complex interplay of the vagally mediated sensing of acidemia and the Bezold-Jarisch reflex, i.e., also vagally mediated intermittent hypotensive ABP responses to UCO-triggered FHR decelerations. This hypothesis needs to be validated in specifically designed animal experiments, for example by repeating the experiments underlying the present study with the variation of performing cervical bilateral, left or right vagotomies. This would allow evaluating the contribution of the vagus nerve to the dynamic interplay between the progressive systemic acidemia, the ensuing systemic inflammatory response, accounting for vagus nerve's lateral asymmetry, to the evolution of FHR decelerations and ABP responses over the period of worsening UCOs comparable in duration to stages 1 and 2 of pushing (18).
We show the relevance of timescales ranging in [2.5 − 8] seconds (equivalently [0.125 − 0.4] Hz in frequencies) for early detection of both acidemia and CVD, matching the scales classically used in FHR analysis and referred to as longterm (57,63). We observed that reduced pH closely relates to larger m LT which may be interpreted as an increase of baseline FHR (62,63), and lower entropy rate h LT , in agreement with earlier findings reported in the literature (34). More importantly, a per-individual distance metric was constructed from these four (population-normalized) features to quantify a self-referencing departure from a healthy state for each subject independently. Such a definition raises two issues. Firstly, it requires, as is often the case, that monitoring is started early enough while the fetus is still in a healthy condition, so as to create a self-reference to normal on a per-individual basis. If fetuses are already in distress when monitoring is initiated, the distance, albeit increasing with distress, may fail to detect CVD correctly. Secondly, the definition of the vector, and hence the distance, requires a normalization, which is performed in the present work at the population level, i.e., using an average across subjects. Although such an average should converge rapidly with the population size, this dependence requires further investigations.
It has been documented that sheep fetuses have an individual cardiovascular phenotype in their responses to increasing acidemia due to repetitive intermittent hypoxia (15). Chronically hypoxic fetuses have diminished cardiovascular defenses to hypotensive stress (78). Studying the same dataset, we demonstrated that under the conditions of repetitive UCOs and in comparison to the fetuses who were normoxic on the onset of the UCOs, the hypoxic fetuses exhibit accelerated acidosis (21), altered temporal profile of neuroinflammation following UCOs (20) and deceleration reserve (79). In the For each time-window of size T = 20 min, we divide the number of missing data points by the expected number of points (=20×60×4), and we then average this ratio over all time-windows available in a given regime. present study, we could not identify any influence of the phenotype (normoxic/hypoxic) in any of the metrics including the distance D. We are currently exploring whether the initial vector u 0 may contain such information. Conversely, the finding that the presented approach functions well without the consideration of a pre-existing hypoxia or pattern of labor contractions (UCO severity) is an additional bonus from the clinical viewpoint. Lastly, we recognize that the group of chronically hypoxic animals may have failed to recover from surgical instrumentation adequately, i.e., they were already decompensating rather than becoming "spontaneously" hypoxic for reasons of utero-placental dysfunction preceding the surgery. Overall, the constructed distance proves able to detect accurately the occurrence of acidemia and CVD from the analysis of FHR only, and without recourse to pH. This opens the route to investigating the relevance of such metrics in clinical practice, as it is non-invasive and much faster than biochemical measurements like pH. Further, for practical purposes, the present studies show that the computation of features and distance is robust to FHR sampled at 4 Hz and, to some extent, to missing data. Also, the FHR features and distance are computed in sliding-windows, permitting an on-line and quasi-real time analysis of the evolution of the dynamics of FHR, and thus in relation to a local health state of the fetus. The extent to which the 20-min sliding-window size, chosen here for proof-of-concept developments, can be further reduced to 10 or 5-min, is under investigation. In conclusion, we propose a real-time FHR-based metric predicting CVD which should be of a great help for health practitioners managing the delivery.

DATA AVAILABILITY STATEMENT
The data analyzed in this study is subject to the following licenses/restrictions: experimental data may be made available on request (to MF). Requests to access these datasets should be directed to mfrasch@uw.edu.

ETHICS STATEMENT
The animal study was reviewed and approved by University of Western Ontario Council on Animal Care/Canadian Council on Animal Care.

AUTHOR CONTRIBUTIONS
SR designed the signal processing tools, ensured their practical implementation, conducted their application to data, and the analysis and interpretation of the results. He also prepared the figures and tables reported in the article. NBG contributed to the design of the signal processing tools, to the analysis and interpretation of the results, and to the writing of the manuscript. PA contributed to the interdisciplinary connections between signal processing and medical doctor teams, to the interpretation of the results, and to the writing of the manuscript. NG contributed to the measurements and the manuscript. MF designed the experiment, conducted the surgeries, and measurements. He also performed the expert visual detection of the cardiovascular decompensation used as ground truth here. He contributed to writing the article and to the interpretation of the results. All authors contributed to the article and approved the submitted version.