Impact Factor 3.367 | CiteScore 4.3
More on impact ›

Original Research ARTICLE

Front. Physiol., 30 October 2020 |

Modeling Lung Derecruitment in VILI Due to Fluid-Occlusion: The Role of Emergent Behavior

  • 1Department of Medicine, Vermont Lung Center, Larner College of Medicine, The University of Vermont, Burlington, VT, United States
  • 2Department of Telecommunications and Control Engineering, University of São Paulo, São Paulo, Brazil
  • 3Department of Bioengineering, College of Engineering, Design & Computing, University of Colorado Denver, Aurora, CO, United States
  • 4Department of Biomedical Engineering, Boston University, Boston, MA, United States

Ventilator-induced lung injury (VILI) is driven by the processes of volutrauma and atelectrauma, which can act synergistically to compromise the blood-gas barrier. We have postulated that this synergy arises through a rich-get-richer mechanism whereby atelectrauma causes holes to form in the blood-gas barrier while concomitant volutrauma causes susceptible holes to progressively enlarge as VILI worsens. We previously developed an analytical model based on this idea that accurately predicts the progressive increases in lung elastance seen immediately following a recruitment maneuver as VILI progresses over the course of hours. In the present study we extend this model to account for the rate of change of elastance, due to closure of lung units, in the minutes following a recruitment maneuver. We found that the distribution of unit closing velocities throughout the lung can be described by a power law with an exponent of −2 that matches previously published power laws associated with the dynamics of lung recruitment. Our model thus reveals lung collapse as an example of emergent complex behavior and links the dynamics of altered function in the injured lung to structural damage in a way that explains the mechanisms of injury progression arising from the ongoing stresses and strains applied by mechanical ventilation.


Mechanical ventilation plays a major role in critical care, providing life support to patients in respiratory failure. However, because of the stresses and strains it imposes on the tissues of the lung, mechanical ventilation can also cause ventilator-induced lung injury (VILI) (Dos Santos and Slutsky, 2000; Gattinoni et al., 2003; Slutsky and Ranieri, 2013). The two principle injury mechanisms of VILI are over-distension of the tissues that gives rise to volutrauma, and cyclic recruitment and derecruitment of lung units that gives rise to atelectrauma (Parker et al., 1993; Dreyfuss and Saumon, 1998). Volutrauma and atelectrauma are both involved in compromising the blood-bas barrier of the lung, which allows plasma-derived fluid and proteins to leak into the airspaces. This disrupts surfactant function, which then decreases lung compliance and thus increases the tissue stresses associated with mechanical ventilation, making the leak worse. At the same time, the lung becomes progressively more derecruited, which places the remaining open lung regions at greater risk of over-distension (Gattinoni and Pesenti, 2006; Gattinoni et al., 2016). Patients with pre-existing lung injury, particularly acute respiratory system distress syndrome (ARDS), are highly susceptible to VILI (Slutsky and Ranieri, 2000; Dreyfuss and Hubmayr, 2016). Various ventilation strategies have been proposed to manage ARDS and prevent VILI, but once underway it is difficult to break out of the vicious cycle of progressing injury (Slutsky and Ranieri, 2013).

We have shown in a number of studies in mice (Seah et al., 2011; Smith et al., 2013, 2017) that the nature and degree of lung injury is reflected in what we term the derecruitability of the lung. Derecruitability is measured by applying a recruitment maneuver (deep inflation) to maximally open closed airspaces followed immediately by a period (typically 3 min) of mechanical ventilation at a designated tidal volume (Vt) and level of positive end-expiratory pressure (PEEP). Lung stiffness (H) rises progressively during the ventilation period to a degree that is modest (perhaps 20%) in a normal lung, but is dramatically enhanced in the injured lung to an extent that reflects the degree of injury (Smith et al., 2013). Some of the post-recruitment increase in H in the normal lung likely reflects changes in surface tension arising from the dynamics of surfactant at the air-liquid interface and viscoelastic adaptation of tissue stress, but the highly exaggerated increases in H seen in the injured lung are largely the result of progressive derecruitment of airspaces (Massa et al., 2008; Smith et al., 2017; Hamlington et al., 2018b).

We denote as Drate the mean rate of increase in H over the 3 min following recruitment, while H1 is the first measurement of H made immediately after recruitment (Smith and Bates, 2013; Smith et al., 2013). Both Drate and H1 increase with the severity of lung injury, but they do not mirror each other exactly. Drate quantifies derecruitment dynamics that manifest over a timescale of minutes, probably reflecting instabilities in the fluid layer lining the small airways and alveoli (Smith et al., 2013). H1, on the other hand, quantifies either much more rapid derecruitment phenomena or an inability to recruit at all, more likely explicable on the basis of rapid alveolar flooding by accumulated edema fluid and/or unstable collapse of lung units caused by high surface tension in the air-liquid interface (Smith et al., 2013). Drate and H1 together thus comprise a pair of sensitive biomarkers of lung injury (Smith et al., 2013).

