Classification of patients undergoing weaning from mechanical ventilation using the coherence between heart rate variability and respiratory flow signal

Mechanical ventilation remains the main medical treatment for acute respiratory failure, and is one of the defining interventions of intensive care medicine as a specialty. Around 50% of patients in intensive care units require mechanical ventilation (Tobin 2004, Boles et al 2007). Weaning is the process of transferring the work of breathing from the ventilator back to the patient, and is still one of the most challenging problems in clinical practice (Blackwood et al 2011, Burns et al 2013). Withdrawal of mechanical ventilation should be performed as soon as patients can breathe on their own. Both unnecessary delay and premature weaning may have adverse effects on patient outcome, prolonging mechanical ventilation and length of stay in the intensive care unit (MacIntyre 2004, Frutos-Vivar et al 2006).

Several studies have been performed to identify physiological variables that are useful indicators for undertaking a weaning trial (Yang and Tobin 1991, Meade et al 2001, Casaseca-de-la Higuera et al 2006, Stawicki 2007, Santos Lima 2013). McConville et al recently published a review of strategies to reduce the duration of mechanical ventilation, and a list of risk factors for unsuccessful discontinuation (McConville and Kress 2012). They found that most patients who receive mechanical ventilation have acute respiratory failure in the postoperative period, pneumonia, congestive heart failure, sepsis, trauma, or acute respiratory distress syndrome. Esteban et al explored whether mortality in mechanically ventilated patients has changed over time (Esteban et al 2013). Their results show that the proportion of patients successfully liberated from mechanical ventilation at the first attempt has increased over time, while the rate of reintubation after scheduled extubation has remained similar.

Tobin, in an editorial (Tobin 2012), noted that the removal of endotracheal tubes is welcomed by critically ill patients, not just because the tube itself is uncomfortable but also because it is a barrier to communication and extubation indicates that they are on the road to recovery. Unfortunately, as many as 20% of patients require reintubation shortly after the extubation process. Further, in a small proportion of these patients, that require rapid reintubation, this process is the start of a lethal chain of events. While there are several approaches to weaning, it remains unclear what is the optimal strategy.

Analysis of the heart rate variability (HRV) signal is a non-invasive tool widely used to assess autonomic nervous system activity (Task Force 1996). HRV analysis is based on the study of temporal oscillations between heartbeats. It is well known that respiratory activity affects cardiovascular regulation (Orini et al 2012). Respiratory sinus arrhythmia is an index of cardiac vagal excitation, and it can be considered to be a phenomenon that results directly from the interaction between the cardiovascular and respiratory systems (Grossman and Taylor 2007, Ritz 2009). Cardiorespiratory interdependence during weaning trials does seem to be related to specific aspects of the coordination of dynamic autonomic function (Pinsky 2000, Caminal et al 2010).

Various studies have been performed in order to determine parameters related to cardiorespiratory interaction and their influence on the weaning process (Shen et al 2003, Garde et al 2010, Langley et al 2010, Huang et al 2014). To date, however, it is not clear whether the functional relations between breaths and heart beats are more stable in patients with successful trials. On the other hand, spectral analysis of heart rate, respiration and blood pressure signals is a non-invasive approach that is widely used to investigate cardiovascular and cardiorespiratory control mechanisms (Task Force 1996, Meste et al 2005, Mainardi 2009).

Our group have previously assessed the characteristics of patients undergoing weaning using spectral parameters and mutual spectral analysis of the RR interval and respiratory flow signal (Arcentales et al 2010). In other studies (Caminal et al 2010, Garde et al 2010), we found that using parameters extracted from joint symbolic dynamics applied to time series of heart rate and respiratory frequency, and a support vector machine classifier was a suitable approach for describing cardiorespiratory interactions of patients during the weaning process.

The aim of this study was to assess the contribution of the spectral components of the HRV and respiratory flow signals and their coherence to identifying patients with successful spontaneous breathing trials, unsuccessful trials, or initially successfully trials but who are unable to maintain spontaneous breathing and require the reinstitution of mechanical ventilation within 48 h. We considered two methods for classification these patients, applying principal component analysis (PCA) to the spectral signals (S-Class), and sequential floating feature selection (SFFS) to extracted features (F-Class). A fuzzy K-nearest neighbour classifier was used to determine the best spectral signals and parameters. At the outset, we hypothesized that spectral coherence of the signals would be able to improve the classification of patients due to undergo weaning, in particular, distinguishing between those in whom trials would be successful and those who would require reintubation.

