key: cord-1013727-px6wtop3 authors: Minucci, Sarah; Heise, Rebecca L.; Valentine, Michael S.; Kamga Gninzeko, Franck J.; Reynolds, Angela M. title: Mathematical modeling of ventilator-induced lung inflammation date: 2020-11-17 journal: bioRxiv DOI: 10.1101/2020.06.03.132258 sha: c42a3b90d40601c9bea1db6a909e5679533c9a07 doc_id: 1013727 cord_uid: px6wtop3 Respiratory infections, such as the novel coronavirus (SARS-COV-2) and other lung injuries, damage the pulmonary epithelium. In the most severe cases this leads to acute respiratory distress syndrome (ARDS). Due to respiratory failure associated with ARDS, the clinical intervention is the use of mechanical ventilation. Despite the benefits of mechanical ventilators, prolonged or misuse of these ventilators may lead to ventilation-associated/ventilation-induced lung injury (VILI). Damage caused to epithelial cells within the alveoli can lead to various types of complications and increased mortality rates. A key component of the immune response is recruitment of macrophages, immune cells that differentiate into phenotypes with unique pro- and/or anti-inflammatory roles based on the surrounding environment. An imbalance in pro- and anti-inflammatory responses can have deleterious effects on the individual’s health. To gain a greater understanding of the mechanisms of the immune response to VILI and post-ventilation outcomes, we develop a mathematical model of interactions between the immune system and site of damage while accounting for macrophage polarization. Through Latin hypercube sampling we generate a virtual cohort of patients with biologically feasible dynamics. We use a variety of methods to analyze the results, including a random forest decision tree algorithm and parameter sensitivity with eFAST. Analysis shows that parameters and properties of transients related to epithelial repair and M1 activation and de-activation best predicted outcome. Using this new information, we hypothesize inter-ventions and use these treatment strategies to modulate damage in select virtual cases. 1 Inflammation occurs in the lungs when an immune response is initiated 2 to eliminate an insult. Types of insults include inhaled pathogens, such 3 pneumonia, tuberculosis, SARS-COV-2, or other harmful particles. In the 4 most severe cases this leads to acute respiratory distress syndrome (ARDS). Due to respiratory failure associated with ARDS, the clinical intervention 6 is the use of mechanical ventilation. When individuals have a severe form 7 of COVID-19, the disease caused by SARS-COV-2, the disease can lead to 8 respiratory failure and death of the patients. In a recent study, two-thirds of 9 patients admitted for COVID-19 required mechanical ventilation [1] . 10 Despite the benefits of mechanical ventilators, prolonged or misuse of 11 these ventilators may lead to ventilation-induced lung injury (VILI). In this 12 work we will focus on the tissue damage associated with mechanical venti-13 lation and resulting immune cell recruitment. The damage caused to alve-14 olar sacs (clusters of alveolar cells) during mechanical ventilation can lead 15 to volutrauma (extreme stress/strain), barotrauma (air leaks), atelectrauma 16 (repeated opening and closing of alveoli), and biotrauma (general severe in-17 flammatory response). If the trauma increases, it can lead to multi-system 18 organ failure [2, 3] . 19 It has also been shown that the inflammatory response of the elderly is 20 altered in the lungs and other areas [4, 5] . As compared to younger indi- Mortality or discharge to extended-care facilities increased for each decade of 26 age greater than 65 years [7] . Additionally, the case fatality rate for COVID-27 19 patients over 70 years old and over 80 years old was around 50.8% and 28 14.8% of the total number of deaths, respectively [8] . This is in agreement 29 with other studies reporting higher rates of severe outcomes in patients with 30 COVID-19 aged 65 or older [9] . The change in the inflammatory response 31 with patient age combined with the increased need for ventilation and in-32 creased mortality rate among the elderly stresses the need to investigate the 33 influence of aging in VILI. The framework we have built here addresses VILI 34 with various parameters and initial conditions that can be narrowed in future 35 studies with data from different age groups and/or insults to explore dynam-36 ics and driving factors in various diseases related to age and/or outcome. 37 We used mathematical modeling to investigate the role of the pulmonary 38 immune response and treatments in ventilator-induced damage. We adapted The primary focus of this model is to examine the effects of damage on 154 the alveolar epithelium, in particular alveolar type II cells, since they are 155 responsible for restoration of the epithelium. In this section we begin with a 156 simple model, concentrating on the novel aspect of incorporating epithelial 157 cells and relative damage due to inflammation. We then add variables to 158 more accurately model the dynamics within this system. 159 We begin with a small three-dimensional system of differential equations, 160 shown in Eqs (1)- (3) , where E h is the proportion of the local space filled by 161 healthy cells, E d is the proportion of the local space filled by damaged cells, 162 and E e represents dead cells or empty "space" that can be replaced/filled 163 with healthy cells. Each term represents a biological event explained by 164 the brackets above the term. This first model includes only the baseline 165 abilities of epithelial cells to proliferate and repair themselves in the pres-166 ence of sustained damage. We do not explicitly model proliferating and 167 non-proliferating cells; the parameter p is modulated to reflect the general 168 mechanism by which neighboring epithelial cells renew surrounding "space" (4)-(5). In the presence of sustained stretch (s > 0), the E d nullcline switches from a vertical line to a line with slope (r + b)/s. The second equilibrium point changes from (0, 1) to Figure 1 : The epithelial subsystem generates a transcritical bifurcation for the parameter p. Bifurcation diagram for the proliferation parameter p for the epithelial system with stretch and no immune response. Other parameters are set to r = 2.6, s = 0.22, and b = 0.74. The unstable equilibrium below p < p * = 0.497 is not included in the figure, since it is not biologically relevant. in Eq (7) represents this loss of damaged cells. The stability analysis is similar to that from the model without the immune response, with additional parameters m, n that can shift steepness of the nullcline or the speed at which the system approaches or diverges from an equilibrium. The parameter p once again plays an important role in the stability of the two critical points, (0, 0) and There is a transcritical bifurcation when the value of p is varied; given interactions and events involved in VILI which we will explore in the next 254 section. 255 Figure 3 : Schematic describes interactions between immune system components. Immune system components shown in this schematic are macrophages, neutrophils, various pro-and anti-inflammatory mediators, and epithelial cells. Green boxes represent various types of neutrophils, colored circles represent naive, M0, M1, and M2 macrophages. White boxes represent healthy, damaged, and dead epithelial cells/empty space (E h , E d , and E e , respectively). Other colored boxes represent various types of mediators that are produced by epithelial and/or immune cells and signal to immune cells. Unactivated immune cells become activated by various mediators (p b , a b in the blood and p, a locally) and perform either pro-inflammatory or anti-inflammatory roles which are meant to remove debris (E e ) and promote repair of damaged epithelial cells. Dynamics between cells and mediators in the blood (not shown) are similar to the detailed dynamics shown for local inflammation. By adding variables to the two-dimensional system proposed above, we 257 developed a system of coupled ordinary differential equations to model the in- this is not always the case. One way to include aspects of the spatial hetero-268 geneity without explicitly modeling space is to use a compartmental model. Each compartment represents a well mixed environment and, when biologi- Proliferation of healthy cells, upregulated by PIM term s p , and we also model natural decay of these mediators. Anti-inflammatory mediators, such as the anti-inflammatory signaling these AIMs to inhibit this process is controlled by a b∞ . We also account for 456 intrinsic decay of neutrophils in the last term of Eqs (21) through (24). Neutrophils go through a multi-step process of rolling along and subse- Table 2 for ranges used for each parameter. Using LHS with these ranges we generated 100,000 parameter sets. Future Many of these sets had initial conditions associated with a severely in-548 flamed lung without ventilation, which did not seem biologically realistic. To correct for this we eliminated sets based on their initial condition for 550 E e (empty/dead cells). We performed all of the analysis below with three 551 different thresholds to see whether the exclusion of these parameter sets af-552 fected the results. In this paper we focus on the 23,086 parameters sets that 553 had E e (0) < 50% and show a summary of all results for E e (0) < 25% and 554 E e (0) < 75% in the supplementary materials. We did not find any major 555 differences when varying this inclusion threshold. Simulations were separated into three categories of disease progression: can be used. We included supplementary predictors calculated from the 628 transients, described in Table 3 . Adding these predictors allowed for the Our aim is to understand how recruitment of the immune response and We generated 100,000 parameter sets using LHS with parameter ranges 668 given in Table 2 Table 4 shows a summary of 735 the results from the various methods used to examine predictors' significance 736 in determining model output. Column 1 of Table 4 shows the predictors in 737 which all three groups were different from one another, as determined by the Random Forest eFAST (Ordered) (Not ordered) (Ordered output) 0.5h 2h 6h Table 4 : Summary of three different methods used to determine the most influential predictors, including parameters and other factors. Columns 1 & 2 show results for all 23,086 parameter sets. Column 1: significance testing results for predictors in which all three outcome groups are statistically different (p-value < 0.01). For ease of comparison between columns, the predictor is listed next to its counterpart in the ordered random forest list, if listed in that column. Column 2: average importance values determined by random forest decision trees. The top ten are ordered from highest to lowest importance. Table 4 ). The consistency of 772 the importance of these parameters and predictors using different methods 773 supports the idea that they play a significant role in the sensitivity of model Starting with a parameter set that gives rise to an E h transient that starts healthy and ends in a moderate inflammation state, we applied various treatment strategies by changing three key parameters, b r (rate at which healthy epithelial cells self-repair), k mne (rate of collateral damage to epithelial cells by macrophages and neutrophils), and x mne (Hill-type constant which regulates the effectiveness of macrophages and neutrophils in damaging epithelial cells We intervened in a case that starts healthy and ends in moderate inflam-794 mation. Note in Fig 11, we used Latin hypercube sampling to find biologically meaningful parameter 883 sets, producing a total of 23,086 acceptable parameter sets. This "virtual 884 cohort" produces a variety of dynamics that can be generated by the model. 885 We classified parameter sets into categories of healthy, moderate inflamma- transients that are linked to outcome and to determine candidate treatments. 890 We utilized several methods to determine the most important parameters 891 for model output, particularly epithelial health. Using eFAST, a sensitivity 892 analysis method for non-linear, non-monotonic ODEs, we found parameters 893 that, when fluctuated, caused a statistically significant difference in out-894 put than that generated by a dummy parameter. We then compared these 895 results with more non-conventional and less computationally intensive meth-896 ods. The random forest decision tree algorithm generated values denoting 897 the importance of parameters and other predictors on epithelial health and 898 is particularly useful for large data sets, such as the parameter sets in our 899 virtual cohort. Additionally, significance testing determined statistically sig-900 nificant differences in parameters grouped by outcome. 901 We were able to not only include parameter values in this analysis but 902 also other predictors later found to be important, including the M1 peak ratio Analysis showed that properties and parameters related to epithelial re-911 pair and M1 activation and de-activation were especially predictive of out-912 come. We used b r , k mne , x mne , and k en to simulate treatments for a parameter 913 set in the virtual cohort that started healthy and ended in a moderate inflam-914 mation disease progression. We found that modulating b r is effective in most 915 cases, and the other four can be helpful in some. The chosen case responded 916 differently to treatments and these were paired with varied M1 activation 917 dynamics, indicating that macrophage activation is tied to epithelial health 918 in VILI. Our approach of developing a virtual cohort and selecting important pa-920 rameters is a first step in identifying the driving mechanisms behind VILI 921 and how they contribute to outcomes. However, experimental data will be 922 necessary to better understand the immune response to VILI and identify bio- collected, which can be explored in future work. Another area of further study is determining why some virtual cases can 928 recover with a short intervention time while others need indefinite treatment. 929 We hypothesize that this has to do with patient-specific initial conditions and Parameters and other predictors that show a statistically significant difference (p-value<0.01) between all three disease progression groups, using Kruskal-Wallis and Wilcoxon tests. Criteria: Ee(0)<75% Ee(0)<50% Ee(0)<25% Covid-19: most patients require mechani-950 cal ventilation in first 24 hours of critical care Cytokines and biotrauma in ventilator-induced lung injury: a critical 956 review of the literature Ventilator-Induced Lung In-959 jury Inflammation, chronic ob-963 structive pulmonary disease and aging, Current Opinion in Pulmonary 964 Diabetes, Hyperglycemia, 968 and Inflammation in Older Individuals: The Health, Aging and 969 Body Composition study Characterization of Lung 974 Inflammation and its Impact on Macrophage Function in Ag-975 ing Age, duration of mechanical ventilation, and outcomes of 980 patients who are critically ill Characteristics of and Important Lessons 982 From the Coronavirus Disease 2019 (COVID-19) Outbreak in China: 983 Disease Control and Prevention Short-term 988 outcomes in individuals aged 75 or older with severe coronavirus 989 disease (COVID-19): First observations from an Infectious Dis-990 eases Unit in Southern Italy Neutrophil Recruitment and Func-1028 tion in Health and Inflammation Contribution of Neutrophils to Haemophilus 1037 Influenzae Induces Neutrophil Necrosis Ageing and the Immune 1042 System: Focus on Macrophages Ad-1047 vanced Age Impairs Macrophage Polarization A Biomath-1053 ematical Model of Pneumococcal Lung Infection and Antibi-1054 otic Treatment in Mice A Mathematical Model of Intrahost Pneu-1059 mococcal Pneumonia Infection Dynamics in Murine Strains Mathematical Model 1064 of a Three-Stage Innate Immune Response to a Pneumococcal Lung 1065 Modeling the Immune Rheostat 1069 of Macrophages in the Lung in Response to Infection A Systems Perspective 1074 of Host-Pathogen Interactions: Predicting Disease Outcome in 1075 Tuberculosis Identi-1079 fying Control Mechanisms of Granuloma Formation Dur tuberculosis Infection Using an Agent-Based Model Within-Host Influenza Dynamics: A Small-Scale Math-1086 ematical Modeling Approach Boolean Mod-1090 eling of Cellular and Molecular Pathways Involved in Influenza Infection A Dynamical Model 1095 of Human Immune Response to Influenza A Virus Infec-1096 tion An Agent-Based 1101 Model of Inflammation and Fibrosis Following Particulate Exposure 1102 in the Lung Kimp-1106 ton The Role of Inflammation Resolution 1108 Speed in Airway Smooth Muscle Mass Accumulation in Asthma: 1109 Insight from a Theoretical Model Multiscale CT-Based Computational Modeling of Alveolar 1114 Gas Exchange during Artificial Lung Ventilation, Cluster (Biot) and Pe-1115 riodic A Computa-1119 tional Model of Unresolved Allergic Inflammation in Chronic 1120 Asthma Strain-induced 1125 inflammation in pulmonary alveolar tissue due to mechanical ven-1126 tilation Time dependence of recruit-1130 ment and derecruitment in the lung: a theoretical model Predicting ventilator-induced lung injury using a lung injury cost 1137 function The Pressure-Volume Curve Is Greatly 1142 Modified by Recruitment, American Journal of Respira-1143 tory and Critical Care Medicine Pida-1147 parti, Quantification of Age-Related Lung Tissue Me-1148 chanics under Mechanical Ventilation Determinants and Limits 1153 of Pressure-Preset Ventilation: a Mathematical Model of Pressure 1154 Modeling the dy-1157 namics of recruitment and derecruitment in mice with acute lung 1158 injury Anal-1163 ysis for Stress Environment in the Alveolar Sac Model A Mathemat-1168 ical Model of Pulmonary Gas Exchange Under Inflammatory 1169 Cytokine Re-1173 sponse Is Determined by Duration of Receptor and Sig-1174 nal Activation Mathematical Modeling of Pro-and Anti-Inflammatory Sig-1180 naling in Macrophages, Processes Epithelial Repair Mechanisms 1183 in the Lung Lung Epithelial Wound 1188 Healing in Health and Disease Acute Lung Injury How Macrophages Orchestrate Resolution of Inflammation 1193 and Tissue Repair Simulating, Analyzing, and Animating Dy-1197 namical Systems, Software, Environments and Tools, Soci-1198 ety for Industrial and Applied Mathematics Alternative Activation of Macrophages Key Mechanisms Governing Resolution of Lung Inflamma-1209 tion Neutrophils: Cinderella of Innate Immune 1213 Pulmonary Macrophage Subpopulations in the In-1218 duction and Resolution of Acute Lung Injury, American Jour-1219 nal of Respiratory Cell and Molecular Biology Neutrophil Kinetics in Health and 1224 Disease Anti-inflammatory cytokines Hub-1230 mayr, Stretch Induces Cytokine Release by Alveolar Epithe-1231 lial Cells in Vitro M2 Macrophages Induced by Prostaglandin E2 and IL-6 from Cervical Carcinoma Are Switched to Activated M1 Macrophages 1239 by CD4+ Th1 Cells Phagocyte partnership during the on-1243 set and resolution of inflammation Comparison of Three 1247 Methods for Selecting Values of Input Variables in the Analysis of 1248 Output from a Computer Code A methodology 1251 for performing global uncertainty and sensitivity analysis in systems 1252 biology Matlab Functions for PRCC and eFAST Effects 1258 of Age on the Synergistic Interactions between Lipopolysaccha-1259 ride and Mechanical Ventilation in Mice Mechanical ventilation using non-injurious 1266 ventilation settings causes lung injury in the absence of pre-existing 1267 lung injury in healthy mice Sensitivity Analysis 1270 in Practice: A Guide to Assessing Scientific Models An alternative way to com pute Fourier amplitude sensitivity test (FAST), Computa-1274 tional Statistics & Data Analysis A Quantitative 1278 Model-Independent Method for Global Sensitivity Analysis 1279 of Model Output Study of the sensitivity of coupled reac-1284 tion systems to uncertainties in rate coefficients. I Theory Study of the sensitivity of coupled 1289 reaction systems to uncertainties in rate coefficients. II Applica-1290 tions An Evaluation with the Fourier Ampli-1294 tude Sensitivity Test (FAST) of Which Land-Surface Parameters Are 1295 of Greatest Importance in Atmospheric Modeling Classification and Regression by randomForest The Corsini En-1302 cyclopedia of Psychology Controlling the false discovery rate: a 1306 practical and powerful approach to multiple testing Total number of sets that reached steady-state Severe inflammation ES: 2659 (208) 2314 (208) Numbers in parentheses are the number of sets that leave the state and enter the state at the end of the simulation for initial condition (IC) and ending state (ES) [10] M. Torres, J. Wang, P. J. Yannie, S. Ghosh, R. A. Segal, A. M. 994 Reynolds, Identifying important parameters in the inflammatory pro-