We have previously shown in mice that, as VILI develops, the time-course of H1 can be accurately accounted for in terms of the accumulated leak of plasma-derived material through holes in the blood-gas barrier (Mori et al., 2018). These holes appear to be generated by both atelectrauma and volutrauma acting synergistically in a rich-get-richer fashion (Hamlington et al., 2018a). The rate of accumulation of edema in the airspaces, and thus the displacement of air that results, is a function of the number and sizes of these holes, both of which increase with time during injurious mechanical ventilation. In the present study, we extend this analysis to Drate in order to establish a consistent theory of derecruitability of the injured lung that links both H1 and Drate to the biophysical processes taking place during the development of VILI.

Materials and Methods

Model Development

In the following model development we take a similar analytical approach to our previous study of the dependence of H1 on degree of lung injury, this time focusing on the complementary parameter Drate. To account for the progressive increase in Drate that accompanies the development of VILI, it is necessary to model events taking place over two different scales of time. The longer time scale, on the order to minutes to hours, concerns the progression of VILI itself. The shorter time scale, on the order of seconds to minutes, involves the dynamic recruitment and derecruitment of alveoli and small airways that are reflected in our periodic assessments of VILI. In order to avoid ambiguity between these two time scales in the following mathematical development, events taking place over the shorter time scale (seconds to a few minutes) will be referenced to time denoted by τ while events taking place over the longer time scale (many minutes to hours) will be referenced to time denoted by t.

We first consider the short time-scale assessment of VILI, which is based on the computational model of recruitment/derecruitment dynamics proposed by Bates and Irvin (2002) consisting of a parallel distribution of lung units. Each unit has its own airway to which is applied the same pressure. The units open and close according to their randomly assigned critical opening and closing airway pressures (Po and Pc, respectively). These pressures are drawn from probability distribution fPo(Po) and fPc(Pc), respectively. However, a unit does not necessarily open or close as soon as it traverses its particular value of either Po or Pc. Rather, a latency is built into the process by associating with each unit a virtual trajectory variable x that can assume any value between 0 (corresponding to the unit being definitively closed) and 1 (definitively open). An open unit remains open until x=0, while a closed unit remains closed until x=1. Movement of x along the virtual trajectory occurs to the right (i.e., toward increasing values of x) when P > Po, and to the left when P < Pc, with the rate of change x being proportional to the amount by which P either exceeds Po or is less than Pc. The constants of proportionality that determine the velocity of x rightward and leftward are s0 and sc, respectively. That is,

d x d τ _ = ( P - P o ) s o ; for P P o
= ( P - P c ) s c ; for P P c
= 0 ; for P c < P < P o (1)

where 0≤x <  1. The velocity constants s0 and sc for each unit are drawn from probability distributions fso(s0) and fsc(sc), respectively.

Following the original model (Bates and Irvin, 2002), we assume that Po and Pc are distributed according to Gaussians with variances σ0 and σc, and means μ0 and μc, respectively. In accord with the findings of Massa et al. (2008) in acid-injured mice, we set the two variances to be equal (i.e., σ0 = σc = σ) and the two means separated by δP such that the opening pressures are greater than the closing pressures (i.e., μ0 = μc + δP)

f P c ( P c ) = 1 σ 2 π e - 1 2 ( P c - μ c σ ) 2 (2)


f P o ( P o ) = 1 σ 2 π e - 1 2 ( P o - μ c - δ P σ ) 2 (3)

We assume that σ and δP do not change with increasing injury, so the progression of VILI is reflected entirely in the way that μc increases with time.

The velocity distributions fso(s0) and fsc(sc) were originally specified on purely empirical grounds as conforming to hyperbolic functions (Bates and Irvin, 2002), motivated by the ubiquitous appearance of such functions in naturally occurring complex systems (West, 1990; Robertson et al., 2000; West et al., 2000). Solving the model analytically using such functions is not possible, however, because hyperbolas are not finite-integrable over s ∈ (0,∞), but finite integrals can be obtained by having the exponent of the power law be different from unity and by specifying a small positive lower limit, sm, for integration. Accordingly, we set

f s c ( s c ) = α - 1 s m ( s c s m ) - α (4)

where 0 < sm < sc < ∞ and α > 1.

