BY-NC-ND 4.0 license Open Access Published by De Gruyter September 30, 2016

Resistance in a non-linear autoregressive model of pulmonary mechanics

Ruby Langdon, Paul D. Docherty, Bernhard Laufer and Knut Möller


Respiratory system modelling can enable patient-specific mechanical ventilator settings to be found, and can thus reduce the incidence of ventilator induced lung injury in the intensive care unit. The resistance of a simple first order model (FOM) of pulmonary mechanics was compared with a flow dependent term of a non-linear autoregressive (NARX) model. Model parameters were identified for consecutive non-overlapping windows of length 20 breaths. The analysis was performed over recruitment manoeuvres for 25 sedated mechanically ventilated patients. The NARX model term, b1, consistently decreased as positive end expiratory pressure (PEEP) increased, while the FOM resistance behaviour varied. Overall the NARX b1 behaviour is more in-line with expected trends in airway resistance as PEEP increases. This work has further verified the physiologically descriptive capability of the NARX model.

1 Introduction

Mechanical ventilation (MV) is essential for patients with acute respiratory distress syndrome (ARDS) in the intensive care unit (ICU). ARDS can be caused by a variety of injuries or illnesses, and it generally involves inflammation in the lungs, and increased pulmonary elastance [1]. Mechanical ventilation pushes air into the lungs and ensures adequate gas exchange is maintained. However, suboptimal ventilator settings can result in ventilator induced lung injury (VILI) [2]. Respiratory system modelling can enable patient-specific pulmonary mechanics to be captured, and aids clinicians in setting optimal ventilator settings for each patient, reducing the incidence of VILI [3].

There are a wide range of physiologically and clinically relevant models. However, simple models have been limited in their abilities to describe all relevant behaviour [4], [5], and more complex models have had issues with non-identifiability or practical usability [6], [7]. In particular, many models use a single term to represent the effects of airway resistance, which may be an unreasonable assumption [8].

A nonlinear autoregressive (NARX) model of respiratory mechanics has been proposed by Langdon et al. [9]. The NARX model has successfully captured airway pressure waveforms in ARDS patients during recruitment manoeuvres. The model uses basis functions to capture elastance across pressure, and multiple flow-dependent, time-dependent coefficients to capture resistance, passive lung relaxation, and viscoelastic effects. In particular, the NARX model has been more successful than the first order model (FOM) at capturing the expiratory relaxation, and both the pressure relaxation and oscillation during an end-inspiratory pause.

In this paper, we aim to evaluate the changes in resistance captured by the NARX model in response to PEEP changes. When PEEP is increased, resistance is expected to decrease as the higher pressure causes recruitment of alveoli, and a widening of the bronchial path vessels [10]. The ability of the NARX model to capture the reduction in resistance at higher PEEP levels will confirm its ability to describe expected behaviour and will validate this element of the modelling approach.

2 Material and methods

2.1 Data

The data was obtained from a study conducted in the ICUs of eight German hospitals between 2000 and 2002. Full details of this study are available in Stahl et al. [11]. Measurements were taken from 28 patients who suffered from ARDS. A wide range of patient conditions were captured, as the cause of ARDS varied, the length of ventilation was between two and 19 days, and the patient age ranged from 17 to 77 years. The protocol was approved by the local ethics committee of each participating institution, and informed consent was obtained from the patient, or his or her legally authorised representative.

The volume controlled method of ventilation was used, with a tidal volume of 8 ± 2 ml/kg body weight. An end-inspiratory pause of ≥ 0.2 s was used. Sedatives were titrated to achieve a Ramsay sedation score of 4–5, and neuromuscular blocking drugs were administered as needed. Therefore patient breathing efforts did not occur.

Airway pressure and flow were measured at the airway opening, using a pneumotachometer and piezoresistive transducer. The volume was calculated from continuous integration of the flow, with adjustment for volume drift. The sampling rate was 125 Hz. The sampling rate was later reduced to 62.5 Hz for analysis, in order to reduce noise in the data and improve the speed of simulations.

During the study, patients underwent recruitment manoeuvres (RMs). We have selected a data set from each patient that spans 8–10 min. During this period, patients were ventilated at zero PEEP (ZEEP) for approximately 5 min, before PEEP was increased in steps of 2 cmH2O after every 5–10 breathing cycles. This continued until a peak inspiratory pressure of approximately 50 cmH2O was reached.

The analysis in this paper was performed on 25 patients. One patient of the 28 was excluded because a recruitment manoeuvre was not performed. Two patients were excluded due to unusually high intrinsix PEEP, causing nonlinear behaviour that the NARX model was unable to capture well [9].

2.2 Respiratory models

The first order model is the foundation of the NARX model:


where: Paw is the measured airway pressure (cmH2O), R is the airway resistance (cmH2Os/l), V˙ is the airway flow rate (l/s), E is the pulmonary elastance (cmH2O/l), V is the inspired volume (l), and P0 is the end-expiratory pressure (cmH2O).

