I Introduction
The respiratory management of extremely preterm infants (gestational age 28 weeks) is specialized and challenging due to their immature lungs and inability to maintain functional breathing. Despite advances in neonatal care and noninvasive ventilation, these infants are a highrisk population that often requires invasive mechanical ventilation (IMV) during the beginning of life [1, 2].
Prolonged IMV is associated with several shortterm and longterm morbidities, such as pneumonia, lung injury, bronchopulmonary dysplasia, and neurodevelopmental problems [3, 4]. Therefore, clinicians emphasize the need for prompt extubation using ventilatory settings, blood gases, hemodynamic stability and clinical experience to determine readiness [4].
Unfortunately, extubation failure rates in this population are high and variable, ranging from 20% up to 70% depending on the failure criteria and the time frame used for definition [5]. Reintubation may be challenging and has been identified as an independent risk factor for increased morbidities and death [6]. The Automated Prediction of Extubation Readiness (APEX) research project aims to quantify cardiorespiratory behavior and evaluate its association with extubation outcomes, with the goal of assisting clinicians in the extubation decision.
As part of the APEX study, we have previously examined the ability of cardiorespiratory features [7] and clinical parameters [8] to predict extubation readiness. Cardiorespiratory behavior can also be characterized by segmenting the physiological time series data into patterns of respiratory behavior that are mutually exclusive. The resulting pattern sequence provides a highlevel description of respiratory behavior. This paper describes an approach to model this sequence using different types of Markov models. The objective was to identify distinguishing characteristics of the infants who succeeded extubation compared to those who failed.
Markov chain modeling provides a robust probabilistic mechanism for characterizing sequential data. It has been applied with great success in several domains including in natural language processing
[9] and in gene pathway analysis [10]. In this paper, we first examined the distribution of the respiratory patterns. We then modeled the patterntopattern transitions over time as discretetime Markov chains. Because usual discretetime Markov chains can be sensitive to the sampling rate of the signal, we also used semiMarkov models, which allow separate modeling of the time spent in each state, as well as the crosspattern transitions. Lastly, we fit multichain Markov models to examine the potential nonstationarity of the respiratory patterns over time.The rest of the paper is organized as follows: Section II describes the APEX study; Section III describes methodology and modeling details; Section IV reports the results of the experiments; and Section V discusses the conclusions of this work.
Ii APEX Study Design
Iia Infant Population
Eligible infants were recruited from Canada (Royal Victoria Hospital, Montreal Children’s Hospital, Jewish General Hospital, Quebec) and USA (Detroit Medical Center, MI, and Women and Infants Hospital of Rhode Island, RI) with Birth Weight (BW) 1250 g. Infants were recruited while undergoing IMV prior to their first extubation.
Infants were excluded if they had any major congenital anomalies such as heart disease, or were receiving any vasopressor or sedative drugs at the time of extubation. Ethics approval was obtained from the ethics board at each institution and data was collected after informed parental consent was obtained.
IiB Data Acquisition
Cardiorespiratory signals were acquired from infants in supine position in their incubators. Signals were collected for 60 min under Invasive Mechanical Ventilation (IMV), followed by a 5 min period of EndoTracheal Tube Continuous Positive Airway Pressure (ETTCPAP) immediately prior to extubation. Respiratory signals were measured using Respiratory Inductance Plethysmography (RIP) bands placed around the infant’s ribcage (RCG), at the level of the nipple line, and around the infant’s abdomen (ABD), 0.5 cm above the umbilicus (Viasys^{©}, Healthcare, USA).
Signals were sampled at 1000 Hz using the PowerLab 16/30 analogdigital data acquisition system with a 16bit analogtodigital resolution (ADInstruments, Bella Vista, Australia, 2009^{©}). Signals were antialias filtered at 500 Hz before acquisition. The signals were analyzed using MATLAB^{TM} (The MathWorks Inc.).
The outcome of the extubation was also recorded. A failed extubation outcome was defined as the need for reintubation within 7 days of extubation.
IiC Respiratory Patterns
RIP signals were analyzed using Automated Unsupervised Respiratory Event Analysis (AUREA), which extracts samplebysample metrics of respiratory power, synchrony between RCG and ABD, and movement artifact [11]
. AUREA uses kmeans clustering to assign each sample to one of the following respiratory patterns:

Pause (PAU): A cessation of breathing.

Synchronous Breathing (SYB): RCG and ABD are in phase.

Asynchronous Breathing (ASB): RCG and ABD are out of phase.

Movement Artifact (MVT): Associated with infant moving or nurse handling.