When modeling how the fraction of closed lung, Fclosed, changes during a derecruitability test, we avoid having to consider fso(s0) specifically by assuming, as did Albert et al. (2009), that it is only derecruitment that takes place to any significant extent during the 3 min following a recruitment maneuver. That is, we assume that the fully recruited lung derecruits progressively during these 3 min as if exposed to an effective pressure, Pderecruit, that is fixed throughout the breath and that reflects the mean airway pressure applied by the mechanical ventilator. Pderecruit is determined by the applied level of PEEP plus some positive ΔP that depends on the ventilation regimen that is applied on top of the PEEP. As the lung becomes progressively more injured, and thus stiffer, mean airway pressure increases. This pressure increase is offset to some extent, however, by the compliance of the gas in the piston and connecting tubing of the flexiVent ventilator (Scireq, Montreal, Canada) that was used to collect the experimental data. While the effect of this compliance on calculation of H was compensated for digitally (Schuessler and Bates, 1995), it still absorbed an increasing fraction of the pressure increases as the lung stiffened (Smith et al., 2013). Mean airway pressure also increases somewhat over the course of a single derecruitability test. However, in order to obtain an analytical solution to our model we assumed that Pderecruit remains unchanged over each of the 3 min derecruitment maneuvers. This is clearly not precisely the case, but we take it to be a reasonable first-order approximation since tidal volume is quite low during the tests.

For a unit to have the possibility of closing during a derecruitability test it must therefore satisfy the condition Pc > Pderecruit. Whether or not a particular unit actually closes during the test, however, depends on whether it moves along its virtual trajectory rapidly enough for x to reach a value of 0 before the 3 min test is complete and that there is no significant retrograde motion of x during the inspiratory portions of the test. If the duration of the test were to continue indefinitely, the asymptotic value to which Fclosed would tend is that which satisfies the above condition, which is

F closed ( ) = P derecruit [ 1 σ 2 π e - 1 2 ( P c - μ c σ ) 2 ] d P c (5)

This situation is never actually achieved, of course, because the derecruitability test is terminated after 3 min.

We assume that the recruitment maneuver establishes an initial condition for which all units open and are positioned at x=1 along their respective virtual trajectories. As derecruitment proceeds during the test, the various units in the model close sequentially. The time, τ, at which a given unit closes is determined by how quickly it traverses the length of its virtual trajectory from x=1 to x=0, which is given by

τ = 1 s c | P derecruit - P c | (6)

which means that the units that have closed by time τ are those that have closing pressures satisfying 1τ(Pc-P)<sc. Assuming, as did Albert et al. (2009), that fPc(Pc) and fsc(sc) are statistically independent, Fclosed(P,τ) is given by

F closed ( P , τ ) = P [ ( P c - P ) τ ] - 1 1 σ 2 π e - 1 2 ( P c - μ c σ ) 2 α - 1 s m ( s c s m ) - α d P c d s c (7)

A detailed step-by-step solution of Eq.7 is presented in the Supplementary Material.

Solving Eq. 7 we find that α = 2 and Eq. 4 of thus becomes

f s c ( s c ) = 1 s m ( s c s m ) - 2 (8)

Massa et al. (2008) found σ = 3 in mice with acute lung injury caused by hydrochloric acid instillation. Smith et al. (2013) found μc(t) = kt, where k is a constant, in injuriously ventilated mice. Since PEEP is zero and H is measured using small-amplitude perturbations in lung volume, we assume that P=0 to a first approximation (Mori et al., 2018). Drate can thus be written as:

D rate ( t ) = H 1 ( t ) s m 2 π [ π 2 3 k t ( e r f ( k t 3 2 ) + 1 ) + 9 e - 1 2 ( k t 3 ) 2 ] (9)

From our previous proposed model (Mori et al., 2018) we can relate Eq. 9 to H1 when t is sufficiently large:

lim t D rate ( t ) = 3 k s m 1 a ( H 1 ( t ) - H 0 2 ) (10)

Experimental Data

We tested the predictions of the above model against a set of data collected in mice and published previously (Smith et al., 2013). A detailed description of the experimental method used, as well as animals ethics approval, are given in this previous publication. Briefly, healthy 8- to 10-week-old female BALB/c mice (18.3–24.4 g) were anesthetized with 90 mg/kg intraperitoneal sodium pentobarbital and then connected to a computer-controlled mechanical ventilator (flexiVent, Scireq, Montreal, Canada). Following paralysis (0.5 ml/kg IP injection of pancuronium bromide) animals received a 21 min protocol starting with 16.5 min of injurious ventilation delivered at zero PEEP with a very large Vt (1.0, 1.1, 1.2, and 1.3 ml in 4 different groups of animals) followed by a derecruitability test (4.5 min protocol consisting of a deep inflation followed by ventilation at zero PEEP with Vt = 0.25 ml during which elastance was measured every 20 s). This 21 min protocol was repeated continuously for 4 h or until the animal died. The way that the lungs derecruit over time during the derecruitment test is reflected in the way that H increases with time (Allen and Bates, 2004). Finally, the resulting estimations of Drate from the derecruitability tests were fitted to Eq. 15, and the 95% joint confidence regions for sm and k were determined as described by Mori et al. (2018).