The NARX model builds upon the structure of the FOM:


where: M is the number of b-spline basis-functions to be used, i is the index of a particular basis function of degree d, ai is the coefficient for a given basis function (cmH2O/l), and i,d(Paw(t)) is the basis function value for a given pressure measurement. The sum of the basis functions multiplied by their ai coefficients defines a pressure dependent elastance.

There are j = 1 … L ⋅ bj coefficients (cmH2Os/l) that capture the pressure responses that occur due to flow and changes in flow. The subscript −j in the second term refers to the previous time sample. Thus, each Paw(t) is calculated from information from the previous L data points.

The b1 coefficient corresponds most closely to the resistance term R in the FOM. The coefficients b2bL capture other effects related to the pressure response to flow. These may include the effects of inertance, relaxation of the lung during expiration, and viscoelastic effects occurring during the end-inspiratory pause. The b allow the NARX model to capture both the end-inspiratory pause and the expiration curve much more accurately than the FOM.

2.3 Analysis

The NARX model was identified on a moving window of 20 breaths that shifted across the recruitment manoeuvre of each data set (Figure 1). Coefficients ai and bj were identified independently for each of these windows, allowing the trend of the resistance term b1 to be observed over time as PEEP increased.

Figure 1 Three pressure windows of length 20 breaths.

Figure 1

Three pressure windows of length 20 breaths.

Previous work determining the optimal parameters of the NARX model for this cohort found that d = 1, L = 350, and M = 5 provided a good fit to the data, when identified on the entire recruitment manoeuvre with a pressure range of 0–50 cmH2O [9]. When identified over a smaller pressure range, a smaller value of M is appropriate to provide a robust result. Therefore, M = 3 was used in this analysis, as the NARX model was identified on breaths at one to three adjacent PEEP levels only.

The NARX model used in this analysis has 353 parameters. The number of data points per breath is in the range 190–380 for different patients, due to varied breathing rates of patients, and a constant sampling frequency of 62.5 Hz. Therefore the window length must be at least two breaths long to avoid non-identifiability of the NARX model. A window length of 20 breaths was chosen as it enabled robust and stable b1 values that exhibited the underlying behaviour of the patients.

To observe how b1 varied with PEEP, boxplots for the 25 data sets were created. The average PEEP for each window was calculated, and the boxplots were plotted using data from windows where the average PEEP was 0.1–5, 5.1–10, 10.1–15, and 15.1–20 cmH2O.

The same analysis was performed using the R value of the FOM. All analysis was undertaken on an i7 quad core PC with 16GB RAM using MATLAB 2014a 64 bit functions and the statistical toolbox (Mathworks, Natick, MA).

3 Results

Figure 2 shows how the NARX b1 term and the FOM R term changed over time, for three data sets. Separate linear trend lines have been plotted over the ZEEP portion and over the recruitment manoeuvre portion. These three cases show the range of behaviours observed in this cohort. In the recruitment manoeuvre, the NARX b1 tended to decrease. In contrast, the FOM resistance exhibited three different types of behaviour. The FOM resistance increased, decreased, or had a quadratic shape as PEEP increased.

Figure 2 The NARX b1 (left) and the FOM R (right) for three patients. T = time at the start of the window.

Figure 2

The NARX b1 (left) and the FOM R (right) for three patients. T = time at the start of the window.

Table 1 shows the mean gradient of the linear trend line applied to the ZEEP and RM portions for b1 and R. Within the 90% confidence interval, the RM gradient of b1 was negative, as b1 consistently decreased with PEEP. In fact there was only one patient for which the RM b1 gradient was positive. In contrast, the varied behaviour of R as PEEP increased meant that the 90% confidence interval for the R gradient in the RM surrounded zero. The t-test revealed a significant difference between the mean gradient of b1 and R during the recruitment manoeuvres. There was no significant difference in the behaviour of R and b1 at ZEEP.

Table 1

The mean gradient of the NARX b1 and FOM R terms for the ZEEP and RM portions.

Mean ZEEP gradientMean RM gradient
and 90% CIand 90% CI
NARX b1−0.011 [−0.021, −0.0012]−0.031 [−0.041, −0.021]
FOM R−0.0012 [−0.0027, 0.0003]−0.0035 [−0.011, 0.0036]

The boxplots for the b1 and R terms (Figure 3) shows that the median of b1 decreased as PEEP increased. The median of R also decreased with PEEP, but with a much shallower slope. The boxplot at PEEP = 10.1–15 cmH2O contains 23 data sets out of the total of 25, because the highest PEEP for two patients was less than 10 cmH2O. Similarly, the boxplot at PEEP = 15.1–20.1 cm H2O contains data from 17 patients only.

Figure 3 Boxlots for the NARX b1 and FOM R parameters, normalised to the average value at ZEEP. The box limits are the 25th and 75th percentiles, and the whiskers show the range limited to data points that are within 1.5 IQR.

Figure 3