Unknown (UNK): Ambiguous patterns not belonging to any other pattern category.
AUREA provides repeatable results with no human intervention. An example of RIP signals and corresponding patterns assigned by AUREA to the different samples is shown in Fig. 1.
Iii Methodology
The objective of this work was to model the sequence of respiratory patterns for the 5 min period of ETTCPAP. Separate models were developed for infants who succeeded or failed extubation.
In the following subsections, we describe our initial analysis of the pattern distributions, which motivated the choice of Markov modeling. We then describe the three different types of Markov models we applied, based on the observed data characteristics.
Iiia Analysis of Pattern Distributions
We examined the distribution of pattern durations (for PAU, SYB, ASB, MVT, and UNK) to determine if there are marked differences over time during the ETTCPAP period. To do this, the ETTCPAP respiratory pattern signal is split into two equaltime segments for each patient. The length of each segment of a particular pattern was recorded. The total occurrences for each outcome group were used to plot a probability histogram. The result is a probability density histogram of the pattern durations within the first or second half of the ETTCPAP, for either failure or success patients.
IiiB Fitting a Discretetime Markov Chain
Problems which can be phrased as transitions among discrete states over time can be modeled using Markov chains [12]. A Markov chain model consists of a set of states (in our case, the 5 respiratory patterns); the distribution of initial states; and the transition matrix , which holds the probabilities of all possible statetostate transitions. At each times step the system is in some state, and transitions to a next states according to . In a homogeneous Markov chain, the probability of transitioning from state to state is independent of the time step , as well as of the history preceding state :
(1) 
This is known as the Markov property [12], which makes such models tractable by removing the need to compute full conditionals
. Empirically, the transition probabilities can be computed using maximum likelihood estimation:
(2) 
where is number of times the transition from state to state was observed in the data.
Note that in our data, the distribution of initial states is unreliable due to infant and device handling at the beginning of data collection episodes, so in fact we did not include it in the model.
IiiC SemiMarkov Model
SemiMarkov models are different from standard Markov chains in that selftransitions, i.e. transitions from a state to itself, are collapsed. Instead, each state has a dwell time distribution
, which models the duration spent in the state until a transition out of the state occurs, coupled with a probability distribution of transitioning to
other states. The latter results in a chain where every transition results in a state that differs from its predecessor.Using this framework is useful for two reasons. First, it increases the resolution of the offdiagonal elements of the transition matrix , especially when selftransitions are extremely likely (which is the case in our data, as will be shown in Sec. IV). Secondly, if the sampling rate of the data were to change, the simple Markov chain model would drastically change, whereas the semiMarkov model would not be significantly affected. For example, if the sampling rate were to double, assuming that this does not affect the state labeling, the selftransitions of a Markov chain that works at the sampling rate would roughly double; in contrast, the semiMarkov chain would still have the same dwell time distribution and same probability of transitioning from a state to its successors.
IiiD MultiChain SemiMarkov Model
Multichain Markov models can be used to examine nonstationarity over time in sequential data. Separate Markov chains are fit to different, nonoverlapping time segments of the data. The differences in the models can be used to identify and quantify nonstationarity. We fit two semiMarkov models to the first and second half of the 5minute respiratory sequence of patterns, partly as a soft start step, but also because modeling on a finer temporal scale (e.g. minutebyminute) would leave too few samples for robust modeling, especially given that certain patterns (states) have low overall frequency.
IiiE KullbackLeibler (KL) Divergence
The KullbackLeibler (KL) divergence is a measure of the deviation of one probability distribution from another. The KLdivergence of a distribution P from another Q is defined as:
(3) 
Th KL divergence is 0 if the two distributions are identical, and it goes to if does not have support for the whole domain of .
Note that the KLdivergence is nonsymmetric: , which is not desirable in our application. Hence, we apply symmetrized KLdivergence to compare distributions over pattern transitions, defined as:
(4)  
Iv Experiments and Results
At the time of this work, there were a total of 186 patients in our database: 136 extubation successes and 50 extubation failures.
Iva Analysis of Pattern Durations
Fig. 2
shows the mean fraction of ETTCPAP time spent in each of the five respiratory patterns for the success and failure groups, with error bars to display the standard deviation. Note that failure patients spent significantly longer time in the Pause pattern, and less time in Synchronous Breathing.
We also examined and fit parametric models to the histograms of durations of different respiratory patterns. Fig.
3shows the density histograms of the Pause pattern in the success and failure groups, fit separately from the two parts of the ETTCPAP. We also show an approximate fit using exponential distributions. The means of these distributions are different (2.2s in the success population, 3s and 2.7s in the failure population).
IvB Discretetime Markov chains
Discretetime Markov chains were fit to the sequence of respiratory patterns obtained during ETTCPAP for the success and failure populations. The transition matrices (not shown) had strong and equal diagonal elements (or selftransitions). In particular, the likelihood of going from every state to itself in both the success and failure population was the same: 0.99. As a result, the numerical resolution in the offdiagonal elements (crosspattern transitions) was very low (close to 0). Hence, the symmetric KLdivergence between the transition probabilities was close to 0 (0.001).
IvC Stationary SemiMarkov Model
We fit semiMarkov models to better visualize crossstate transitions. Tables I and II show the transition matrices for the success and failure populations, respectively. Given each start state (each row), the transition that is most likely is highlighted in bold font. It was observed that in both populations, the most probable transitions were the same for all except the Pause and Movement Artifact states. It should be noted, however, that even in the cases (Asynchrony, Synchrony and Unknown) where the the most probable transitions were same in both groups, the actual probabilities of these transitions differed. The symmetric KLdivergence of 0.27 also indicated a much larger difference than when using standard Markov models.
The dwell time distributions of the semiMarkov models were fit and are summarised in Table III. In each respiratory pattern, the distribution type which best fit the dwell times based on the Bayesian Information Criterion BIC [13] was found to be same for the success and failure population. However, the values of the parameters of these distributions differed.
IvD Experiment 3: Multichain SemiMarkov model
The 5 min ETTCPAP data was split into 2 halves for each patient. One semiMarkov chain was fit to the first 2.5 min and another to the second. The symmetric KL divergence between the models for the success and failure groups showed an increase from 0.30 in the first half to 0.37 in the second half.
PAU  ASB  MVT  SYB  UNK  