It should be pointed out that we also measured parameters related to the resistance of both the conducting airways and the respiratory system tissues in these mice. We focus here exclusively on the elastance parameter H in our model because it scales inversely with the fraction of lung that remains open, provided intrinsic tissue stiffness does not change, thus serving as a convenient biomarker of derecruitment. Lung tissue resistance scales similarly to elastance (Smith et al., 2020) and therefore adds no new information, so we did not include it in our model. Airway resistance tends to correlate less well with lung derecruitment because the proximal branches of the airway tree, which contribute much to airway resistance, mostly remain patent even when significant portions of the lung periphery derecruit (Smith et al., 2020).

Model Fitting

The model described by Eq. 9 was fit to the elastance data versus time collected within each derecruitment maneuver with MATLAB (The Mathworks, Natick, MA, United States) using a non-linear least squares method (trust-region algorithm) for each animal individually in all 4 groups, resulting in a pair of values of sm and k for each animal. The maximum number of evaluations and iterations were set to 600 and 400, respectively, and the termination tolerance value of the cost function (mean squared residual) was 10–6. Finally, the 95% joint confidence regions for sm and k were determined as described by Mori et al. (2018).


Figure 1 shows the fits of Eq. 9 to experimental data from four groups of over-ventilated mice (Vt = 1.0, 1.1, 1.2, and 1.3 ml, with zero PEEP), together with the experimental data. The model fits pass within one standard error either side of the data points in most cases, and are very close to these ranges in the other cases. The R2 values for all fits were greater than 0.85. This demonstrates that the model accurately accounts for the way that injury, as reflected by H, accelerates over time and how this acceleration increases with increasing Vt.


Figure 1. Drate (average ± SE) versus time as VILI developed in the 4 groups of mice (symbols) together with the model fits to Eq. 9 (lines). Black: V_t = 1.0 ml (RMSE = 0.35); Red: V_t = 1.1 ml (RMSE = 0.64); Green: V_t = 1.2 ml (RMSE = 0.36); Blue: V_t = 1.3 ml (RMSE = 0.65). The data were collected in a previous study (Smith et al., 2013) (RMSE – root mean squared residual.).

Figure 2 shows the two model parameters in Eq. 9, sm and k, versus Vt. There was no dependence of sm on Vt by ANOVA (p = 0.11, Figure 2A), whereas k was significantly greater for the lungs exposed to the two higher Vt compared to the two lower Vt (p < 0.001, Figure 2B). Normality of the distributions of all parameter groups was supported by the Shapiro-Wilk test (p > 0.05). The Bonferroni post hoc test was performed to assess inter-groups differences for k. All statistical tests were performed in MATLAB 2019b (The Mathworks, Natick, MA, United States).


Figure 2. Mean ± 2 SEM for the two key model parameters s_m (A) and k (B) and for the 4 groups, 1.0 ml (n = 6), 1.1 ml (n = 6), 1.2 ml (n = 5) and 1.3 ml (n = 4). Statistically significant differences between groups are indicated by: *P < 0.05), **P < 0.01.

The joint 95% confidence intervals show clean separation between the parameter pairs for the four groups of mice (Figure 3).


Figure 3. Joint paired 95% confidence regions for the model parameters. The optimum parameter values (those that minimize the root mean squared difference between measured and modeled Drate) are located at the centers of the ellipses. Black: V_t = 1.0 ml; Red: V_t = 1.1 ml; Green: V_t = 1.2 ml; Blue: V_t = 1.3 ml.

Figure 4 shows plots of Drate(t) versus H1 for all animals studied. Also shown are the model predictions of Drate(t) from Eq. 9 versus H1 calculated according to Mori et al. (2018). Lastly, Figure 5 shows, in a Bland-Altman plot, the difference between modeled and measured Drate values as a function of H1.


Figure 4. Rate of change of lung stiffness (Drate) versus initial measurement of post-recruitment stiffness (H_1) measured in the 4 groups of mice (symbols) together with the model fits (lines). Black: V_t = 1.0 ml; Red: V_t = 1.1 ml; Green: V_t = 1.2 ml; Blue: V_t = 1.3 ml.


Figure 5. Difference between modeled and measured Drate as a function of H_1 (open dots). Mean is represented by the full line whereas 95% interval (± 1.96 SD) is represented by dashed lines.