Electrocardiogram (ECG - lead II) and respiratory flow (FLW) signals were measured in 121 patients undergoing weaning from mechanical ventilation (74 males, aged 63.9  ±  17 years; 47 females, aged 66.3  ±  15 years). These patients were recorded between January 2003 and April 2006, in the Intensive Care Department at the Santa Creu i Sant Pau Hospital, Barcelona, Spain and the Getafe Hospital, Getafe, Spain. All subjects were studied according to a protocol approved by the local ethics committee. ECG signals were acquired using a SpaceLabs Medical monitor (now Spacelabs Healthcare, Snoqualmie, WA, USA). Respiratory flow signals were recorded using a pneumotachograph, connected to an endotracheal tube, and consisting of a Datex-Ohmeda monitor (GE Healthcare, Milwaukee, WI, USA) and a Validyne MP45-1-871 Variable-Reluctance Transducer (Validyne Corp., Northridge, CA, USA). The signals were recorded for 30 min, at a 250 Hz sampling rate and 12 bit resolution.

Using clinical criteria based on a T-tube trial, patients were included randomly in this study, according to standard indices used to assess readiness for a spontaneous breathing trial (SBT): treatment of acute respiratory failure (with inspired oxygen fraction [FiO₂]  ⩽  0.4, oxygen saturation [SO₂]  ⩾  90% and the need for positive end-expiratory pressure [PEEP  ⩽  5 cm to H₂O]), hemodynamic stability (absence of myocardial ischemia and/or heart failure, cardiac frequency  ⩽  140 bpm, and stable arterial blood pressure with tolerance to a reduction in inotropic support), and adequate respiratory muscle function (acceptable respiratory rate,⩽  30–35 breaths per minute). The patients were ventilated in assist-control (AC) mode or pressure support (PS) mode.

The patients undergoing an SBT were disconnected from the ventilator just at the moment of the test, and left to breathe through an endotracheal tube for 30 min with close monitoring. Any patients who were not able to breathe spontaneously were reconnected, while those who were able to maintain spontaneous breathing were extubated. When a patient was still able to maintain spontaneous breathing after 48 h, the weaning trial process was considered successful. If not, the patient was reintubated.

According to the outcome of the SBT, the patients were classified into three groups: GSucc, 73 patients with successful weaning trials who could maintain spontaneous breathing after 48 h; GFail, 33 patients who failed to maintain spontaneous breathing and were reconnected to mechanical ventilation after 30 min; and GRein, 15 patients who successfully passed the weaning trial, but had to be reintubated within 48 h (figure 1). The baseline characteristics of the patients such as age, sex, respiratory rate, ventilation mode, days of mechanical ventilation and the main diagnosis are summarized in table 1.

Zoom In Zoom Out Reset image size

Both the ECG and the respiratory flow signals were preprocessed to reduce artifacts and remove linear trends. Time series of the cardiac interbeat interval were extracted automatically from the ECG signal using an algorithm based on wavelet analysis (Martínez et al 2004). Ectopic beats were identified and removed using an algorithm based on local variance estimation (BailÓn et al 2006). The HRV signal derived from the interbeat interval was resampled at 5 Hz using a cubic spline function. The respiratory flow signal was decimated at the same frequency as the HRV signal, and both signals were filtered and synchronized.

The study of the HRV and the respiratory flow signals and their interaction involves the following processing steps: spectral estimation and calculation of the coherence between these signals, feature extraction from them, dimensionality reduction using PCA and SFFS, and classification through a fuzzy K-nearest neighbour method (figure 2).

Zoom In Zoom Out Reset image size

We propose classify the groups of patients (GSucc, GFail and GRein) on weaning process, through the principal components analysis applied to the spectral signals of HRV and FLW and the magnitude squared coherence between them. This classification is compared with the one made using the features extracted of these signals.

3.1. Spectral estimation