PAU  0  0.27  0.09  0.26  0.38 
ASB  0.10  0  0.16  0.29  0.45 
MVT  0.12  0.32  0  0.43  0.14 
SYB  0.06  0.25  0.15  0  0.54 
UNK  0.13  0.28  0.04  0.55  0 
PAU  ASB  MVT  SYB  UNK  

PAU  0  0.28  0.06  0.39  0.28 
ASB  0.12  0  0.21  0.28  0.40 
MVT  0.17  0.41  0  0.32  0.09 
SYB  0.14  0.21  0.14  0  0.52 
UNK  0.15  0.30  0.03  0.52  0 
Success  Failure  
Pause  Exponential  Exponential 
Asynchrony  GeneralizedExtremeValue  GeneralizedExtremeValue 
k=0.63, = 1.30, =1.85  k=0.65, = 1.36, =1.81  
Movement  GeneralizedPareto  GeneralizedPareto 
k=0.22, = 3.62  k=0.11, = 3.31  
Synchrony  InverseGaussian  InverseGaussian 
=8.61, = 3.61  =7.83, = 3.41  
Unknown  GeneralizedPareto  GeneralizedPareto 
k=0.07, = 2.07  k=0.10, = 2.05 
V Discussion
The experiments presented highlight interesting differences between extremely premature infants that succeed and those who fail extubation, which are apparent during ETTCPAP periods prior to extubation. Infants who go on to fail extubation spend a longer portion of the ETTCPAP period in Pause and less time in Synchronous Breathing than successfully extubated patients, as seen in Figs. 2 and 3.
The work also highlights the modeling power of Markov models for this type of data. While the standard Markov chains provided models with too many selfstate transitions, the use of semiMarkov models allowed modeling of crosspattern state transitions, showing that the transition probabilities between success and failure groups were indeed different. This was seen in Tables I and II. We also observed similarities between the populations: the most likely next respiratory state from the breathing (Synchronous and Asynchronous) and Unknown patterns were the same in both groups.
Modeling of the dwell time distributions also showed important similarities as well as differences between the two groups. Patients who fail extubation had longer occurrences of Pauses. The dwell time in both groups followed the same distribution type in each pattern, suggesting an underlying consistency in infant respiratory behavior that is unaltered by the extubation outcome.
The multichain semiMarkov model presents evidence of nonstationarity in the distribution of estimated transition probabilities over time: the transition probabilities of the two groups diverged from the first to the second half of ETTCPAP, as previously demonstrated in adults [14].
A Markov chain model encodes knowledge about state transitions and how frequently they occur. It makes the simplifying assumption that likelihoods of state transitions are independent of time, but as we have shown, nonstationarity in state transitions can be captured by fitting separate models to different time segments. The Markov models we developed based on patterns identified from RIP signals show potential for further use in probabilistic prediction of extubation readiness. We are currently carrying out work in this direction.
The semiMarkov models can also be used to simulate breathing behavior of a patient, by sampling respiratory states from the fitted transition model, and dwell times from the timing distribution associated with each state. This allows understanding possible trends of breathing state evolution for the two populations, which can also be used to provide a measure of variability in the trajectories.
The high probability of transitioning from SYB and ASB to UNK (Tables I, II) raises questions about the Unknown pattern. In particular, AUREA was validated in older infants who had surgery, during the postoperative period. In extremely preterm infants, the rib cage is highly compliant and as the diaphragm shortens, there is very little expansion of the rib cage with substantial increase in the motion of the abdominal wall [15]. Thus, we suspect that some of the UNK patterns found in the APEX study data may be related to predominant abdominal breathing with lowamplitude ribcage movements, which is characteristic of this population.
In conclusion, modeling results from this work give novel insights into the respiratory behavior of extremely preterm newborns in the period prior to extubation. The ultimate goal of this research is to build classifiers that can effectively determine whether or not an infant is ready for extubation. This work highlights several relevant methods that constitute a good basis for the future work of building these classifiers.
References
 [1] Michele C. Walsh, Brenda H. Morris, Lisa A. Wrage, Betty R. Vohr, W. Kenneth Poole, Jon E. Tyson, Linda L. Wright, Richard A. Ehrenkranz, Barbara J. Stoll, and Avroy A. Fanaroff. Extremely low birthweight neonates with protracted ventilation: Mortality and 18month neurodevelopmental outcomes. The Journal of Pediatrics, 146(6):798 – 804, 2005.
 [2] Martin Keszler and Guilherme Sant’Anna. Mechanical ventilation and bronchopulmonary dysplasia. Clinics in Perinatology, 42(4):781 – 796, 2015. Bronchopulmonary Dysplasia: An Update.
 [3] Rachel A. Joseph. Prolonged mechanical ventilation: Challenges to nurses and outcome in extremely preterm babies. Critical Care Nurse, 35(4):58–66, 2015.
 [4] H AlMandari, W Shalish, E Dempsey, M Keszler, P G Davis, and G Sant’Anna. International survey on periextubation practices in extremely preterm infants. Archives of Disease in Childhood  Fetal and Neonatal Edition, 2015.
 [5] Annie Giaccone, Erik Jensen, Peter Davis, and Barbara Schmidt. Definitions of extubation success in very premature infants: a systematic review. Archives of Disease in Childhood  Fetal and Neonatal Edition, 99(2):F124–F127, 2014.
 [6] Brett J. Manley, Lex W. Doyle, Louise S. Owen, and Peter G. Davis. Extubating extremely preterm infants: Predictors of success and outcomes following failure. The Journal of Pediatrics, 173:45 – 49, 2016.
 [7] D. Precup, C. A. RoblesRubio, K. A. Brown, L. Kanbar, J. Kaczmarek, S. Chawla, G. M. Sant’Anna, and R. E. Kearney. Prediction of extubation readiness in extreme preterm infants based on measures of cardiorespiratory variability. In 2012 Annual International Conference of the IEEE Engineering in Medicine and Biology Society, pages 5630–5633, Aug 2012.
 [8] P. Gourdeau, L. Kanbar, W. Shalish, G. Sant’Anna, R. Kearney, and D. Precup. Feature selection and oversampling in analysis of clinical data for extubation readiness in extreme preterm infants. In 2015 37th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), pages 4427–4430, Aug 2015.