We have developed an analytical model of the time-course of derecruitment in the injured lung that mimics observations of how lung derecruitability evolves during the development of VILI. This analytical model is based on our previous numerical approach of ascribing recruitment and derecruitment dynamics to the behavior of a virtual trajectory associated with various lung units (Bates and Irvin, 2002; Massa et al., 2008). This model uses a similar analytic approach to our recent study that examined how H measured immediately after a recruitment maneuver (H1) increases as a function of worsening lung injury (Mori et al., 2018). In this sense, the present study complements our previous model analysis since together they provide a complete analytic description of how a deep inflation of the injured lung recruits closed lung units and then how these units derecruit again over time as governed by the values of the three fixed parameters α, sm, and k, where α and k depend on Vt. The model of the present study accurately predicts the manner in which the temporal evolution of Drate varies with both time and injurious Vt (Figure 1), dependencies that are embodied in a single parameter, k, thought to largely reflect how surface tension increases with time (Figure 2). Surface tension elevation is a consequence of surfactant dysfunction, itself caused either by direct damage to the surfactant or by the presence of material that has leaked from the vasculature into the airspaces of the lung. H can also increase as a result of accumulated fluid that floods alveoli. This model thus links developments in physical injury of the parenchyma tissue to the dynamics of lung function derangement via the physiologic mechanism of airspace derecruitment.

Analytical models generally do not have the same flexibility as numerical models to describe the complex details of experimental data. On the other hand, when the model equations are based on putative underlying mechanisms, as in the present case, analytical models may offer insight of a general nature that more empirical nature of numerical models cannot. In this regard, the model in the present study supports the suitability of a power-law function as describing the velocity distribution along virtual trajectories (Eq. 4). The parameter sm sets the minimum closing velocity, which corresponds to those units that close so slowly they remain open throughout the entire derecruitability test. It is possible that these units would eventually close if the derecruitability tested lasted long enough, but they make no contribution to the data when the tests last only 3 min. Interestingly, the exponent of the power law was found to be −2 purely on the basis of the empirically observed linear increase in lung stiffness following a recruitment maneuver at zero PEEP. In other words, this exponent value comes from the qualitative behavior of the data rather than the fit of a particular model. This is significant because power laws have been invoked in a previous model of lung recruitment during lung inflation to account for avalanches of sequential opening events along the airway tree that are triggered when critical opening pressures are surpassed (Suki et al., 1994, 2000). Moreover, both inter-crackle intervals and terminal air space sizes that are linked by avalanches throughout the airway tree are characterized by power-law distributions having an exponent of −2 (Suki et al., 2000; Alencar et al., 2001). The similarity between the probability distribution for sc and these other avalanche-based mechanisms in the lung might not be a coincidence, but rather might follow in some way from the fractal branching structure of the airway tree. Thus, although our model consists of a parallel collection of respiratory units that cannot open in the kinds of cascades observed in the airway tree (Suki, 2002), the fractal structure of the lung is embedded in our model through the power-law form of the distribution of sc (Eq. 4) and thus conforms to previous work by Suki et al. (1994, 2000), Alencar et al. (2001), and Suki (2002).

Having the geometry of the airway tree be implicit in our model in the functional form for sc could also explain the findings of Massa et al. (2008) that fso(s0) and fsc(sc) are unchanged following lung injury. That is, although lung injury would be expected to damage the alveolar parenchyma and the airway epithelium, there is no reason to suppose it would alter the branching structure of the airway tree, which would leave the functional form for sc preserved. On the other hand, we do expect that injury would increase the mean values of the distributions for Pc and Po (Eqs 2 and 3, respectively) on the basis of increases in surface tension at the air-liquid interface within the lungs as a result of blood-gas barrier perforation and the subsequent airspace accumulation of fluid and protein (Massa et al., 2008; Smith et al., 2013). This conveniently allows us to account for the progression of lung injury largely in terms of the evolution of the value of a single parameter, k in Eq. 9 (Figures 2, 3).

We previously reported an apparent master relationship between and Drate and H1 (Smith et al., 2013). The model developed in the present study, however, suggests this is not quite the case because rate of change of Drate with injury (Eq. 10) is predicted to depend on k, which depends on Vt (Figure 2B), even though the initial value of Drate (Eq. 10 of the Supplementary Material) does not (Figure 2A). On the other hand, the dependence of dDrate/dH0 on Vt is rather subtle (Figure 4), and is easily masked by noise in the data, so the notion of a single master relationship may still be a useful approximation. Also, both H1 (Mori et al., 2018) and Drate (Figure 4) exhibit exponential growth over the time-scale of our experiment (up to 4 h), but they are unlikely to continue to increase indefinitely beyond this point because this would correspond to derecruitment of the entire lung. Even if this was a possibility, however, forcing a sizeable Vt into an ever shrinking lung would eventually lead to pressures that would cause rupture of the tissues. In any case, it appears that neither eventually would ever be realized because the mice died when about 75% of the lung became derecruited, at which point H1 and Drate were still increasing exponentially.