Boxlots for the NARX b1 and FOM R parameters, normalised to the average value at ZEEP. The box limits are the 25th and 75th percentiles, and the whiskers show the range limited to data points that are within 1.5 IQR.

4 Discussion

The resistance to flow is primarily described by the NARX model via variance in the b1 parameter. Figure 2 and Figure 3 show that the resistance captured by the NARX b1 term consistently reduced as PEEP increased, over 25 recruitment manoeuvres. This behaviour is concomitant with expected behaviour, as high pressures cause widening of airway passages, thus reducing resistance.

In contrast, the behaviour of the FOM R term was inconsistent. In response to the RM, the R term either decreased, increased, or had a quadratic shape (Figure 2). While Figure 3 shows that the general trend of R as a function of PEEP was downwards, the variance exhibited by the FOM was significantly larger than that of the NARX b1 value.

The average trend in both b1 and R during ZEEP was a very shallow negative slope. The 90% confidence interval indicated no consistent trend in R at ZEEP, and the t-test showed no significant difference in the gradient of R and b1 at ZEEP (Table 1). This behaviour was expected, as patient’s resistance at ZEEP would be likely to be roughly constant, unless patient condition changes.

The R term of the FOM had inconsistent behaviour across patients. Many patients experienced an increase in this modelled resistance as PEEP increased, from either the beginning of the RM or from midway through the RM. It is suspected that those patients that had an apparent increase in modelled resistance at higher pressures were actually exhibiting non-linear elastance behaviour that was being compensated for by increased resistance values in the FOM. In contrast, the non-linear elastance behaviour was captured by the more complex NARX model, due to the multiple elastance parameters and the other flow-dependent terms. Thus, the b1 term was robust to this non-linearity.

This work has further validated the descriptive capability of the NARX model for capturing changes in airway resistance over PEEP steps. Hence the NARX may offer a unique methodology for observing changes in ARDS patients, most notably, those suffering from chronic obstructive pulmonary disease (COPD). Future work will assess the physiological relevance and numerical robustness of the b2 – bL coefficients that are necessary to capture the viscoelastic effects observed in the end expiratory pause, and expiratory phase [9].

Author’s Statement

Research funding: The author state no funding involved. Conflict of interest: Authors state no conflict of interest. Ethical approval: The protocol was approved by the local ethics committee of each participating institution, and informed consent was obtained from the patient, or his or her legally authorised representative.


[1] Slutsky AS, Ranieri VM. Mechanical ventilation: lessons from the ARDSNet trial. Respir Res. 2000;1:73–7.Search in Google Scholar

[2] Dreyfuss D, Saumon G. Ventilator-induced lung injury. Lessons from experimental studies. Am J Respir Crit Care Med. 1998;157:294–323.Search in Google Scholar

[3] Chiew YS, Chase JG, Shaw G, Sundaresan A, Desaive T. Model-based PEEP optimisation in mechanical ventilation. BioMed Eng Online.2011;10:111.Search in Google Scholar

[4] Mount LE. The ventilation flow-resistance and compliance of rat lungs. J Physiol. 1955;127:157–67.Search in Google Scholar

[5] van Drunen E, Chiew YS, Chase J, Shaw G, Lambermont B, Janssen N, et al. Expiratory model-based method to monitor ARDS disease state. BioMed Eng Online. 2013;12:57.Search in Google Scholar

[6] Sundaresan A, Yuta T, Hann CE, Chase JG, Shaw GM. A minimal model of lung mechanics and model-based markers for optimizing ventilator treatment in ARDS patients. Comput Methods Programs Biomed. 2009;95:166–80.Search in Google Scholar

[7] Schranz C, Docherty PD, Chiew YS, Möller K, Chase JG. Identifiability analysis of a pressure-depending alveolar recruitment model, in IFAC Proceedings Volumes; 2012. p. 137–42.Search in Google Scholar

[8] Mols G, Kessler V, Benzing A, Lichtwarck-Aschoff M, Geiger K, Guttmann J. Is pulmonary resistance constant, within the range of tidal volume ventilation, in patients with ARDS? Bri J Anaesth. 2001;86:176–82.Search in Google Scholar

[9] Langdon R, Docherty PD, Chiew Y-S, Möller K, Chase JG. Use of basis functions within a non-linear autoregressive model of pulmonary mechanics. Biomed Sig Process Contr. 2016;27:44–50.Search in Google Scholar

[10] Damanhuri NS, Docherty PD, Chiew YS, van Drunen EJ, Desaive T, Chase JG. A patient-specific airway branching model for mechanically ventilated patients. Comput Math Methods Med. 2014.Search in Google Scholar

[11] Stahl CA, Moller K, Schumann S, Kuhlen R, Sydow M, Putensen C, et al. Dynamic versus static respiratory mechanics in acute lung injury and acute respiratory distress syndrome. Crit Care Med. 2006;34:2090–8.Search in Google Scholar

Published Online: 2016-9-30
Published in Print: 2016-9-1

©2016 Ruby Langdon et al., licensee De Gruyter.

This work is licensed under the Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 License.