|Home | About | Journals | Submit | Contact Us | Français|
Heart rate variability (HRV) is an important index of cardiovascular state . In healthy adults, heart rate generally increases during inspiration and decreases during expiration. Such influence of respiration on the heart beat variation through mechanical and autonomic pathways is called “respiratory sinus arrhythmia” (RSA) . Adults have a breathing pattern with a mean respiratory rate of approximately 12 breaths per minute (0.2Hz) and, assuming a linear influence, these oscillations are directly reflected in the variability of the heartbeat. The RSA can be evidenced as a peak in the HRV spectrum in the normal respiratory frequency range. Such a well defined peak in the HRV spectrum reflecting RSA indicates the presence of cardio-respiratory interactions in adults .
In infants, cardio-respiratory interactions are considered as an important indicator of the level of maturation of vagally-mediated autonomic influence on the heart , although the exact relationship between HRV and respiration in preterm infants remains elusive [5–11]. The normal respiratory rate of infants is approximately 60 breaths per minute (1Hz). However most preterm infants have irregular breathing patterns with periodic breathing and pauses in breathing (apnea) that introduces frequencies lower than the normal range . As a result, the respiratory modulation of heart rate, if any, will be occurring at various ranges of frequencies, from ~1Hz (normal breathing) and below. Hence, in preterm infants, the signature of RSA, (i.e. the peak in the power spectrum of HRV at the normal breathing frequency of ~1Hz) may not be observed in the HRV spectrum due to irregularity in breathing. Notably, a further complication in relying on traditional spectral analysis is that heart rate fluctuations may exist at the respiratory frequencies even in the absence of respiration .
A wide collection of mathematical tools have been successful in quantifying important cardiovascular control mechanisms in adults [14–17]. In particular, application of frequency domain methods to peak-to-peak (RR) series detected from the electrocardiogram (ECG) signal alone and together with cardiovascular covariates such as respiration and blood pressure have led to highly refined models and analysis tools, as well as successful efforts in defining specific standards, all now traditionally classified as “heart rate variability studies” . Despite the valuable research devoted to the study of cardiovascular mechanisms, there is still need for specifically target standards and methods to assess the cardio-respiratory functions in the early stages of life.
As reported in the literature , the standard low frequency (LF: 0.04–0.15 Hz) and high frequency (HF: 0.15–0.4 Hz) ranges classified for adult HRV analysis do not apply to the analysis of HRV in infants. Attempts have been made to adapt new standards to newborn physiology [5–11] . Different frequency ranges within 0.01 Hz to 1.5 Hz are used for frequency domain analysis of HRV in preterm infants. Generally, any frequency above 0.2 Hz has been classified as HF in the case of infant HRV . However, a standard classification for infants has not yet been established.
The coherence function is a quantitative measure reflecting the strength of the linear interaction between HRV and respiration in adults. Coherence is traditionally calculated as the cross- spectral density between HRV and respiration normalized by the corresponding auto-spectral density functions . The coherence function takes values between zero, indicating absence of linear interactions, and one, indicating exclusive linear interactions. However it has been pointed out that the estimation of coherence using cross-spectral density does not account for causality of the two time signals considered . The method assumes that the two signals interact in an open loop in which respiration has a unidirectional influence on the HRV whereas there may be information flow from HRV to respiration in preterm infants. For example, fluctuations in ventilation results in fluctuations in arterial pH and pCO2, which in turn affect ventilatory drive via the central and peripheral chemoreceptors. The time delay and dynamics in this feedback system is modified by fluctuations in systemic and cerebral circulations, which are influenced by HRV. Correlated fluctuations in heart rate, arterial pressure cerebral circulation have been recorded in infants , and have been related to variability in breathing patterns of preterm infants . Therefore information flow from HRV to respiration could be important in some infants making the interactions between HRV and respiration bi-directional and closed loop.
Another possible source of information flow from HRV to the respiratory control of ventilation could relate to neural reflexive baroreceptor influences on central neural respiratory activity (c.f. ). As a result, even with a non-negligible information flow in the opposite direction any interactions estimated using open loop paradigm can significantly differ from the actual interactions occurring in the closed loop. Such effect has been demonstrated in the interactions between variations of the heart beat and systolic blood pressure in adults .
To determine the significant interactions between HRV and respiration, a threshold level of coherence has been generally set either arbitrarily [25–27] or based on statistical criteria derived from the sampling distribution . Any value of coherence above the threshold is considered as significant. Methods that are able to compute more appropriate coherence thresholds, theoretically or experimentally derived by the knowledge of estimator and signals under investigation, can avert discretionary use of a threshold not based on theoretical model or empirical approach. In preterm infants, as cardio-respiratory interactions are weaker (or even absent) than adults, it is even more important to establish a solid statistical criterion to assess a reliable significance threshold for the coherence function.
To address all of these issues, we propose a new framework for studying the interaction between HRV and respiration in preterm infants. We present a bivariate model specifically tailored for the preterm infants in which auto- and cross- spectral density functions of the RR and respiration are evaluated using autoregressive model parameters. Using the procedure , we estimated causal coherence and causal gain at different frequency bands. We then determined the significance of the coherence by application of surrogate data analyses .
The combined methodologies, bivariate modeling and surrogate data analysis, were applied on data collected from 11 preterm infants to test whether our methods can detect differences in RSA during different known physiological perturbations, and whether characteristic dynamic changes occur as a consequence of such perturbations.
The data employed in our analysis are from an experimental protocol (see details below) designed to study the effect of mild vibro-tactile stimulation in preterm infants. ECG signals were recorded at a sampling rate of 200 Hz using surface electrodes placed below the left clavicle and at right abdomen. Respiratory inductance plethysmography (Somnostar PT; Viasys Healthcare, Yorbalinda, CA) was used at a sampling rate of 100 Hz to record abdominal respiratory movements. The recorded ECG and respiratory signals were displayed in real-time during the studies and stored on hard disk using a PC compatible A–D data acquisition system (Embla S700, Broomfield, CO). Further details of the recording procedure can be found in .
The R-wave peaks were detected using a derivative and threshold algorithm. The peak to peak (RR) series were visually inspected and corrected for any artifacts or erroneously detected peaks. To improve the detection, we applied a fourth order band-pass zero-phase Butterworth filter. RR series, together with the corresponding respiratory signals obtained were interpolated and re-sampled at 3Hz.
As a preliminary analysis to understand the variation in heart rate, we computed the following time domain measures using the standard approach : 1) The standard deviation of all RR intervals (SDNN). 2) The standard deviation of the average of all RR intervals in all the 5-minute segments (SDANN). 3) The square root of the mean of the sum of the squares of differences between adjacent RR intervals (RMSSD). 4) The mean of the standard deviation of all RR intervals for all 5-minute segments (SDNNi). 5) The standard deviation of differences between adjacent RR intervals (SDSD).
As frequency domain indices have been established for adults by defining specific frequency ranges of interest, a critical point in this study was to establish a similar classification of frequency ranges for infants. Starting from the observation that our preterm infants have a predominant breathing frequency ~1Hz, we introduce a more refined characterization of the respiratory range, thereby classifying four different frequency ranges: Low Frequency (LF: 0.01–0.15Hz) and three High Frequency ranges (HF1: 0.15–0.45Hz, HF2: 0.45–0.7 Hz and HF3: 0.7–1.5 Hz). HF3 is generally the frequency range corresponding to the eupneic respiratory rhythm of preterm infants.
To obtain the frequency domain measures, both the RR and respiration are considered as output variables of a multivariate autoregressive model . The coefficients of the model are determined by solving the extended Yule-Walker equations  and autospectra, coherence and gain are derived in the frequency domain from these coefficients. The statistical significance of the coherence for each infant is then determined by surrogate data analysis . Altogether, this approach provides a quantification of the linear relationship between RR and respiration, as well as its significance, defined along the entire range of frequencies. Consequently, specific indices for each of the frequency bands defined above can be derived from the spectral estimates.
A bivariate autoregressive model is employed to study the interaction between RR and respiration (RP). The model is defined as
where M is the order and is set at 32, N is the total number of data points,
w(n) represents the white noise and aij(k) represents the autoregressive coefficients. We used a recursive algorithm to determine the coefficients of the autoregressive model; spectral components are determined from these coefficients . Thus in the frequency domain the model is represented as,
where with i,j =1,2 and , a complex quantity
The Eqn (2) can be reformulated as
The coherence γ2 at a specific frequency f is evaluated using the classical definition as
where PRR(f) and PRP(f) are the auto-spectral density functions of RR and RP respectively. PCROSS(f) is the cross spectral density between RR and RP. The auto-and cross-spectral density functions are evaluated as
The causal coherence is calculated using the same Eqn (3) with the corresponding loop set to zero, thus for the respiration to RR causal coherence (RP → RR), we set h21 = 0 and for RR to respiration causal coherence (RR → RP) h21 = 0. Similarly, corresponding gains are calculated as
It has to be noted that the estimated causal gains (Eq 6) are relevant only if the corresponding causal coherence values are significant. Hence detecting the significant coherence values between RR and respiration is an important step in establishing the presence of interactions in preterm infants.
As both coherence and gain are estimated along the entire frequency range up to the Nyquist frequency, we can further compute the maximum coherence and the corresponding gain in each of the defined frequency band (LF, HF1, HF2 or HF3).
Using this approach for each original signal, a set of surrogate series is generated that retain properties with respect to the spectral power of the original signal but are uncoupled with respect to the other original signal. In our case these two signals are RR and RP. Thus for each segment of RR (or RP) signal we generated an ensemble of 100 pairs of Fourier Transform (FT) surrogate time series.
The FT surrogates were constructed by computing the Fourier Transform of the RR (RP) and then randomizing the phases, while keeping the magnitude unchanged. By this procedure the surrogate series of RR (RP) has the same magnitude for the power spectrum as the original RR (RP) but are uncoupled with respect to RP (RR). We set the threshold for zero coherence at the 95th percentile of the coherence sampling distribution with the condition that any coherence value above this threshold is significant. Thus this approach allows determining the significance of coherence spectra on an individual subject basis.
We illustrate our methods for the analysis of RR interval and respiration time series from data collected as part of a comprehensive study designed to investigate the effect of vibro-tactile stimulation on the stability of breathing in preterm infants . The study was conducted at the Newborn Intensive Care Unit, University of Massachusetts Memorial Healthcare and approved by the Committee for the Protection of Human Subjects in Research at the University of Massachusetts Medical School. Results from similar recordings using a different sample of subjects have previously shown that mild vibro-tactile stimulation applied continuously for 10-minutes to the preterm infants improves the stability of eupneic breathing . The data for the analysis in the present paper were obtained from a different group of infants in which stimulation occurred only when there was an apneic pause in breathing.
Infants included in this analysis were: 1) preterm with gestational age (GA) <36 wks and post-conceptional age (PCA) >30 wks at time of study; 2) spontaneously breathing room air or receiving supplemental oxygen through nasal cannulae at a fixed flow rate. As part of the comprehensive study, each of the 11 infants participated in an experimental protocol that evaluated the effects of apnea-triggered stimulation. A mattress vibrotactile stimulation was actuated only when the infant stopped breathing for >5 sec and was terminated automatically upon spontaneous resumption of breathing or via clinical intervention if resulting apnea was >20sec. Apnea- triggered stimulation was controlled by a customized program (Labview v.6.1, National Instruments, Austin TX) using an apnea detection algorithm. Apnea-triggered stimulation was implemented for 10-minute periods (On periods), alternated with 10-minute intervals during which no stimulation was given (Off periods).
The ECG and respiratory data considered in our analysis for each of the 11 infants consists of two 10 minute selected control periods, two 10 minute periods that immediately followed the control periods when the mattress was in Off condition, and two 10- minute segments when the mattress was in On condition. This allowed a total of 20 minute segments for each of the control, Off and On conditions.
We considered two 10 minute segments from 11 preterm infants for three conditions: 1) control, 2) Off and 3) On mattress vibration: a total of 66 segments. For each segment we computed the time and frequency domain indices, and each index is then averaged among subjects for each of the three conditions. Figure 1 shows an example of the intervention protocol with two 10 minute segments of Off and On conditions.
Figure 2 reports the group averages of the time domain measures during control, Off and On conditions for the 11 preterm infants. Three of the measures, RMSSD, SDNNi and SDSD appear to be reduced during intervention indicating a reduction in the total heart rate variability. However separate ANOVA revealed no main effect of intervention condition for any of the time domain indices; however, RMSSD approached significance [F(2,30)=3.062, p=0.062).
In this section we report the frequency domain measures computed using bivariate algorithm as well as significant levels of coherence determined using surrogate data analysis.
After preliminary optimal order analysis (Akaike Information Criterion) using data from control conditions, an autoregressive order of 32 was fixed for uniformity for all selected 10 minute segments using a “maximum denominator” criterion . In particular, we tested our model with different orders ranging from 8 to 36 in steps of four and found that for those infants with irregular breathing, the model was generally higher, with orders approaching the upper bound of 32 more often than those with regular breathing. In fact, this high order (as compared to orders commonly used by applications to adult populations) is justified both by the infant’s faster heart beat  and by the use of interpolated 3 Hz time series, which together increase the Nyquist frequency up to at least 1.5 Hz. We estimated the bivariate autoregressive parameters as well as the noise variances, to compute the coherence, the RR and RP power spectral densities, and the causal gains (Eq. 6). The maximum coherence as well as the gain at the frequency corresponding to maximum coherence is calculated in each of the defined frequency bands.
Figure 3 presents the frequency domain indices of a 200s segment during control conditions from an infant. The infant (subject 4) exhibited a regular breathing pattern during the 200s segment as shown in (A), confirmed by a respiratory power spectrum concentrated in HF3. Furthermore, a considerable amount of power density is observed in the HF3 range of the RR spectrum (as shown in (B)). The causal gains are also presented. The classical coherence value shows a peak value in the HF3 range which is the outcome of the causal interaction of respiration to RR which has been demonstrated as a peak in the HF3 in the causal coherence (RP → RR). Such a peak is absent in the other direction as evident from the causal coherence (RR → RP).
Figure 4 represents the power spectral density measures during control, Off and On conditions calculated for each of the 10 minute segment for all the subjects. Due to irregularity in breathing, we can see large power in the LF range of respiration. Figure 5 represents the cardio-respiratory interactions using classical coherence and causal coherences (both RP → RR and RR → RP). We found that the coherence values during Off and On conditions were greater than the control periods.
Using repeated measure ANOVA with three factors, (i) conditions: control, Off and On, (ii) frequency bands: LF, HF1, HF2 and HF3, (iii) measures: classical and causal, we find a main effect for frequency bands (F(3,30)=7.134, p=0.004; Greenhouse Geiser epsilon=0.678) and for causality (F(1,10)=229.65, p<0.001, Greenhouse Geiser epsilon=1.0). However, no main effect is observed for conditions (p= 0.094). Post hoc analysis to determine which frequencies were significantly different from each other (with Bonferroni correction for multiple comparisons), independent of conditions and measures, revealed HF3 significantly greater than HF2 (p=0.003) and HF2 smaller than LF (p=0.028). Independent of frequency bands or conditions, classical coherence (estimated marginal mean = 0.184, SE=0.18) was significantly larger than causal coherence (estimated marginal mean =0.094, SE =0.009). As the coherence values are very low (<0.5), the significance of the coherence was further assessed using surrogate data analysis.
Figure 6 shows the coherence measures along with the threshold for zero coherence function between 0.01Hz and 1.5 Hz for a 10 minute segment of control conditions from the two infants (subject 8 and subject 9). Focusing on the HF3 range, it can be seen that for subject 9 (shown in (A)) classical coherence as well as causal coherence (RP → RR) values are significant based on surrogate data analysis, as the maximum value of coherence in this range is greater than the threshold for zero coherence in the same range. On the other hand, RR → RP the causal coherence is less than the threshold for zero coherence. For subject 8, however the classical coherence and RR → RP causal coherence are significant. These observations suggest that for subject 9 the interaction from RP to RR predominates, whereas in subject 8 the heart beat variations present a more important anticipatory phase with respect to respiration.
Additional results from the surrogate data analysis are compiled and presented in Table 1. If the coherence values are above the threshold for zero coherence in the respective frequency bands, it is marked as “Y”, and as “N” otherwise. The classical coherence values are represented in the column marked “I” whereas both RR → RP and RP → RR causal coherence are presented in the columns marked “II” and “III” respectively. The number of subjects with significant coherence values in each of the frequency bands as well as in each measure is shown below each sub-panels of the table.
Since causal gains are significant only if their corresponding causal coherence values are significant, we considered only those significant values as represented in Figure 7. Any non-significant values are excluded.
We established the presence of RSA in preterm infants by detecting significant value of coherence using a new methodological framework for the analysis of cardiovascular variables. Our approach is adapted from a set of mathematical tools that have been successful in quantifying important cardiovascular control mechanisms in adult humans, to specifically reflect the physiology of the developing cardiovascular system.
In contrast to the traditional spectral techniques , our bivariate approach considers both the RR and respiratory signal in a single modeling framework, thereby providing frequency measures in a computationally efficient way. The bivariate algorithm is able to separate the mutually directional causal interactions between RR and respiration as opposed to traditional spectral methods whereby the relationship between RR and respiration is assumed as unidirectional. Our approach also provides a more realistic linear model of RR and respiration which enables to quantify the relationship with a transfer function approach and, at the same time, to model the heart rate variations as composed of a discrete number of oscillatory dynamics that can be more clearly associated to other interconnected periodic physiological rhythms, such as the respiratory cycle. Furthermore, the gain is determined at specific oscillations where coherence indicates maximum coupling, with significant levels of coherence established using surrogate data analysis. Thus, the bivariate model along with the surrogate data analysis provides a meaningful statistical framework for establishing a quantitative assessment of cardio-respiratory coupling in each preterm infant.
To explore the variation of HRV in preterm infants, we used classical time domain analysis and compared to vibro-tactile stimulation in three conditions. Classical time domain indices did not detect any variation in these indices among different conditions.
We then evaluated the frequency domain measures by performing a bivariate analysis of RR and respiration using multivariate autoregressive model. We found that the coherence during the respiratory frequency range (HF3) is significantly different from other high frequency ranges. These findings suggest presence of RSA in preterm infants. This outcome is important because it was achieved without limiting the analysis to segments to only regular breathing patterns; frequency ranges that also included pauses in breathing as well as irregular breathing patterns are also incorporated, without resorting to filtering techniques which might bias the identification. Such variability in breathing is common in preterm infants.
We performed surrogate data analysis to determine the significance of coherence value in each infant. Using this approach we found that only 5 out of 11 infants (< 50%) have significant coherence during control conditions in the HF3 range suggesting that the interactions in the remaining infants are indistinguishable from fictitious interactions. However during vibro-tactile stimulation (Table 1), we found that 10 out of 11 (~90%) infants have significant coherence values at HF3.
The closed-loop causal formulation quantifies features otherwise not identifiable by standard measures of heart rate variability, and allows for separation of the interactions between RR and respiration occurring in both directions, as evident from the significant causal coherence values. In particular we found that the number of subjects with significant coherence value from RR to respiration during HF3 is higher during intervention conditions. Considering the increase in significant number of causal gain values (Figure 7) during intervention, we conclude that the mild vibro-tactile stimulation might enhance or establish cardio-respiratory coupling in preterm infants. This observation may reflect the stabilizing effect of stochastic mechanosensory stimulation on eupneic breathing . The mechanism of this putative enhancement in cardio-respiratory coupling is uncertain in preterm infants. Stretch sensitive somatic or visceral mechanoreceptors could enhance the activity of cardio-respiratory activities at multiple central integrative sites .
Among the different frequency bands, we found large percentage of power in LF range. It has been pointed out that a pause in breathing introduces frequency in the low frequency ranges  and since most of our preterm infants have pauses greater than 5 seconds, these apneas are indeed reflected in the LF. Of note, in some infants we also observed some genuine slow frequency patterns in the respiratory signal, further contributing to LF powers. The mechanism responsible for LF fluctuations in infant respiration is uncertain. A number of factors have been proposed to contribute to respiratory instabilities in neonates [33–35], which might contribute to LF fluctuations in respiration. Proposed mechanisms include an increased loop gain of chemical feedback , the modulatory effects on respiration of sleep , and irregularities related to intrinsic nonlinear properties of the central respiratory oscillator [39–41]. Our finding of a bidirectional interaction between HRV and respiration raises the interesting possibility that fluctuations in cardiac dynamics might be an important factor for causing pathological fluctuations in breathing, including apnea episodes.
Although it might appear that these low frequency oscillations can be eliminated by a filter to improve the computational efficiency of the model, the core physiological mechanism related to apnea can be adversely concealed. Hence we included all the low frequency components observed in the data for our analysis and estimated the original power distribution at the expense of only a marginal increase in complexity for the bivariate model. Our choice of the highest optimal order model as indicated among the considered population was cautionary, in the sense that we wanted to achieve the highest frequency resolution for an accurate study of the distribution among the newly defined frequency ranges. Of course, use of lower orders would be admissible on a case by case basis, but always with supporting Akaike results. Of note, a further undersampling of the signal would also decrease the optimal order; however one has to make sure that all power is within the Nyquist frequency.
The advantage of our approach is that the parametric linear regression defined by our bivariate modeling framework can compute frequency domain measures in a parsimonious way, and allows for modeling different oscillations underlying cardiovascular control. In addition, the bivariate analysis enables the determination of frequency domain measures at predefined frequency ranges, such as the respiratory frequency ranges. Surrogate data analysis provides a powerful method for computing coherence thresholds using empirical distributions derived from each infant. The threshold for coherence determined at a specific percentile of the empirical sampling distribution allows us for a highly sensitive and specific statistical quantification, critical for reliable detection of weak interactions in a single infant.
The limitation of our approach is that we employed a bivariate model that is suitable for stationary or quasi stationary data sets, thus assuming that each of the epochs that we analyzed satisfies such criteria. However there are non-stationary events in these segments particularly associated with apnea that may affect the stationary nature of the respiratory signals. In an immature system, both the cardiovascular and respiratory mechanisms may exhibit transient interactions, but due to system vulnerability, the interactions may be fast changing and non persistent. To understand such changes we need a model which estimates the measures in a time varying manner.
Another limitation is that we derived a continuous representation of RR using an interpolation method. The low HRV in these preterm infants require a fine time resolution, and a simple interpolation may not be accurately representing the real dynamic process associated with the RR variations. A point process model, possibly including respiration as one of the covariates  may address both limitations as well as provide a more accurate determination of time varying representation of cardio-respiratory interactions (work is in progress).
The data employed in our analysis is obtained from an experimental protocol designed to study the effect of vibro-tactile stimulation on preterm infants. In many species, maternal licking of newborn at birth is considered as an early somatosensory stimulation to affect breathing . Although our analysis suggests the presence of RSA in preterm infants and identifies changes in the coupling due to vibro-tactile stimulation, it is not known whether such mild interventions improve breathing in these preterm infants. It has been shown that the vibro-tactile stimulation applied continuously for a specific duration of time improves breathing in preterm infants , however we have not studied the effect of intermittent stimulation on breathing patterns of these preterm infants which is beyond the scope of this work and will be presented elsewhere as our investigation further progresses.
This work was supported by CIMIT under U.S. Army Medical Research Acquisition Activity Cooperative Agreement W81XWH-07-2-0011 (R.B.), by National Institutes of Health (NIH) Grants R01-HL084502 (R.B), R01-DA015644 (E.N.B), DP1-OD003646 (E.N.B), R01-HL49848 (D.P.), R01-HL071884 (D.P.), by an American Heart Association Scientist Development Grant (E.B.S.), intramural Pilot Project Programs of the University of Massachusetts Medical School and Mental Retardation Developmental Disabilities Research Center (D.P.) and the Hansjörg Wyss Institute for Biologically Inspired Engineering at Harvard University (D.P.). The information contained herein does not necessarily reflect the position or policy of the Government, and no official endorsement should be inferred.
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.