Another feature of the relationship between Drate and H1 is the early transient seen in Figure 4, which shows that Drate to begin to increase early in the development of lung injury before H1 changes noticeably. We have previously postulated (Smith et al., 2013) that increases in Drate are due to increases in surface tension that give rise to derecruitment dynamics manifesting over a timescale of minutes, while increases in H1 reflect flooded alveoli that have become permanently derecruited. This being the case, the relationship between Drate and H1 in Figure 4 suggests that surfactant deactivation by plasma-derived fluid and proteins starts to occur as soon as this material arrives in the airspaces whereas a certain volume of fluid as to accumulate in the airspaces before it can cause noticeable alveolar flooding. It also suggests that Drate is more directly related to the “rich-get-richer” hypothesis of VILI that we have previously advanced (Hamlington et al., 2018a; Mori et al., 2018) and thus may be the more sensitive marker of the onset of lung injury. On the other hand, even though the model errors for Drate do not appear to depend on H1 (Figure 5), at large t, Drate, and H1 become linearly dependent on each other (Eq. 10) suggesting that both variables ultimately embody the same pathophysiological information, although the fits of the present study are not as good at this end of the data range compared to those of our previous study (Smith et al., 2013) which curve over and follow the data better. Regardless, both Drate and H1 are increased by fluid accumulation in the airspaces, so some kind of monotonic relationship between the two parameters is to be expected.

Our model has a number of limitations. For example, although we have characterized derecruitability of the lung in terms of only two parameters, H1 and Drate, the increase in H with time during a derecruitability test is only approximately linear in its early stages; eventually, H begins to plateau (Allen et al., 2002; Allen and Bates, 2004; Massa et al., 2008). Such a plateau is inevitable as only those lung units for which P < Pc will eventually close as τ→∞. Drate is thus an empirical characterization of only part of a more complex relationship. We also neglected intra-breath variations in the pressures applied to the lungs by mechanical ventilation, as well as the way they depend on the degree of lung injury. We assumed a constant pressure for derecruitment in order to be able to solve the model equations analytically (Eq. 7) and so gain the benefits of insight that a closed-form solution provides. In reality of course, these pressures are not constant because mean airway pressure rises for a given Vt as the lung becomes progressively more consolidated, so more quantitatively accurate predictions would require the numerical solution of a computational model in which airway pressures are allowed to vary in a realistic manner. This may partly explain why the ellipse areas of the joint parameter confidence regions in Figure 3 tend to be greater for larger Vt (i.e., greater VILI); we would generally expect a simple analytical model to fit better to data from normal or mildly injured lungs compared to those that are severely injured. We also waived details concerning alveolar micromechanics such as interactions between neighboring respiratory units, although we recognize the potentially important role of such interactions in lung injury (Knudsen and Ochs, 2018). Finally, we have assumed that the dynamics of derecruitability in the lung are accurate indicators of the degree of lung injury, but they are not direct measures of such. Nevertheless, numerous prior studies from our laboratory have established a tight link between derecruitability and direct measures of damage to the blood-gas barrier (Allen and Bates, 2004; Allen et al., 2005; Smith et al., 2017; Hamlington et al., 2018a,b).

In summary, we have developed an analytical version of our previously developed computational model of the dynamics of lung derecruitment following a recruitment maneuver. We found that the distribution of velocity constants defining rates of movement along the virtual trajectories associated with the pulmonary airways follows a power-law with the same exponent of −2 as found in other investigations into the dynamics of lung recruitment and derecruitment. We speculate that this may reflect the fractal-like structure of the airway tree throughout which recruitment and derecruitment may take place as cascades of open and closing events. Our study also suggests that the rate of change in lung stiffness represented by the parameter Drate, to the extent that this remains constant during a 3 min derecruitability test, is a particularly sensitive marker of the degree of lung injury.

Data Availability Statement

The datasets generated for this study are available on request to the corresponding author.

Ethics Statement

The animal study was reviewed and approved by Institutional Animal Care and Use Committee of the University of Vermont.

Author Contributions

VM developed the model analysis, coded the simulations, and drafted the manuscript. BJS, BS, and JB contributed to the development of the model concept, checked the mathematics, and edited and approved the manuscript. All authors contributed to the article and approved the submitted version.


This study was funded by NIH grants R01 HL-124052, R00 HL-128944, and R01 HL-142702 and Coordenação de Aperfeiçoamento de Pessoal de Nível Superior – Brazil (CAPES) – Finance Code 88881.135413/2016-01. These funding agencies contributed to the salary support of the authors of this study. R01 HL-142702 will provide the publication fees for this manuscript.

Conflict of Interest

BJS and JB filed the following patent that is potentially related to the subject matter of the manuscript: “Variable ventilation as a diagnostic tool for assessing lung mechanical function” JB and BJS. PCT Application WO2015127377 A1, Filed on February 23, 2014 (C538).

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at:


Albert, S. P., Dirocco, J., Allen, G. B., Bates, J. H., Lafollette, R., Kubiak, B. D., et al. (2009). The role of time and pressure on alveolar recruitment. J. Appl. Physiol. 106, 757–765. doi: 10.1152/japplphysiol.90735.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Alencar, A. M., Buldyrev, S. V., Majumdar, A., Stanley, H. E., and Suki, B. (2001). Avalanche dynamics of crackle sound in the lung. Phys. Rev. Lett. 87: 088101.

Google Scholar

Allen, G., and Bates, J. H. (2004). Dynamic mechanical consequences of deep inflation in mice depend on type and degree of lung injury. J. Appl. Physiol. 96, 293–300. doi: 10.1152/japplphysiol.00270.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

Allen, G., Lundblad, L. K., Parsons, P., and Bates, J. H. (2002). Transient mechanical benefits of a deep inflation in the injured mouse lung. J. Appl. Physiol. 93, 1709–1715. doi: 10.1152/japplphysiol.00473.2002

PubMed Abstract | CrossRef Full Text | Google Scholar

Allen, G. B., Pavone, L. A., Dirocco, J. D., Bates, J. H., and Nieman, G. F. (2005). Pulmonary impedance and alveolar instability during injurious ventilation in rats. J. Appl. Physiol. 99, 723–730. doi: 10.1152/japplphysiol.01339.2004

PubMed Abstract | CrossRef Full Text | Google Scholar

Bates, J. H., and Irvin, C. G. (2002). Time dependence of recruitment and derecruitment in the lung: a theoretical model. J. Appl. Physiol. 93, 705–713. doi: 10.1152/japplphysiol.01274.2001

PubMed Abstract | CrossRef Full Text | Google Scholar

Dos Santos, C., and Slutsky, A. (2000). Invited review: mechanisms of ventilator-induced lung injury: a perspective. J. Appl. Physiol. 89, 1645–1655. doi: 10.1152/jappl.2000.89.4.1645

PubMed Abstract | CrossRef Full Text | Google Scholar

Dreyfuss, D., and Hubmayr, R. (2016). What the concept of VILI has taught us about ARDS management. Intens. Care Med. 42, 811–813. doi: 10.1007/s00134-016-4287-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Dreyfuss, D., and Saumon, G. (1998). Ventilator-induced lung injury: lessons from experimental studies. Am. J. Respir. Crit. Care Med. 157, 294–323. doi: 10.1164/ajrccm.157.1.9604014

PubMed Abstract | CrossRef Full Text | Google Scholar

Gattinoni, L., Carlesso, E., Cadringher, P., Valenza, F., Vagginelli, F., and Chiumello, D. (2003). Physical and biological triggers of ventilator-induced lung injury and its prevention. Eur. Respir. J. 22, 15s–25s.

Google Scholar

Gattinoni, L., Marini, J. J., Pesenti, A., Quintel, M., Mancebo, J., and Brochard, L. (2016). The “baby lung” became an adult. Intens. Care Med. 42, 663–673. doi: 10.1007/s00134-015-4200-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Gattinoni, L., and Pesenti, A. (2006). “The concept of baby lung,” in Applied Physiology in Intensive Care Medicine, (Berlin: Springer), 303–311. doi: 10.1007/3-540-37363-2_44

CrossRef Full Text | Google Scholar

Hamlington, K. L., Bates, J. H., Roy, G. S., Julianelle, A. J., Charlebois, C., Suki, B., et al. (2018a). Alveolar leak develops by a rich-get-richer process in ventilator-induced lung injury. PLoS One 13:e0193934. doi: 10.1371/journal.pone.0193934

PubMed Abstract | CrossRef Full Text | Google Scholar

Hamlington, K. L., Smith, B. J., Dunn, C. M., Charlebois, C. M., Roy, G. S., and Bates, J. H. T. (2018b). Linking lung function to structural damage of alveolar epithelium in ventilator-induced lung injury. Respir. Physiol. Neurobiol. 255, 22–29. doi: 10.1016/j.resp.2018.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Knudsen, L., and Ochs, M. (2018). The micromechanics of lung alveoli: structure and function of surfactant and tissue components. J. Histochem. Cell Biol. 150, 661–676. doi: 10.1007/s00418-018-1747-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Massa, C. B., Allen, G. B., and Bates, J. H. (2008). Modeling the dynamics of recruitment and derecruitment in mice with acute lung injury. J. Appl. Physiol. 105, 1813–1821. doi: 10.1152/japplphysiol.90806.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