Power spectral densities (PSDs) of HRV and the respiratory flow signals were estimated at frequency f applying Welch's averaged modified periodogram method (Welch 1967, Proakis and Manolakis 2008), defined by

$$\begin{eqnarray}&&{{S}_{xx}}(f)=\frac{1}{LMU}\underset{i=0}{\overset{L-1}{\sum}} {\left| \underset{n=0}{\overset{M-1}{\sum}} x\left(n+\text{i}D\right)w(n)\left({{\text{e}}^{-j2\pi fn}}\right)\right|}^{2}\end{eqnarray} \tag{ 1 }$$

where x(n) denotes the autocorrelation function of each signal, L is the number of segments, M the length of the segments, D the overlap between segments, and w(n) is a Hamming window. Here, U is a normalization factor related to the characteristic of the window function w(n), which removes the energy bias introduced by the windowing,

$$\begin{eqnarray}&&U=\frac{1}{M}\underset{n=0}{\overset{M-1}{\sum}} {\left\vert w(n)\right\vert}^{2}.\end{eqnarray} \tag{ 2 }$$

Cross-power spectral density was also estimated by applying Welch's periodogram to the cross-correlation function between the HRV and respiratory flow signals. Its normalization is known as the coherence function, and is given by

$$\begin{eqnarray}&&{{\Gamma}_{xy}}(f)=\frac{{{S}_{xy}}(f)}{\sqrt{{{S}_{xx}}(f)}\sqrt{{{S}_{yy}}(f)}}\end{eqnarray} \tag{ 3 }$$

where x and y are the signals, and S_xy(f) the cross-spectral density at frequency f between these signals. The magnitude squared coherence (MSC) is defined by Clifford Carter et al (1973)

$$\begin{eqnarray}&&\mid {{\Gamma}_{xy}}(f){{\mid}^{2}}=\frac{\mid {{S}_{xy}}(f){{\mid}^{2}}}{{{S}_{xx}}(f){{S}_{yy}}(f)}\end{eqnarray} \tag{ 4 }$$

and normalized such that

$$\begin{eqnarray}&&0\leqslant \mid {{\Gamma}_{xy}}(f){{\mid}^{2}}\leqslant 1.\end{eqnarray} \tag{ 5 }$$

Different values of L and M were tested to maximize the frequency resolution and reduce the noise in the spectral estimation. The best results were obtained with M  =  2 min, an overlap of 50% (D  =  M/2).

Since the spectral estimation method yields non-zero values of the magnitude squared coherence where there is no spectral coherence (Faes et al 2004), a threshold for real coherence was calculated using two 30 min white noise signals. These signals were analyzed in the range of 0 to 0.8 Hz. For each iteration, the maximum coherence value was recorded. A threshold of 0.25 was set for the average value of the maximum coherence and the magnitude square coherence that did not reach this threshold was considered null.

3.2. Feature extraction

Commonly, the power spectral density of the HRV is analyzed considering the following spectral bands (Task Force 1996): very low frequency (VLF: 0–0.04 Hz), low frequency (LF: 0.04–0.15 Hz), and high frequency (HF: 0.15–0.4 Hz). In addition, as signals from respiratory activity sometimes extended to frequencies beyond 0.4 Hz, we considered the HF band up 0.8 Hz (HF: 0.15–0.8 Hz). Therefore, the HRV and MSC signals were obtained for these frequency bands, and the FLW signal was obtained in the range of 0 to 0.8 Hz.

The power spectral densities of the HRV, FLW, and MSC signals were characterized by the following frequency domain parameters: power peak (P_p), frequency peak (f _p), central requency (f _c), total power (P_T) and the area under spectral signal (P_a). For HRV and MSC signals, P_a was estimated for each frequency band (VLF, LF and HF), as well as the power ratio between VLF and HF (R_VLF/HF), LF and HF (R_LF/HF), VLF and P_T $\left({{R}_{\text{VLF}/{{P}_{T}}}}\right)$, LF and P_T $\left({{R}_{\text{LF}/{{P}_{T}}}}\right)$, HF and P_T $\left({{R}_{\text{HF}/{{P}_{T}}}}\right)$.