[9]
Prerna Rai and Arvind Lal.
Google pagerank algorithm: Markov chain model and hidden markov model.
2016.  [10] D. Li and H. Q. Wang. A markov chain modelbased method for cancer classification. In 2012 8th International Conference on Natural Computation, pages 1064–1068, May 2012.
 [11] C. A. RoblesRubio, K. A. Brown, and R. E. Kearney. Automated unsupervised respiratory event analysis. In 2011 Annual International Conference of the IEEE Engineering in Medicine and Biology Society, pages 3201–3204, Aug 2011.

[12]
Andrei Andreevich Markov.
Extension of the law of large numbers to dependent quantities (in russian).
(2nd Ser.):135–156, 1906.  [13] Gideon Schwarz et al. Estimating the dimension of a model. The annals of statistics, 6(2):461–464, 1978.
 [14] M. Orini, B. F. Giraldo, R. Bailon, M. Vallverdu, L. Mainardi, S. Benito, I. Diaz, and P. Caminal. Timefrequency analysis of cardiac and respiratory parameters for the prediction of ventilator weaning. In 2008 30th Annual International Conference of the IEEE Engineering in Medicine and Biology Society, pages 2793–2796, Aug 2008.
 [15] A. Charles Bryan, Glenn Bowes, and John E. Maloney. Control of Breathing in the Fetus and the Newborn. John Wiley & Sons, Inc., 2011.
Comments
There are no comments yet.