Mori, V., Smith, B. J., Suki, B., and Bates, J. H. (2018). Linking physiological biomarkers of ventilator-induced lung injury to a rich-get-richer mechanism of injury progression. Ann. Biomed. Eng. 47, 638–645. doi: 10.1007/s10439-018-02165-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Parker, J. C., Hernandez, L. A., and Peevy, K. J. (1993). Mechanisms of ventilator-induced lung injury. Crit. Care Med. 21, 131–143. doi: 10.3109/9781420019261-9

CrossRef Full Text | Google Scholar

Robertson, H. T., Altemeier, W. A., and Glenny, R. W. (2000). Physiological implications of the fractal distribution of ventilation and perfusion in the lung. Ann. Biomed. Eng. 28, 1028–1031. doi: 10.1114/1.1310215

CrossRef Full Text | Google Scholar

Schuessler, T. F., and Bates, J. H. (1995). A computer-controlled research ventilator for small animals: design and evaluation. IEEE Trans. Biomed. Eng. 42, 860–866. doi: 10.1109/10.412653

CrossRef Full Text | Google Scholar

Seah, A. S., Grant, K. A., Aliyeva, M., Allen, G. B., and Bates, J. H. (2011). Quantifying the roles of tidal volume and PEEP in the pathogenesis of ventilator-induced lung injury. Ann. Biomed. Eng. 39, 1505–1516. doi: 10.1007/s10439-010-0237-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Slutsky, A. S., and Ranieri, V. M. (2000). Mechanical ventilation: lessons from the ARDSNet trial. Respirat. Res. 1:2.

Google Scholar

Slutsky, A. S., and Ranieri, V. M. (2013). Ventilator-induced lung injury. New Engl. J. Med. 369, 2126–2136.

Google Scholar

Smith, B. J., Bartolak-Suki, E., Suki, B., Roy, G. S., Hamlington, K. L., Charlebois, C. M., et al. (2017). Linking ventilator injury-induced leak across the blood-gas barrier to derangements in murine lung function. Front. Physiol. 8:466. doi: 10.3389/fphys.2017.00466

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, B. J., and Bates, J. H. (2013). Assessing the progression of ventilator-induced lung injury in Mice. IEEE Trans. Biomed. Eng. 60, 3449–3457. doi: 10.1109/tbme.2013.2267151

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, B. J., Grant, K. A., and Bates, J. H. (2013). Linking the development of ventilator-induced injury to mechanical function in the lung. Ann. Biomed. Eng. 41, 527–536. doi: 10.1007/s10439-012-0693-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, B. J., Roy, G. S., Cleveland, A., Mattson, C., Okamura, K., Charlebois, C. M., et al. (2020). Three alveolar phenotypes govern lung function in murine ventilator-induced lung injury. Front. Physiol. 11:660. doi: 10.3389/fphys.2020.00660

CrossRef Full Text | Google Scholar

Suki, B. (2002). Fluctuations and power laws in pulmonary physiology. Am. J. Respirat. Crit. Care Med. 166, 133–137. doi: 10.1164/rccm.200202-152pp

PubMed Abstract | CrossRef Full Text | Google Scholar

Suki, B., Alencar, A. M., Tolnai, J., Asztalos, T., Peták, F., Sujeer, M. K., et al. (2000). Size distribution of recruited alveolar volumes in airway reopening. J. Appl. Physiol. 89, 2030–2040. doi: 10.1152/jappl.2000.89.5.2030

PubMed Abstract | CrossRef Full Text | Google Scholar

Suki, B., Barabási, A.-L., Hantos, Z., Peták, F., and Stanley, H. E. (1994). Avalanches and power-law behaviour in lung inflation. Nature 368:615. doi: 10.1038/368615a0

PubMed Abstract | CrossRef Full Text | Google Scholar

West, B. J. (1990). Physiology in fractal dimensions: error tolerance. Ann. Biomed. Eng. 18, 135–149. doi: 10.1007/bf02368426

PubMed Abstract | CrossRef Full Text | Google Scholar

West, G. B., Brown, J. H., and Enquist, B. J. (2000). “The origin of universal scaling laws in biology,” in Scaling in Biology, eds J. H. Brown and G. B. West (Oxford: Oxford University Press), 87–112.

Google Scholar

Keywords: acute lung injury, analytical model, alveolar flooding, surface tension, lung elastance

Citation: Mori V, Smith BJ, Suki B and Bates JHT (2020) Modeling Lung Derecruitment in VILI Due to Fluid-Occlusion: The Role of Emergent Behavior. Front. Physiol. 11:542744. doi: 10.3389/fphys.2020.542744

Received: 13 March 2020; Accepted: 08 October 2020;
Published: 30 October 2020.

Edited by:

Jeremy Andrew Simpson, University of Guelph, Canada

Reviewed by:

Robert Huhle, Dresden University of Technology, Germany
Martin Kneyber, University Medical Center Groningen, Netherlands

Copyright © 2020 Mori, Smith, Suki and Bates. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Vitor Mori,