Additionally, the PSD of the FLW signal was characterized through the discriminant band (DB) defined by the power between the frequencies corresponding to 10% of the value of P_p (right and left of f _p). For this parameter (DB), P_a and the ratio between DB power and the total power of the signal $\left({{R}_{\text{DB}/{{P}_{T}}}}\right)$ were also estimated.

Summarizing, a total of 33 features were extracted from the signals of each patient: 12 parameters for HRV and MSC signals, respectively, and 9 for the respiratory flow signal (see table 2). An example of the spectral analysis for each signal and the extracted parameters is illustrated in figure 3.

Zoom In Zoom Out Reset image size

Table 2. Frequency domain parameters extracted from the HRV, FLW and MSC signals.

Feature	Description
HRV and MSC
Pp [n.u.]	Peak power
f p [Hz]	Peak frequency
f c [Hz]	Central frequency
PT [n.u.]	Total power (0 to 0.8 Hz)
Pa [n.u.]	Area under spectral signal for each frequency band (VLF, LF, HF)
RVLF/HF	Power ratio between VLF and HF
RLF/HF	Power ratio between LF and HF
${{R}_{\text{VLF}/{{P}_{T}}}}$	Power ratio between VLF and total power
${{R}_{\text{LF}/{{P}_{T}}}}$	Power ratio between LF and total power
${{R}_{\text{HF}/{{P}_{T}}}}$	Power ratio between HF and total power
FLW
Pp [n.u.]	Peak power
f p [Hz]	Peak frequency
f c [Hz]	Central frequency
PT [n.u.]	Total power (0 to 0.8 Hz)
DB [n.u.]	Discriminant band (power between right and left of f p
f l [Hz]	Frequency at the left of DB
f r [Hz]	Frequency at the right of DB
Pa [n.u.]	Area under the curve of PSD (DB)
${{R}_{\text{DB}/{{P}_{T}}}}$	Power ratio between DB and total power

n.u.: normalized units. The signals were normalized to the total power.

3.3. Dimensionality reduction

Principal component analysis (PCA). A data matrix X contains in each row p_i the spectral signal of each patient

$$\begin{eqnarray}&&\mathbf{X}={{\left[{{\mathbf{S}}_{{{p}_{1}}}}\;\;\;{{\mathbf{S}}_{{{p}_{2}}}}\;\;\ldots \;\;{{\mathbf{S}}_{{{p}_{np}}}}\right]}^{T}}\end{eqnarray} \tag{ 6 }$$

were ${{S}_{{{p}_{i}}}}$ contains the values of the frequency range between 0 to 0.8 Hz.

According to singular value decomposition theory (Jolliffe 2002), X can be written as

$$\begin{eqnarray}&&\mathbf{X}=\mathbf{UL}{{\mathbf{A}}^{\prime}}\end{eqnarray} \tag{ 7 }$$

where U and A are the orthonormal columns, and L is a diagonal matrix. The principal component (PC) scores are given by

$$\begin{eqnarray}&&\mathbf{Z}=\mathbf{UL}.\end{eqnarray} \tag{ 8 }$$

Therefore, the spectral signals are re-defined in the new dimensional space by the PCs. The variance captured by each PC is defined by the eigenvalues on the diagonal of L (Jolliffe 2002).

Sequential floating feature selection (SFFS) is a feature selection method to obtain the optimal characteristics and optimize the classification rate (Pudil and Novovi 1994). The criterion for the feature selection was maximizing the percentage of the patients belonging to each class with an optimal accuracy.

3.4. Classification

A fuzzy K-nearest neighbour method was used to classify patients into each group. Specifically, the degree of membership of each patient to its class was determined based on the Euclidean distance function for each K-nearest neighbour.

Let W={x₁, x₂, ..., x_n} be a training set formed by n labeled samples. The degree of membership for a validation sample to the each class is calculated based on a subset of W formed by the K-nearest neighbour to these sample set, and is defined by Keller et al (1985):

$$\begin{eqnarray}&&{{u}_{i}}(x)=\frac{{\sum}_{j=1}^{k}{{u}_{ij}}\left(\frac{1}{\|x-{{x}_{j}}{{\|}^{2/(m-1)}}}\right)}{{\sum}_{j=1}^{k}\left(\frac{1}{\|x-{{x}_{j}}{{\|}^{2/(m-1)}}}\right)}\end{eqnarray} \tag{ 9 }$$

were u_ij represents membership in the ith class of the jth labeled sample. The m variable determines the contribution of each neighbour to the membership value according to its distance. With m increasing, the neighbours are more evenly weighted, and the effect of the relative distance is lower. When m approaches to one, the closer neighbors are weighted far more heavily than those farther away. The optimal value selected was m  =  3.

The minimum membership percentage for each class was considered in 65%. Samples with values below this percentage were considered misclassified.

3.5. Statistical analysis

In order to evaluate the results of this work, we decided to classify the weaning candidates considering the following comparisons between patients:

• Successful versus failed patients (GSucc versus GFail).
• Successful versus reintubated patients (GSucc versus GRein).
• Failed versus reintubated patients (GFail versus GRein).

Normally, GSucc group is bigger than GFail and foremost GRein groups, and consequently the sample set is unbalanced. Therefore, a random undersampling process was applied to select the training set (He and Garcia 2009). In order to avoid any bias due to the random selection, the same process was evaluated 50 times for each validation patient. The final percentage of membership is the average value of these 50 iterations.

4.1. Performance of S-Class

4.2. Performance of F-Class

Spectral signals were characterized through the features described in section 3.2. Applying the sequential floating feature selection method the best parameters to classify the groups of patients were obtained. Table 6 presents mean, standard deviation and p-value of the most statistically significant parameters.

Table 6. Mean and standard deviation of the most statistically significant parameters when comparing GSucc, GFail and GRein patients, using the sequential floating feature selection method.

GSucc		GRein	GFail	p-value GSucc versus GFail	p-value GSucc versus GRein	p-value GFail versus GRein
Feature PSD HRV
Pp	20.59  ±  15.56	15.59  ±  11.17	25.22  ±  20.94	n.s.	n.s.	0.04
Feature PSD FLW
f p	0.38  ±  0.10	0.52  ±  0.13	0.40  ±  0.14	<  0.0001	n.s.	0.05
${{R}_{\text{DB}/{{P}_{T}}}}$	0.73  ±  0.19	0.91  ±  0.33	0.76  ±  0.12	<  0.0001	n.s.	0.01
Feature MSC
f p	0.35  ±  0.18	0.43  ±  0.23	0.25  ±  0.25	0.02	n.s.	0.05
f c	0.37  ±  0.09	0.44  ±  0.10	0.38  ±  0.10	0.004	n.s.	n.s.
${{R}_{\text{HF}/{{P}_{T}}}}$	0.77  ±  0.32	0.80  ±  0.28	0.66  ±  0.32	n.s.	0.04	0.05

Features of PSD of HRV signal: the mean value of P_p is lower in GFail than in GRein group of patients (p  =  0.04). In general, P_p value is shorter in GFail patients than in the other groups.

Features of PSD of FLW signal: the mean value of f _p and ${{R}_{\text{DB}/{{P}_{T}}}}$ are lower in GSucc than in GRein and GFail, respectively. The best discriminant values were obtained when comparing GSucc versus GFail (p  <  0.0001). When comparing GFail versus GRein, the best result was obtained with the mean value of ${{R}_{\text{DB}/{{P}_{T}}}}$ (p  =  0.01).

Features of MSC signal: the results show that f _p(p  =  0.02) and f _c(p  =  0.004) are shorter in GSucc than in GFail patients. Furthermore, it is shown that ${{R}_{\text{HF}/{{P}_{T}}}}$ is lower in GRein than in GSucc (p  =  0.04) and GFail (p  =  0.05) patients.

The accuracy, sensitivity and specificity of these classifications show that the best results is obtained with the features extracted from FLW signal when comparing GS versus GF an accuracy of 74%.

In this study, we have explored the classification of patients who are candidates for weaning from mechanical ventilation using the spectral components of the HRV and respiratory flow signals, and their coherence. A relevant issue in clinical practice is to distinguish between patients who successfully pass the spontaneous trial but later need reintubation and those who pass the test and can maintain spontaneous breathing.

Two classification methods have been studied, applying PCA to the spectral signal (S-Class), and the sequential floating feature selection to parameters extracted from these spectral signals (F-Class).

The S-Class method was found to be well suited to classifying the different groups of patients undergoing weaning trials, particularly using MSC signal when comparing successful patients group. As well as a good balance between sensitivity and specificity, it provided a high accuracy of 92% classifying GSucc versus GFail patients (sensitivity of 92%, specificity of 94%), and 86% classifying GSucc versus GRein patients (sensitivity of 85%, specificity of 93%). In contrast, when GFail and GRein groups were compared, the best parameters were obtained with the respiratory flow signal, with and accuracy of 79% (sensitivity of 82%, specificity of 73%), whereas with the MSC parameters the accuracy is 67%.

Applying the F-Class method the best classification was obtained comparing GSucc versus GFail groups with the respiratory flow signal, and an accuracy of 74%. This result is not so good as the obtained with S-Class method (an accuracy of 88%).

The results with HRV signal showing a similar behavior in all classifications, and with a low discriminatory capability. Their accuracies are between 33% and 46%, with a good balance among sensitivity and specificity.

Overall, the performance of the S-Class method is better than that of F-Class method. This could be related to the dimensionality reduction of the spectral signal, compared to the analysis of features, it being possible to capture more information through the spectral signals.

The most complex group to analyze is that of reintubated patients, because their initial response to a spontaneous breathing test is similar to that of successful patients, but within 48 h, the evolution of the respiratory pattern is closer to that of the patients who failed the test. In this study, among the 88 patients who had passed the spontaneous breathing trial, 15 (17%) patients needed reintubation within 48 h after extubation.

Considering the clinical point of view, the patients that had passed the spontaneous breathing trial but finally were classified into the reintubated group, has failed in the factors that predict the evidence-based extubation process. Consequently, our hypothesis is that the behavior of the spectral coherence is similar in reintubated and failed patients, while the respiratory pattern is different. The results suggest a difference behavior in the MSC signal between group of successful patients compared with failed and reintubated patients, respectively. However, this signal is similar when comparing failed group of patients with reintubated (accuracy of 67%), while the differences have been introduced by the respiratory flow signal (accuracy of 79%).

From a clinical perspective, recent studies show that the proportion of patients successfully weaned from mechanical ventilation at the first attempt is increasing over time, but that the reintubation rate has remained similar (Esteban et al 2013). The percentage of patients that need to be reintubated within 48 h (after successful weaning trials) ranges from 6 to 25% in different populations. These patients tend to require mechanical ventilation for significantly longer, and have longer ICU and hospital stays, and higher mortality (McConville and Kress 2012).

In our previous studies, we proposed respiratory pattern characterization based on statistical analysis of time series extracted from the respiratory and cardiac signals (Caminal et al 2010, Garde et al 2010, Garde et al 2013). In general, we obtained good results classifying GSucc and GFail patients, but the major challenge has continued to be classifying the reintubated patients.

Recently, authors as Huang et al (2014) and Seely et al (2014) are studying changes of heart and respiratory rate variability during the weaning process. They hypothesized that measuring their changes during this process may help clinicians to predict weaning results. The second paper propose to investigate the added value of these measurements for the prediction of extubation outcomes, both individually and in combination, using a predictive modeling.

Despite the improvements in the weaning process of patients undergoing mechanical ventilation, it is difficult to define the risk for extubation failure. With our results, the highlight conclusion is introducing new indices that help to discriminate especially between successful and reintubated patients.

An advantage of the method proposed in this paper, applying PCA analysis to the spectral signal, is the enhanced classification of different groups of patients, especially those requiring reintubation. Besides, the spectral coherence between HRV and the respiratory flow signal, in patients on weaning trial process, can contribute to the extubation decision. Nevertheless, these results should be evaluated on a more number of patients. Additional features and clinical information about the patients should be considered before weaning trials to increase the power of discrimination among these groups.

This work was supported in part by the Spanish Ministry of Economy and Competitiveness under grants TEC2010-21703-C03-01 and TEC2010-20886-C02-01. The authors would like to thank the staff of the Intensive Care Departments at Santa Creu i Sant Pau Hospital, Barcelona, and at Getafe University Hospital, Getafe, Spain, for their collaboration with the database acquisition.