category: literaturenote citekey: navasortegadatadrivenmodelingircu2025 title: Data-Driven Modeling of IRCU Patient Flow in the COVID-19 Pandemic authors: "Navas-Ortega, Ana Carmen; Sánchez-Martínez, José Antonio; García-Flores, Paula; Morales-García, Concepción; Fabregas, Rene" year: 2025 date: 2025-00-00 2025 doi: 10.1016/j.csbj.2025.10.017 url: "https://arxiv.org/abs/2506.09751v1" zotero_key: AI6QUPCU zotero_storage: W8TZV3HF collections: transpordi_planeerimine folder: 001_artiklid firstAuthor: "Navas-Ortega, Ana Carmen"
Ana Carmen Navas-Ortega‡1 , José Antonio Sánchez-Martínez‡1,2, Paula García-Flores1,2, Concepción Morales-García1,2, and Rene Fabregas3
1Department of Pneumology, University Hospital Virgen de Las Nieves, Granada, Spain. 2Biosanitary Research Institute of Granada-Ibs, Granada, Spain.
3Department of Applied Mathematics and Modeling Nature (MNat) Research Unit, Faculty of Sciences, University of Granada, Granada, Spain.
Corresponding Author: Concepción Morales-García, concepcion.morales.sspa@juntadeandalucia.es and Rene Fabregas, rfabregas@ugr.es.
Thursday 12th June, 2025
Background: Intermediate Respiratory Care Units (IRCUs) play a pivotal role between general wards and Intensive Care Units (ICUs), particularly during crises like the COVID-19 pandemic. Their effectiveness is influenced by structure, protocols, and clinical staffing expertise. This study aimed to evaluate the clinical outcomes and operational dynamics of a newly established IRCU featuring a specialized staffing model during the COVID-19 pandemic in Spain.
Methods: A prospective cohort study was conducted at the UHVN IRCU (Granada, Spain) between April and August 2021. Adult patients (n=249) admitted primarily for COVID-19 respiratory failure were included. Data on demographics, Non-Invasive Ventilation (NIV) use, length of stay (LOS), and outcomes (ICU transfer, exitus, recovery) were collected. Statistical analyses compared outcomes stratified by NIV status. A calibrated Ordinary Differential Equation (ODE) compartmental model and an empirical LOS-based convolution model were developed to simulate patient flow dynamics under various scenarios, including admission surges and varying care efficiency.
Results: The median age was 51 years; 31% (n=77) required NIV. Patients requiring NIV were significantly older (median 61 vs 42 years, p<0.001). Overall, 8% (n=20) required subsequent ICU transfer and 3% (n=7) experienced in-IRCU exitus. Notably, among the 172 patients managed without NIV, no ICU transfers or deaths occurred. Of the 77 high-risk patients requiring NIV, 68% recovered within the IRCU without ICU escalation. ODE modeling accurately reflected aggregate outcomes and demonstrated significant system strain during simulated admission surges, partially mitigated by enhanced recovery efficiency. LOS-based modeling yielded consistent peak occupancy estimates.
Conclusion: This IRCU, characterized by specialized clinical staffing including dedicated physicians and trained nurses, demonstrated effective management of severe COVID-19 respiratory failure. High recovery rates were achieved, particularly among NIV patients, potentially alleviating ICU pressure. Dynamic modeling confirmed the unit's vulnerability to admission surges but highlighted the positive impact of care efficiency. These findings highlight the positive outcomes observed in this appropriately structured and staffed IRCU, supporting the potential value of such units in pandemic response and respiratory care.
Keywords: Intermediate Respiratory Care Unit, COVID-19, Non-Invasive Ventilation, Respiratory Failure, Clinical Outcomes, Patient Flow Dynamics, Mathematical Modeling, Nursing/medical staff, Critical Care, Healthcare Operations.
‡ These authors contributed equally.
The effective management of severe respiratory failure represents a cornerstone of modern critical care medicine [1], a challenge dramatically underscored by the unprecedented pressures of the recent COVID-19 pandemic [2, 3, 4]. Within the complex architecture of tiered hospital care, Intermediate Respiratory Care Units (IRCUs) have emerged as a critical component [5, 6], strategically positioned to bridge the gap between the intensive care unit (ICU) and the general ward [7, 8, 9, 10, 11, 12, 13, 14]. Historically developed to provide focused monitoring and non-invasive respiratory support (NIV) [5, 6, 15, 16, 17], such as bilevel or continuous positive airway pressure (BiPAP/CPAP) and high-flow nasal oxygen (HFNO) [1, 18, 19], IRCUs proved indispensable during the pandemic surge [20, 21], functioning both "upstream" to prevent ICU admissions and "downstream" to facilitate timely ICU discharge [22, 23, 24, 25]. The potential benefits, including substantial reductions in the need for invasive mechanical ventilation [26, 27, 28] and associated healthcare costs [29, 17], are well-recognized [7, 30, 31].
Despite the accelerated proliferation of IRCUs and accumulating evidence of their utility [32, 14, 33], a critical knowledge gap persists regarding the specific determinants of their success [12, 34]. While structural factors and patient selection criteria are often discussed [8, 7, 35], the potentially significant impact of dedicated and specialized clinical staffing models, including physicians and nurses with advanced respiratory care competencies [36, 37, 38], remains inadequately quantified. Extensive evidence links general nurse staffing levels and education to improved patient outcomes [39, 40, 41, 42], and specialized physician involvement in critical care is known to be beneficial [43]. However, the specific impact of targeted multidisciplinary clinical expertise in the complex, high-stakes environment of an IRCU—managing sophisticated respiratory support [18, 44], identifying subtle signs of deterioration [45, 46], and enacting timely interventions through coordinated care [8, 47]—has not been systematically evaluated in terms of its contribution relative to other factors [14, 36, 48], representing a significant lacuna in our understanding of optimal intermediate care delivery.
Central to this investigation is the hypothesis that the operationalization of an IRCU, characterized by dedicated staffing including clinicians (both physicians and nurses) possessing specialized training in advanced respiratory support techniques [18, 38], is associated with superior clinical outcomes [26, 49] and optimized patient flow [50, 51, 52], even amidst the exacerbated demands of a pandemic [22, 23]. This study, therefore, embarks on a pioneering investigation to evaluate this premise. Our primary objective is to analyze the clinical effectiveness and operational dynamics of a newly instituted IRCU at the University Hospital Virgen de las Nieves (UHVN), Granada, during the initial five months of its operation amidst the COVID-19 pandemic (April-August 2021). This unit featured 24/7 respiratory physician coverage and nursing staff who received specific training (Fig. 6, 7). We aimed to quantify patient outcomes, including length of stay and ICU transfer rates [29, 24], and to elucidate the performance characteristics of this unit [11, 29], thereby providing crucial, data-driven insights into the outcomes observed in a setting implementing this model of specialized clinical staffing [36, 53, 54, 43].
The principal contribution of this work lies in its multifaceted approach, representing a significant methodological advance over prior descriptive accounts [15, 22, 55]. We present not only a detailed statistical analysis of patient characteristics and outcomes (Section 2.1), stratified by key risk factors such as the requirement for NIV, but also introduce a novel quantitative framework based on dynamic systems modeling [56, 57, 58]. This dual approach offers mechanistic insights into the temporal evolution of patient cohorts within the IRCU [59, 60], allowing for the simulation of system behavior under varying conditions [61, 62] and the identification of critical parameters influencing operational efficiency and patient fate [63, 64, 65]. By explicitly linking a specific care delivery model—characterized by dedicated staffing and enhanced clinical competency (exemplified by targeted nursing training)—to observed clinical results and modeled system dynamics, our study provides a systematic evaluation within this specific context, directly addressing the aforementioned knowledge gap and offering valuable evidence regarding the outcomes associated with implementing such staffing models in intermediate respiratory care [36, 66, 38].
Herein, we first detail the clinical characteristics and outcomes of the 249 patients managed within the UHVN IRCU during the study period, employing comprehensive statistical analysis to identify key determinants of patient trajectory (Section 2.1). Subsequently, we introduce and apply a compartmental Ordinary Differential Equation (ODE) model (Section 2.2) designed to capture the intricate dynamics of patient flow between non-NIV and NIV states and towards final outcomes (ICU transfer, exitus, recovery), exploring system behavior
through simulation under various operational scenarios. This is complemented by an alternative occupancy modeling approach grounded in the empirical Length-of-Stay distribution (Section 2.3). Finally, in the Discussion (Section 3) and Conclusion (Section 4), we elaborate on the implications of our findings for clinical practice, healthcare policy, and future research directions, contextualizing our results within the broader landscape of critical and intermediate care medicine.
A total of 249 patients admitted to the Intermediate Respiratory Care Unit (IRCU) between April and August 2021 were included in the analysis. The cohort exhibited a median age of 51 years (Interquartile Range [IQR]: 38–63 years), with ages spanning from 17 to 98 years (Figure 1A, left panel). Stratification by gender revealed a statistically significant difference in age distribution (Mann-Whitney U = 3405, p = 0.023), with female patients (n=98) presenting a higher median age (55 years, IQR: 40–67) compared to their male counterparts (n=151; median: 50 years, IQR: 37–60) (Fig. 1A, right panel).
Regarding clinical pathways within the entire cohort, 77 patients (31%) utilized Non-Invasive Ventilation (NIV), while 172 (69%) did not. Subsequent transfer to the Intensive Care Unit (ICU) was required for 8% (n=20) of patients, and exitus occurred in 3% (n=7) during the observed IRCU pathway (Fig. 1B, Total rows). While absolute counts stratified by gender in Fig. 1B suggest broadly similar outcome patterns, formal statistical tests confirmed no significant gender association for NIV utilization, ICU transfer, or exitus (see Fig. 2C and D; all p>0.05 [ns]).
The duration of IRCU stay showed considerable variability (Fig. 1C). The median length of stay (LoS) was 5 days (IQR: 1–12 days), with most stays concentrated below 10 days, although a subset required prolonged care exceeding 40 days. Monthly admission data fluctuated, peaking in June 2021 (n=75 admissions), potentially reflecting regional COVID-19 pandemic activity (Fig. 1C inset). Correlation analysis (Fig. 1D) identified expected significant positive associations (p<0.001 [\*\*\*]) between ICU transfer (T-ICU) and exitus, and between NIV utilization and both ICU transfer and exitus, suggesting higher severity among patients requiring NIV or ICU care. Gender was not significantly correlated with these variables (*p*>0.05).
Stratifying the cohort by NIV treatment status revealed critical clinical differences. Patients requiring NIV (n=77) were significantly older (median age 61, IQR: 56–67) compared to those managed without NIV (n=172; median age 42, IQR: 35–49) (Fig. 2A). This age difference strongly predicted NIV use; adjusted logistic regression showed a significant increase in the odds of receiving NIV per year of age (OR = 1.02, 95% CI: 1.01–1.04 [***], Fig. 2C), and crude risk analysis indicated that patients aged ≥70 years had double the risk compared to those <50 (RR = 2.0, 95% CI: 1.3–3.2 [***], Fig. 2D). As noted previously, gender was not significantly associated with NIV utilization (Fig. 2C and D; [ns]).
Final outcome distributions differed markedly between the two groups (Fig. 2B). The 172 patients managed without NIV experienced favorable outcomes, with no ICU transfers or exitus events recorded (100% recovery within these defined categories). In contrast, the 77 higher-risk patients requiring NIV had the following mutually exclusive outcomes: ICU transfer (n=18, 23%), exitus (without prior ICU transfer) (n=7, 9%), and recovery (n=52, 68%). The resulting high crude Risk Ratios comparing NIV users to non-users for ICU transfer (RR = 89.6, 95% CI: 5.5–1464, p<0.001 [***]) and exitus (RR = 31.4, 95% CI: 1.8–546, p<0.001 [***]) (Fig. 2D) must be interpreted in light of significant confounding by indication. Adjusted analyses confirmed that older age remained a significant risk factor for exitus specifically within the NIV group (OR = 1.07 per year, 95% CI: 1.01–1.13, p=0.013 [*], Fig. 2C), while gender showed no significant association with outcomes in this subgroup either (Fig. 2C and D; [ns]). Crucially, the successful recovery of 68% of NIV-treated patients without ICU escalation highlights the IRCU's effectiveness in managing this high-risk population.
Figure 1| Patient characteristics and clinical outcomes in the IRCU (n=249). (A) Age distribution overall (median 51, IQR 38-63) and by gender (Male median 50, IQR 37-60; Female median 55, IQR 40-67). Violin plots show distribution; inner boxes show median/IQRs. (B) Distribution of NIV utilization, ICU transfer, and Exitus for the total cohort, stratified by gender. Bars show percentage (N). Key totals: NIV Yes = 77 (31%), ICU Transfer Yes = 20 (8%), Exitus Yes = 7 (3%). Note potential inconsistency in absolute counts within gender rows for NIV category vs the reported Total N=77 for NIV Yes. (C) Distribution of length of stay (LoS) in the IRCU (days) (median 5, IQR 3-8). Inset: Monthly admissions (Apr-Aug 2021). (D) Correlation matrix heatmap for binary clinical variables. Asterisks denote statistical significance (*** p<0.001, * p<0.05).
Figure 2| Analysis stratified by NIV treatment status and multivariable associations. (A) Age distribution stratified by NIV treatment status (No NIV: n=172, median 42, IQR 35-49; NIV Yes: n=77, median 61, IQR 56-67). (B) Distribution of final, mutually exclusive outcomes stratified by NIV status. For NIV-Yes (n=77): ICU Transfer=23% (n=18), Exitus=9% (n=7), Recovery=68% (n=52). For NIV-No (n=172): Recovery=100% (n=172). Bars show percentage (N). Outcomes defined for mutual exclusivity: ICU Transfer, then Exitus (no prior ICU), then Recovery. (C) Forest plot of adjusted Odds Ratios (OR) with 95% CI from logistic regression (adjusted for Age and Gender where applicable). Ref group for Gender is Female. (D) Forest plot of crude Risk Ratios (RR) with 95% CI comparing risks between groups (Gender [Male vs Female ref.], Age Groups [<50 ref.], NIV Status [Yes vs No ref.]). Note high RRs for NIV vs No NIV reflect confounding by indication. Significance levels: ns (*p*≥0.05), * (p<0.05), *** (p<0.001).
The preceding statistical analysis reveals significant heterogeneity in patient characteristics and clinical trajectories within the UHVN IRCU cohort, highlighting key associations such as older age with Non-Invasive Ventilation (NIV) requirement and markedly different outcomes for patients stratified by NIV status (Section 2.1). While informative, these static analyses offer limited insight into the underlying temporal dynamics governing patient flow, the interplay between admission patterns and care pathways, or the potential impact of changing operational efficiencies. To achieve a more mechanistic understanding and develop a quantitative framework, we adopted a compartmental modeling approach based on Ordinary Differential Equations (ODEs), a technique commonly used in epidemiology and healthcare operations [63].
We model the patient population within the IRCU using a system of coupled ODEs tracking distinct patient states over time t. The state variables represent the number of patients in each compartment: X(t), corresponding to patients admitted to the IRCU not initially requiring NIV; Y (t), representing patients within the IRCU currently receiving NIV; Z(t), denoting the cumulative number of patients transferred from the IRCU to the ICU; W(t), indicating the cumulative number of patients experiencing mortality within the IRCU pathway (without prior ICU transfer from the IRCU); and R(t), which stands for the cumulative number of patients recovered and discharged from the IRCU.
The temporal evolution of these compartments is governed by the following system of ODEs:
$$\frac{dX}{dt} = A(t) - (\alpha + \rho + \eta + \nu)X(t),\tag{1}$$
$$\frac{dY}{dt} = \alpha X(t) - \left[\gamma(1 - \gamma_0) + \theta(t) + \varepsilon(1 - \varepsilon_0)\right] Y(t),\tag{2}$$
$$\frac{dZ}{dt} = \gamma(1 - \gamma_0)Y(t) + \eta X(t),\tag{3}$$
$$\frac{dW}{dt} = \varepsilon(1 - \varepsilon_0)Y(t) + \nu X(t),\tag{4}$$
$$\frac{dR}{dt} = \theta(t)Y(t) + \rho X(t). \tag{5}$$
Here, A(t) is the time-dependent admission rate. For non-NIV patients (X), α is the rate of NIV initiation, ρ is the recovery rate, η is the direct ICU transfer rate, and ν is the mortality rate. For NIV patients (Y), γ is the ICU transfer rate, θ(t) is the potentially time-dependent recovery rate, and ε is the mortality rate (γ0 and ε0 are adjustment factors, assumed zero). Z(t), W(t), and R(t) track cumulative outcomes. We assume Z(t) represents cumulative inflow to ICU, thus setting the outflow rate = 0. The complete derivation of the model, detailed parameter definitions, estimation procedures, and further simulation details are provided in the Supplementary Material.
To ensure the ODE model reflected the specific dynamics of the UHVN cohort, its baseline transition rates were calibrated using observed outcome proportions derived directly from the study data (n=249, Section 2.1). Crucially, based on the clinical finding that patients managed without NIV experienced no ICU transfers or exitus events (Fig. 2B), the corresponding direct transition rates were structurally constrained: η = 0 and ν = 0. For the remaining pathways, the dimensionless outcome proportions were converted into rates [days−1 by incorporating assumed average lengths of stay (LOS) in the non-NIV and active-NIV states (τX = 7 days and τY = 10 days, respectively). These LOS values were predefined parameters used consistently across the simulation analyses described below. This calibration approach provides parameters aligned with the aggregate cohort behavior, a necessary step for developing credible predictive simulations [67]. The detailed methodology for proportion calculation, the LOS assumptions, the resulting rate calculations, and the final parameter values are fully documented in Supplementary Material. Employing this calibrated model structure, we then performed several in silico simulations to explore system dynamics under various critical scenarios, aiming to provide insights beyond the static statistical analysis for resource planning and evaluating the potential impact of changing conditions.
Employing the calibrated ODE model structure (Section 3, Supplementary Material), we performed several in silico simulations to explore the dynamic behavior of the IRCU system under different operational conditions.
First, a baseline autonomous simulation was executed (Simulation A1), utilizing the estimated constant average admission rate (Aavg) and the derived base transition parameters, including the constraint of no direct ICU transfers or mortality from the non-NIV state (η = 0, ν = 0). A constant recovery rate from NIV (θ(t) = *θ*0) was assumed, reflecting the average observed condition. Under these stable parameters, the system evolved towards a steady state (Fig. 3A-C). Panel A illustrates the absolute patient counts, demonstrating that cumulative recoveries (R) represented the predominant outcome pathway. Active patient loads reached maximums of approximately 42 for non-NIV (X) and 18 for NIV (Y) individuals before stabilizing. The total active patient census (X+Y) equilibrated around 60 patients, remaining well below an illustrative capacity threshold of 75 beds (Fig. 3C, dashed red line), with approximately 30.6% requiring NIV in this steady state. At the simulation's end (126 days), the predicted cumulative outcome proportions (relative to total Z+W+R outcomes) were approximately 7.6% for ICU transfer, 1.9% for direct exitus, and 90.6% for recovery (Fig. 3A, text annotations). These simulated proportions exhibit reasonable concordance with the overall cohort outcomes (Fig. 1B), suggesting the baseline parameterization adequately captures the system's aggregate behavior. Correspondingly, the cumulative outcome rates per admission stabilized at approximately 6.9% for ICU, 1.7% for exitus, and 83.2% for recovery (Fig. 3B).
Recognizing the observed fluctuations in monthly admissions (Fig. 1C inset), we subsequently simulated the system's response to a time-varying admission rate, A(t) (Simulation B1, Fig. 3D-F). This scenario utilized a Gaussian function superimposed on the baseline rate Aavg (parameters specified in Fig. 3 caption), while maintaining constant transition parameters (θ(t) = *θ*0). As anticipated, the variable admission rate induced substantial transient dynamics in the active patient populations (Fig. 3D). The peak census of patients requiring NIV (Max Y ≈ 41, Panel D) more than doubled compared to the baseline maximum (Max Y ≈ 18, Panel A). Consequently, the total active patient load reached a markedly higher peak of approximately 138 patients (Fig. 3F), significantly exceeding the illustrative 75-bed capacity threshold and underscoring the potential for resource saturation during admission surges. This simulation vividly illustrates the direct impact of fluctuating admissions on immediate IRCU resource demands, particularly NIV capacity. Notably, while the surge precipitated steeper initial increases in cumulative outcomes (compare Panels D vs. A), the final outcome proportions (7.7% ICU, 1.9% Exitus, 90.4% Recovery, Panel D) remained nearly identical to the baseline scenario. This suggests that, under constant operational efficiency, the primary effect of such an admission surge manifests as transient strain on capacity rather than a fundamental shift in the long-term distribution of patient outcomes.
We then explored the potential impact of temporally varying recovery efficiency for patients receiving NIV (Simulation B2, Fig. 3G-I). Motivated by the IRCU's observed success in managing a majority (68%) of high-risk NIV patients without ICU escalation (Fig. 2B), we modeled a scenario where recovery efficiency improves over time using θ(t) = θ0 + ∆θe−λt, while maintaining a constant average admission rate (Aavg). Compared to the baseline simulation (A1), this transiently enhanced recovery function yielded a discernible improvement in final outcomes (Fig. 3G). Specifically, it resulted in modest reductions in predicted cumulative ICU transfers (approx. 6.8% vs. 7.6% in A1) and direct exitus events (approx. 1.7% vs. 1.9% in A1), with a corresponding increase in cumulative recoveries (approx. 91.5% vs. 90.6% in A1). The cumulative outcome rates per admission (Fig. 3H) reflected this favorable shift (e.g., Final Rate ICU / Adm. ≈ 6.2% vs. 6.9%). Given the constant admission rate, the active patient load dynamics (Fig. 3I) closely resembled the baseline scenario, stabilizing at approximately 60 patients (below the illustrative capacity) with around 30.5% requiring NIV, thereby demonstrating the potential benefits of improved care efficiency under stable admission pressure.
Finally, to synthesize the most plausible dynamic representation incorporating both the observed admission peak and potential efficiency gains, we simulated a scenario combining the variable admission rate A(t) with the variable recovery rate θ(t) (Simulation B3, Fig. 3J-L). This configuration examines the interplay between increased admission pressure and concurrent improvements in care delivery. The results indicate that while the variable A(t) profile continues to drive significant peaks in active patient load (Max X+Y ≈ 136, Panel L) and NIV requirement (Max Y ≈ 39, Panel J), closely mirroring Simulation B1 (Fig. 3F) and substantially surpassing the illustrative capacity limit, the concurrently enhanced recovery efficiency (θ(t)) provides partial mitigation of the adverse downstream consequences. Compared to the surge scenario with constant efficiency (B1, Panel D), this combined scenario (B3, Panel J) resulted in slightly lower final proportions for ICU transfers (approx. 7.0% vs. 7.7%) and direct mortality (approx. 1.8% vs. 1.9%), and higher recoveries (approx. 91.2% vs. 90.4%). This suggests that modeled improvements in care efficiency can modestly counteract the negative shift
in outcome distributions associated with admission surges, although they do not substantially alleviate the peak resource strain (patient volume) imposed by the surge itself.
Figure 3| Model Simulation Results under Different Scenarios. Dynamics of patient compartments (X: No NIV, Y: NIV, Z: Cumulative ICU, W: Cumulative Exitus, R: Cumulative Recovered) over 126 days for four different simulation scenarios based on UHVN IRCU data. Column 1 (A, D, G, J): Absolute patient numbers, final outcome proportions (calculated vs. total Z+W+R outcomes at day 126, shown as percentages), and maximum integer patient counts observed in X and Y states (e.g., Max Y ≈ 18). Column 2 (B, E, H, K): Cumulative outcome rates (Z, W, R) normalized by cumulative admissions (%). Column 3 (C, F, I, L): Total active patient load (X+Y, left axis, purple line) and percentage of active patients requiring NIV (right axis, blue line), with annotations for peak or final/stable values. The dashed horizontal line in panels C, F, I, L indicates an illustrative unit capacity of 75 beds. Base parameters (Aavg, α, ρ, γ, ε, *θ*0) were estimated from cohort data (details in Suppl. Section S5). Rates η (X to ICU) and ν (X to Exitus) were fixed at 0 based on observations. (A-C) Scenario A1 (Baseline): Constant admission rate A(t) = Aavg; constant recovery rate θ(t) = *θ*0. Represents average, stable conditions. (D-F) Scenario B1 (Variable Admissions): Variable A(t) modeled as a Gaussian peak (center = 65, width = 15, scale = 1*.*5) superimposed on Aavg, simulating an admission surge; constant recovery rate θ(t) = *θ*0. (G-I) Scenario B2 (Variable Recovery): Constant admission rate A(t) = Aavg; variable recovery rate θ(t) = θ0 + 0*.*15e −0*.05t* modeling potential efficiency improvement. (J-L) Scenario B3 (Combined): Variable A(t) as in B1; variable θ(t) as in B2. Simulates surge concurrent with efficiency gains.
Finally, sensitivity analyses were performed on the autonomous model (constant Aavg, θ(t) = *θ*0) to assess the relative importance of key transition rates on both final outcomes and peak system load (Fig. 4). We varied the four main transition rates (α, *θ*0, γ, ε) individually by ±30% around their estimated baseline values. These baseline values were α ≈ 0*.*044, θ0 ≈ 0*.*068, γ ≈ 0*.*023, and ε ≈ 0*.*009 days−1 , resulting in reference model predictions of approximately 13 cumulative ICU transfers (Z), 52 cumulative exitus events (W), 622 cumulative recoveries (R), and a peak NIV load (Max. Y) of approximately 18 patients.
Varying the NIV initiation rate from non-NIV patients (α, Transition X to Y) revealed a direct impact on both outcomes and load (Fig. 4A, E). Increasing α led to higher final cumulative counts for both ICU transfers (Z) and exitus events (W) compared to their baseline values (Z ≈ 13, W ≈ 52 at baseline α), alongside a less pronounced decrease in recoveries (R, baseline R ≈ 622) (Panel A). Concurrently, and significantly, the maximum number of patients requiring NIV ('Max. NIV Patients (Y)') increased substantially with α, rising from below 15 to above 19 patients across the tested range, with a baseline value of approximately 18 patients at the baseline α value (Panel E). This suggests that factors accelerating patient progression towards needing NIV (higher α) significantly strain downstream resources (ICU, mortality pathways) and increase the peak demand for NIV beds.
Conversely, varying the base recovery rate from NIV (*θ*0, Transition Y to R) demonstrated a strong influence primarily on outcome distribution and NIV load (Fig. 4B, F). Increasing θ0 markedly decreased the final cumulative counts of ICU transfers (Z, baseline Z ≈ 13) and exitus events (W, baseline W ≈ 52), while substantially increasing recoveries (R, baseline R ≈ 622) (Panel B, shown at baseline θ0 ≈ 0*.*068). Furthermore, a higher recovery rate significantly reduced the maximum number of patients requiring NIV, decreasing from over 22 to below 16 patients across the tested range, with a baseline peak load of approximately 18 patients at the baseline θ0 value (Panel F). This pronounced sensitivity underscores the critical role of effective NIV management and recovery pathways within the IRCU, consistent with the observed high recovery rate (68%) among NIV patients in our actual cohort (Fig. 2B).
Variations in the ICU transfer rate from NIV (γ, Transition Y to Z) primarily affected the ICU outcome, as expected (Fig. 4C, G). Increasing γ led to a strong, near-linear increase in the final cumulative ICU count (Z) starting from the baseline of Z ≈ 13 (at baseline γ ≈ 0*.*023), accompanied by slight decreases in exitus (W, baseline W ≈ 52) and recoveries (R, baseline R ≈ 622) (Panel C). Interestingly, increasing the transfer rate γ slightly decreased the maximum NIV patient load from approximately 20 down to below 18 patients across the tested range, with the baseline value being approximately 18 patients (Panel G), potentially because patients are moved out of the Y state more quickly.
Lastly, varying the direct exitus rate from NIV (ε, Transition Y to W) showed a strong, direct relationship with cumulative mortality (Fig. 4D, H). Increasing ε resulted in a substantial rise in the final exitus count (W) from the baseline of W ≈ 52 (at baseline ε ≈ 0*.*009), with corresponding slight decreases in ICU transfers (Z, baseline Z ≈ 13) and recoveries (R, baseline R ≈ 622) (Panel D). Notably, this parameter had a minimal impact on the maximum NIV patient load across the range tested, which remained close to the baseline value of approximately 18 patients (Panel H), indicating that while crucial for the mortality outcome, its variation within this range doesn't significantly alter the peak bed demand for NIV itself.
Overall, these sensitivity analyses highlight that while all tested parameters influence outcomes, the NIV initiation rate (α) and the NIV recovery rate (*θ*0) appear to have the most substantial impact on the peak operational load (maximum NIV patients, baseline Max. Y ≈ 18). Interventions targeting factors related to α (e.g., earlier treatments preventing deterioration) or θ0 (e.g., enhancing NIV efficacy) could therefore have the largest effects on IRCU resource management, while γ and ε are more directly related to specific downstream outcomes (ICU use, baseline Z ≈ 13, and mortality, baseline W ≈ 52, respectively).
Figure 4| Sensitivity Analysis of the Autonomous Model. Impact of varying key transition rates (α, *θ*0, γ, ε) by ±30% around their baseline values (indicated by markers and dashed lines corresponding to Z ≈ 13, W ≈ 52, R ≈ 622, Max. Y ≈ 18) on model outputs. Baseline parameter values used for reference points are α ≈ 0*.*044, θ0 ≈ 0*.*068, γ ≈ 0*.*023, ε ≈ 0*.*009 days−1 . Top Row (A-D): Absolute final cumulative outcome counts for ICU transfer (Z, yellow, left axis), Exitus (W, green, left axis), and Recovery (R, dark blue/purple, right axis). Bottom Row (E-H): Maximum number of patients simultaneously requiring NIV (Max. NIV Patients (Y)). Parameters varied: (A, E) NIV initiation rate (α, X to Y); (B, F) NIV recovery rate (*θ*0, Y to R); (C, G) ICU transfer rate (γ, Y to Z); (D, H) Direct exitus rate (ε, Y to W). Analyses performed using the autonomous model configuration (constant Aavg, constant θ(t) = θ0 when varying other parameters).
Beyond the state-transition dynamics explored via the ODE model, the empirical distribution of the total Length-of-Stay (LOS) within the IRCU provides an alternative lens through which to analyze system load. Figure 5A-C characterizes this distribution. The observed frequency data (Fig. 5A, points) highlights the variability and right-skewness typical of LOS data, which is captured by a Gaussian Process (GP) fit (Fig. 5A, line and shaded area). The corresponding empirical Probability Mass Function (PMF), f(d) = P(LOS = d), reveals a high probability for short stays, with the mode occurring within the first few days (Fig. 5B). From this distribution, the empirical mean LOS was calculated as 9.0 days. The derived empirical survival function, S(d) = P(LOS ≥ d), representing the probability that a patient remains within the IRCU for at least d days, is shown in Fig. 5C (solid line) alongside smoothed versions (LOWESS and GP-derived). This empirical S(d) inherently captures the overall departure pattern observed in the cohort and serves as the core input for occupancy modeling.
This empirical survival function S(d) enables a direct estimation of the total daily IRCU occupancy, Xtotal(t), via convolution with the daily admission sequence A(t), based on the fundamental relationship Xtotal(t) = Pdmax k=0 A(t − k)S(k) (see Supplementary Material for a detailed mathematical description of the underlying methods). We applied this forward convolution method using the empirical S(d) and illustrative admission scenarios. Figure 5D displays the predicted occupancy (right axis, pink line) resulting from a representative surge admission pattern A(t) (left axis, light green line), which has a baseline rate of 5 admissions/day and a peak surge reaching 15 admissions/day around day 20. The predicted total occupancy demonstrates a characteristic lagged and smoothed response, reaching a peak occupancy of approximately 124 beds around day 24 before declining, considerably exceeding the expected steady-state occupancy of approx. 45 beds under baseline admissions.
To further explore operational implications, Figure 5E compares the predicted occupancy under different admission patterns against an illustrative capacity threshold (approx. 63 beds, red line), derived as 125% of the baseline steady-state occupancy. While the baseline admission rate (5/day, blue dashed line) keeps the occupancy well below capacity, both the stress scenario (sustained 10 admissions/day, yellow dash-dot line) and the surge scenario (red dashed line) clearly breach this threshold for extended periods, highlighting the vulnerability to increased admission pressure given the observed LOS distribution.
Finally, we investigated the sensitivity of occupancy to overall system efficiency, modeled via changes in the Mean LOS achieved (Fig. 5F). Using the baseline constant admission rate (5/day), simulations were run assuming the baseline empirical Mean LOS (9d, solid blue line), a 'Low Efficiency' scenario resulting in a longer Mean LOS (≈11d, dotted orange line), and a 'High Efficiency' scenario achieving a shorter Mean LOS (≈7d, dashed green line). As shown, improved efficiency (shorter LOS) substantially reduces the required bed occupancy, yielding an average reduction of approximately 16 beds (33%) compared to baseline and 20 beds (37%) compared to the low efficiency scenario during the steady-state period (Panel F, text box).
This convolution-based total occupancy modeling complements the ODE model results. While the ODE simulations provide insights into state-specific occupancy (X vs Y ), the convolution model offers a direct prediction of total bed demand grounded in the aggregate observed LOS. Comparing the total peak occupancy predicted here (e.g., approx. 124 in Panel D) with the peak active load from the corresponding ODE simulation (e.g., approx. 136 in Fig. 3L for B3) reveals broadly consistent results between the two distinct modeling approaches, lending confidence to the overall understanding of system dynamics. The empirical PMF f(d) (Panel B) also represents the necessary input for potential future work using deconvolution techniques to estimate historical admissions from occupancy data.
Figure 5| Length-of-Stay Analysis and Occupancy Simulations. (A) Gaussian Process (GP) regression fit to observed IRCU Length-of-Stay (LOS) frequency data (counts of patients discharged per day of stay). Points are observed counts; dark green line is GP mean fit; shaded area is GP mean ± 2 standard deviations. (B) Empirical LOS Probability Mass Function (PMF), P(LOS = d), derived from observed data. (C) Comparison of LOS survival functions, S(d) = P(LOS ≥ d): Empirical (derived from data), LOWESS smoothed (fraction f=0.2), and GP-derived (from the integral of the GP mean fit PMF). The Empirical S(d) is used for convolution simulations. (D) Scenario 1: Predicted total IRCU occupancy (pink, right axis) resulting from forward convolution of an illustrative surge admission pattern (light green, left axis) with the empirical S(d). Annotations indicate baseline and peak admission rates, steady-state occupancy, and peak predicted occupancy with corresponding day. (E) Scenario 3: Comparison of predicted total occupancy under different admission scenarios (Baseline: 5/day; Stress: 10/day sustained; Surge: peak admission rate 15/day) against an illustrative capacity threshold (≈63 beds). Uses empirical S(d). Legend reflects simplified labels used in the plot. (F) Scenario 4: Sensitivity of predicted total occupancy to system efficiency, modeled via different Mean LOS values resulting from varying effectiveness levels (Low Eff: Mean LOS≈11d; Baseline: Mean LOS=9d; High Eff: Mean LOS≈7d). Assumes constant baseline admissions (5/day). Text box quantifies average steady-state occupancy reductions achieved by high efficiency.
Comparing the modeling approaches, the peak total occupancy predictions under surge conditions show broad consistency. The LOS/convolution model predicted a peak of approximately 124 patients (Fig. 5D), while the corresponding ODE simulation (B3) yielded a peak active load (X+Y) of approximately 136 patients (Fig. 3L). This agreement between distinct methodologies lends confidence to the overall assessment of system dynamics under stress. Furthermore, the two models offer complementary perspectives: while the convolution model provides a direct forecast of total occupancy grounded in the empirical LOS distribution, the ODE model elucidates the internal dynamics, revealing the fluctuating distribution between NIV and non-NIV states and the impact of varying transition rates, factors not explicitly visible in the aggregate LOS analysis.
This study provides critical insights into the operational dynamics and clinical effectiveness of an Intermediate Respiratory Care Unit (IRCU) established during the exigencies of the COVID-19 pandemic, with a particular focus on the operational context involving specialized clinical staffing. Our primary findings demonstrate that this newly implemented IRCU successfully managed a substantial cohort of patients with severe COVID-19 associated respiratory failure, many requiring non-invasive ventilation (NIV). Strikingly, while patients requiring NIV were significantly older and represented a higher-risk stratum, 68% of these individuals recovered within the IRCU without needing escalation to intensive care. Furthermore, no ICU transfers or mortality events occurred among patients managed without NIV within the IRCU pathway. These outcomes, achieved during a period of immense healthcare strain, underscore the potential of IRCUs, particularly when structured with key elements such as dedicated, expert clinical teams, to function as highly effective units for delivering complex respiratory care, thereby alleviating critical pressure on ICU resources. Our integrated statistical and dynamic modeling approach further elucidates the quantitative factors governing patient flow and the potential impact of operational changes.
Our findings resonate with the burgeoning literature highlighting the strategic importance of IRCUs as a bridge between the general ward and the ICU, particularly during the COVID-19 crisis [22, 13, 23, 68]. Several studies have reported favorable outcomes and reduced ICU burden associated with IRCU implementation during the pandemic [24, 33]. The observed 8% overall ICU transfer rate in our cohort aligns broadly with figures reported elsewhere, although direct comparisons are complex due to heterogeneity in patient populations, admission criteria, and unit capabilities across centers [22]. A unique aspect of our work is its setting within an IRCU explicitly designed with dedicated 24/7 respiratory physician coverage and nursing staff receiving specialized training (Fig. 6, 7) [36, 53]. While other reports focus primarily on treatment modalities or overall patient outcomes [69], our study provides circumstantial evidence suggesting that this combination of specialized medical and nursing skills, manifested in managing NIV/HFNO, early recognition of deterioration, and protocol adherence, was likely an important factor associated with the successful recovery of the majority of our high-risk NIV patients — recognizing that this observational study establishes association, not direct causation. This is consistent with established links between clinical expertise (both nursing and medical) and improved patient outcomes in critical and intermediate care settings [39, 40, 43].
The specific nurse-to-patient ratio maintained in the UHVN IRCU (approximately 1:4) is pertinent here, aligning with levels often deemed appropriate for intermediate care but subject to considerable variation [70]. Contextualizing this against the robust evidence linking staffing levels to clinical results is crucial. Foundational research and subsequent meta-analyses have repeatedly shown that lower nurse-to-patient ratios are strongly associated with better outcomes, including reduced mortality and fewer adverse events in acute care [39, 40, 42, 41]. Our positive findings, notably the high NIV success rate, appear consonant with this principle, suggesting the dedicated staffing model likely contributed significantly. Conversely, higher patient assignments per nurse are frequently linked to increased risks [71]. Yet, staffing ratios operate within a complex interplay of factors. The educational background and specialized skill set of the nursing staff are also independent predictors of patient outcome [72]. Similarly, the continuous presence of specialized physicians (intensivists or pulmonologists) is associated with better ICU outcomes [43, 38]. It is therefore highly plausible that the targeted respiratory care training received by the UHVN IRCU nurses (Fig. 6, Fig. 7), coupled with the dedicated physician coverage, acted synergistically with the 1:4 nursing ratio, enhancing the team's capacity to manage complex respiratory failure effectively — although disentangling these specific contributions requires further research. Disentangling the specific contributions of staffing numbers versus the specialized competencies of different clinical disciplines within IRCUs remains a key area for future investigation. Furthermore, the high recovery rate among NIV patients reinforces the capacity of IRCUs to prevent substantial numbers of ICU admissions, aligning with pre-pandemic estimates and pandemic observations [7, 73].
The theoretical implications of our study extend to the conceptualization of IRCUs not merely as holding bays, but as dynamic, specialized environments where targeted interventions, potentially facilitated by expert staffing, can significantly alter patient trajectories. Our compartmental modeling (Fig. 3, 4) provides a quantitative framework that moves beyond static outcome reporting. It allows for the dissection of patient flow dynamics, identifying critical bottlenecks and leverage points [74, 58]. Specifically, the high sensitivity of both outcomes and peak NIV occupancy to the NIV initiation rate (α) and, crucially, the NIV recovery rate (*θ*0), offers a mechanistic basis for understanding the impact of clinical efficiency — an efficiency likely enhanced by skilled clinical teams. Interventions enhancing recovery (e.g., optimized NIV protocols, proactive weaning, expert clinical care) demonstrated the potential to substantially reduce adverse outcomes and resource strain, even during admission surges. Practically, our findings support the consideration of wider adoption and adequate
resourcing of IRCUs that incorporate dedicated and specialized clinical staffing models [13, 36, 43]. — While our study cannot prove causality, such efficiency is plausibly supported or enhanced by the presence of skilled clinical teams. — The developed modeling framework (Fig. 3, 5, Fig. S1) offers a valuable tool for hospital administrators and clinical leads for in silico testing of operational strategies, predicting resource needs under varying admission scenarios (e.g., surges), and evaluating the potential return on investment for initiatives aimed at improving care efficiency (e.g., advanced clinical training, staffing adjustments) [61, 57, 75]. The concordance between the ODE and convolution-based occupancy predictions (Section 2.3) lends further credence to the robustness of this modeling approach for operational planning.
Several limitations warrant consideration when interpreting these results. First, as a single-center study conducted over a specific five-month period during the COVID-19 pandemic, the generalizability of our quantitative findings to other settings, patient populations (e.g., non-COVID respiratory failure), or different pandemic phases requires confirmation. Although data were collected prospectively, the analysis is retrospective in nature, limiting our ability to definitively establish causality between the specific IRCU staffing model and patient outcomes; the observed associations, while strong, remain associative. Second, while we captured key outcomes and NIV use, finer-grained data on specific ventilator settings, precise HFNC parameters, daily physiological trajectories, or detailed comorbidity indices were not the focus of this analysis and could offer deeper insights. Third, the modeling approach, while informative, relies on assumptions regarding transition rates and their potential time dependencies. Parameter estimation was based on aggregate cohort data, and while sensitivity analyses explored robustness (Fig. 4), the models simplify complex individual patient pathways. The observed lack of adverse outcomes in the non-NIV group, while empirically accurate for this cohort, necessitated structural assumptions (η = 0, ν = 0) in the model that might differ in populations with higher baseline acuity or different pathologies. Finally, the post hoc power analysis indicated moderate power (63%) to detect medium-sized effects, suggesting that subtle differences, particularly regarding gender or less common outcomes, might have been missed.
Our investigation illuminates the critical operational dynamics and clinical utility of a specialized Intermediate Respiratory Care Unit, particularly highlighting the positive outcomes observed within a model emphasizing dedicated and expert clinical staffing, during the unprecedented challenge of the COVID-19 pandemic. We demonstrated that this newly established unit, featuring specialized physician and nursing support, successfully managed a significant cohort of patients presenting with severe COVID-19 respiratory failure. A central finding is the successful recovery of 68% of patients requiring non-invasive ventilation entirely within the IRCU, obviating the need for ICU escalation in the majority of these high-risk individuals. This outcome, coupled with the notable absence of ICU transfers or mortality among those patients managed without NIV during their IRCU stay, powerfully demonstrates the capacity of intermediate care settings when appropriately staffed and structured to function not merely as observational waypoints, but as active therapeutic environments. These units can effectively absorb significant clinical burden and mitigate pressure on finite critical care resources, a capability of paramount importance during public health crises.
These results significantly extend the growing body of evidence supporting IRCU efficacy [22, 13, 23], particularly in the pandemic context [24]. However, our work distinctively foregrounds the association between the implementation of a specialized staffing model (including targeted nursing competencies and dedicated physician coverage) and these favorable outcomes, aligning with broader findings on the impact of clinical expertise and appropriate staffing levels [39, 36, 40, 43]. Furthermore, the integration of statistical analysis with dynamic systems modeling provides a novel quantitative lens through which to understand patient flow. This approach reveals the high sensitivity of system throughput and, consequently, patient outcomes, to parameters directly influenced by care quality and efficiency—an efficiency likely enhanced by skilled multidisciplinary teams—such as the rate of NIV initiation (α) and the rate of recovery from NIV (*θ*0) (Fig. 4). This combined clinical and modeling strategy offers a more mechanistic understanding of IRCU function than previously available, moving beyond static outcome descriptions towards a framework for operational optimization and in silico evaluation of interventions.
While illuminating, our findings naturally open several promising avenues for future inquiry, building upon the insights gained and addressing the study's limitations. First, large-scale, multi-center prospective studies are imperative to validate the generalizability of our findings across diverse healthcare systems, patient populations (including non-COVID respiratory conditions), and varying phases of pandemic activity. Such studies should explicitly incorporate robust metrics of overall clinical staffing levels, skill mix (for both physicians and nurses), specific training programs, and team functioning to better elucidate the relative contributions of these factors. Second, delving deeper into the physiological mechanisms underpinning successful NIV management and weaning within the IRCU setting, perhaps through enhanced physiological monitoring protocols or comparative effectiveness trials of different NIV modes and interfaces [76, 69], will be crucial for refining clinical practice. Third, the predictive power of our dynamic models (Section 2.2) could be substantially enhanced by incorporating individual patient-level covariates, stochastic elements reflecting real-world variability, and potentially adaptive real-time data assimilation frameworks, perhaps leveraging machine learning approaches [63, 77]. Fourth, rigorous longitudinal studies are needed to track the long-term functional outcomes, health-related quality of life, and healthcare utilization patterns of patients successfully managed within IRCUs compared to those whose trajectory necessitates initial or subsequent ICU admission [78]. Finally, exploring the systematic adaptation and evaluation of IRCUs structured with specialized multidisciplinary teams for the broad spectrum of acute respiratory failure etiologies beyond viral pandemics represents a critical area for future investigation to solidify their role within the permanent healthcare infrastructure. Addressing these questions will further elucidate the full potential and optimal implementation strategies for Intermediate Respiratory Care Units in contemporary healthcare delivery.
This prospective cohort study was conducted within the IRCU at the University Hospital Virgen de las Nieves (UHVN), a tertiary academic medical center in Granada, Spain. The establishment of this IRCU commenced in February 2020, motivated by the anticipated demands of the nascent COVID-19 pandemic (Fig. 6, point 01). A foundational principle of the unit's design was the implementation of a dedicated, specialized clinical staffing model, ensuring continuous expert coverage. Preparatory phases included targeted recruitment and, notably, specialized training workshops for the nursing contingent focusing on advanced respiratory support modalities (Non-Invasive Ventilation [NIV] and High-Flow Nasal Oxygen [HFNO]), initiated in December 2020 (Fig. 6, point 02).
The IRCU, formally accredited by the Spanish Society of Pneumology and Thoracic Surgery (SEPAR) at its highest complexity level, became operational in April 2021 (Fig. 6, point 03). It functioned with dedicated medical staffing, providing 24/7 in-house coverage by respiratory physicians, and specialized nursing personnel. The unit was comprehensively equipped with advanced respiratory support and monitoring technology. Ethical approval for the overarching project involving COVID-19 patients at the institution, under which the data for this specific IRCU analysis were collected, was granted by the Provincial Research Ethics Committee of Granada (CEIM/CEI Provincial de Granada) on June 24, 2020. Written informed consent was obtained from all participants or their legally authorized representatives, in accordance with the Declaration of Helsinki, pertinent data protection regulations, and the provisions of this ethical approval. The period designated for patient inclusion and prospective data collection for the analysis presented herein extended from April 2021 through August 2021 (Fig. 6, points 03-04).
Figure 6| Timeline of UHVN IRCU Implementation and Study Execution. Key milestones include: (01) Feb 2020: Project initiation driven by COVID-19 emergence. (02) Dec 2020: Clinical staff preparation, including nursing staff recruitment and specialized training (NIV/HFNO). (03) Apr 2021: Official IRCU opening and commencement of data collection. (04) Aug 2021: End of the patient inclusion period, encompassing 249 patients who received respiratory support (CPAP, BiPAP, HFNO). (05) Sep 2021: Initiation of data cleaning and statistical analysis.
Eligible participants were adult patients (age ≥ 18 years) admitted to the UHVN IRCU between April and August 2021 (Fig. 6, points 03-04). The primary admission criterion was severe respiratory failure, predominantly associated with confirmed SARS-CoV-2 infection and compatible radiological findings. Patients admitted principally for non-pulmonary primary diagnoses were excluded. While other potential exclusion criteria, such as pre-existing directives limiting therapeutic escalation (e.g., 'do not resuscitate' orders), were considered in individual patient management decisions, they were not systematically applied as absolute exclusion criteria for enrollment in this observational cohort. The final analytic sample comprised 249 consecutive patients receiving care within the IRCU during the defined study interval (Figs. 1, 2).
Upon IRCU admission, patients underwent a standardized initial assessment protocol (Fig. 7, Reception & Evaluation phases), encompassing vital sign measurement and comprehensive clinical evaluation by the dedicated medical and nursing staff. Patients fulfilling criteria for severe respiratory failure (typically defined by an inspired oxygen fraction [FiO2] requirement exceeding 40%) were initiated on NIV or HFNO based upon collaborative clinical judgment (Fig. 7, Treatment phase). The principal respiratory support modalities employed included Continuous Positive Airway Pressure (CPAP) [7, 22, 13], Bilevel Positive Airway Pressure (BiPAP) via Respironics V60 ventilators (Philips) [11, 13], and HFNO [13, 12]. The utilization of HFNO adhered to evolving international guidance during the pandemic, incorporating WHO recommendations for vigilant monitoring [79, 80]. To ensure confidentiality as per ethical approval, patient data were de-identified prior to analysis by the research team.
Continuous cardiorespiratory monitoring informed ongoing clinical management. A key institutional guideline mandated proactive multidisciplinary consultation, particularly involving the Intensive Care Unit (ICU) team, for patients exhibiting sustained high oxygen requirements (e.g., FiO2 > 80% despite NIV optimization) to facilitate timely decisions regarding the appropriate level of care [11, 7]. This collaborative framework guided patient trajectory determinations (Fig. 7, Development phase), encompassing potential ICU escalation, transfer to a general medical ward upon clinical improvement, or direct discharge home. Throughout the IRCU admission, established nurse-to-patient (approximately 1:4) and physician-to-patient (approximately 1:6) ratios were maintained, ensuring continuous specialized clinical oversight. Data were collected prospectively using standardized forms, capturing patient demographics, pertinent comorbidities, respiratory support modalities and key parameters (when available), length of stay (LOS), and definitive IRCU pathway outcomes (ICU transfer, in-IRCU exitus, or recovery/discharge), constituting the dataset for subsequent analysis (Figs. 1, 2).
Figure 7| Patient Management Workflow and Key Clinical Team Actions in the UHVN IRCU. The process encompasses four main stages: (1) Reception: Initial patient identification and proactive recognition. (2) Evaluation: Collection of vital signs and integral clinical assessment (including nursing). (3) Treatment: Selection and administration of respiratory support (NIV: CPAP, BiPAP, HFNO) with continuous multidisciplinary monitoring. (4) Development: Monitoring patient trajectory towards outcomes (Recovery: transfer to General Ward or Discharge; Worsening: transfer to ICU; or Exitus).
Descriptive statistics, including counts (N), percentages (%), medians, and interquartile ranges (IQR), were employed to summarize the cohort's baseline characteristics and clinical outcomes. Distributions of continuous variables were visualized using violin plots (Fig. 1). Intergroup comparisons utilized the Mann-Whitney U test for continuous variables (e.g., Fig. 1A, 2A) and Chi-squared or Fisher's exact tests, as appropriate, for categorical variables (underlying data for Fig. 1B,D, 2B,D). Associations between predictor variables (age, gender, NIV status) and binary clinical outcomes (NIV utilization, ICU transfer, exitus) were evaluated using both crude Risk Ratios (RR) with 95% confidence intervals (CI) (Fig. 2D) and multivariable logistic regression models. The logistic regression analyses yielded adjusted Odds Ratios (OR) with 95% CIs (Fig. 2C), controlling for potential confounding by patient age and gender where applicable. Statistical significance was established at a two-tailed p-value threshold of < 0.05. Data cleaning and formal statistical analysis commenced in September 2021.
The distribution of total Length-of-Stay (LOS) within the IRCU (Fig. 5) was characterized empirically using Probability Mass Functions (PMF) and Kaplan-Meier-derived Survival Functions (S(d)). To facilitate smoothed estimates and potential forecasting applications, Gaussian Process (GP) regression [77], a non-parametric Bayesian technique (methodological details provided in Supplementary Material Section S4), was fitted to the observed LOS frequency data (Fig. 5A), yielding smoothed survival function estimates (Fig. 5C). However, the occupancy simulations presented herein primarily utilized the directly calculated empirical S(d). All statistical computations were performed using Python (version 3.7; libraries: SciPy, Statsmodels, Matplotlib, Scikit-learn) and R (version 4.3.2; e.g., 'pwr' package). Figures, including timeline and workflow diagrams (6, 7) and conceptual model schematics in the Supplement, were prepared using Inkscape (version 1.3) [81]. Given the minimal extent of missing data (<1%), a complete case analysis approach was deemed appropriate.
To contextualize the statistical power of our primary group comparisons (NIV vs. Non-NIV strata), a post hoc sensitivity analysis was conducted using G*Power (version 3.1) [82]. It is imperative to note that this analysis assesses the study's sensitivity to detect a pre-specified, conventionally defined effect size, rather than calculating power based on the observed effect magnitudes. Specifically, we determined the achieved statistical power (1-β) of the study, given the final cohort sample sizes (Total N=249; Non-NIV n1 = 172, NIV n2 = 77), to identify a statistically significant difference corresponding to a Cohen's d of 0.5 between these groups for a continuous outcome. This effect size (d=0.5) was selected based on Cohen's established convention for a 'medium' effect within the biomedical sciences [83]. The analysis assumed a two-tailed independent samples t-test framework, suitable for such group comparisons, employing the standard alpha level (α) of 0.05. The resulting achieved power was calculated to be approximately 0.63 (63%). This signifies moderate sensitivity: under the realized sample distribution, there was a 63% probability of detecting a true difference of medium magnitude (d=0.5), should one exist in the population, as statistically significant at the p<0.05 level. Conversely, this implies a non-trivial risk (Type II error probability, β ≈ 0*.*37) of failing to detect such a medium-sized effect. This moderate statistical power warrants caution in interpreting non-significant findings for comparisons between the NIV and non-NIV groups and underscores the value of future, potentially larger studies to investigate more
The authors affirm that there are no competing interests—financial, personal, or otherwise—that have influenced, or could be perceived to have influenced, the conduct or reporting of the research described in this paper.
The study adhered to the Declaration of Helsinki. Ethical approval for the broader project encompassing the collection of data analyzed in this manuscript was provided by the Provincial Research Ethics Committee of Granada (CEIM/CEI Provincial de Granada) on June 24, 2020. Written informed consent was obtained from all participants or their legal representatives prior to inclusion, following the stipulations of the approved protocol. Participant data confidentiality was maintained throughout the study.
We are truly grateful to every patient who generously dedicated their time and effort to be part of our study. We also want to convey our deep appreciation to the exceptional medical and nursing staff of the IRCU at the UHVN in Granada, Spain. Their unwavering dedication, tireless efforts, and unparalleled expertise have been instrumental in maintaining the seamless functioning of the IRCU and providing outstanding care to our patients. Their remarkable teamwork and exceptional skills have been indispensable in delivering our patients the utmost quality of care. RF acknowledges partial support from the María Zambrano-Senior grant (Spanish Ministerio de Universidades and Next-Generation EU); Grant C-EXP-265-UGR23 funded by Consejería de Universidad, Investigación & Innovación & ERDF/EU Andalusia Program; Grant PID2022-137228OB-I00 funded by the Spanish Ministerio de Ciencia, Innovación y Universidades, MICIU/AEI/10.13039/501100011033 & "ERDF/EU A way of making Europe"; and the Modeling Nature Research Unit, project QUAL21-011.
The Python scripts used for the ODE modeling (Figures 3-4), Length-of-Stay analysis, and occupancy simulations (Figure 5), along with the necessary input data files ("patient_data.csv", "ircu_data.csv"), are publicly available. The specific version of the code and data corresponding to the results presented in this manuscript (Version v1.0.0) has been archived on Zenodo and can be accessed via the following DOI: https://doi.org/10.5281/ zenodo.15286823. The development repository is hosted on GitHub at https://github.com/renee29/IRCU\_ Patient\_Flow\_Modeling\_Scripts.
A.C.N.O., C.M.G., and R.F. conceptualized the study. A.C.N.O., J.A.S.M., and C.M.G. developed the clinical methodology, while R.F. developed the modeling and statistical methodology. Clinical investigation and data curation were performed by A.C.N.O., P.G.F., and J.A.S.M. Formal analysis, software, and visualization were conducted by R.F. Project administration was by A.C.N.O.; resources and clinical supervision by C.M.G. and A.C.N.O.; analytical supervision by R.F. Validation involved A.C.N.O., J.A.S.M., C.M.G., and R.F. The original draft was written by A.C.N.O., J.A.S.M., C.M.G., and R.F. All authors reviewed and edited the final manuscript.
[3] J. S. Whittle, I. Pavlov, A. D. Sacchetti, et al. "Respiratory support for adult patients with COVID-19". In: JACEP Open 1.2 (2020), pp. 95–101. doi: 10.1002/emp2.12071.
[4] M. A. Matthay, R. L. Zemans, G. A. Zimmerman, et al. "Acute respiratory distress syndrome". In: Nature Reviews Disease Primers 5.1 (2019), p. 18. doi: 10.1038/s41572-019-0069-0.
[5] T. L. Petty, D. B. Bigelow, and L. M. Nett. "The Intensive Respiratory Care Unit—An Approach to the Care of Acute Respiratory Failure". In: California Medicine 107.5 (1967), pp. 381–384.
[6] T. L. Petty, S. Lakshminarayan, S. A. Sahn, et al. "Intensive respiratory care unit. Review of ten years' experience". In: JAMA 233.1 (1975), pp. 34–37. doi: 10 . 1001 / jama . 1975 . 03260010036016.
[7] A. Torres, M. S. Niederman, R. Celis, et al. "Intermediate respiratory care units: A European perspective". In: Archivos de Bronconeumología (English Edition) 41 (2005), pp. 505–512. doi: 10.1016/s1579-2129(06)60271-1.
[8] S. A. Nasraway, I. L. Cohen, R. C. Dennis, et al. "Guidelines on admission and discharge for adult intermediate care units". In: Critical Care Medicine 26.3 (1998), pp. 607–610. doi: 10.1097/00003246-199803000-00039.
[9] D. C. Cheng, R. J. Byrick, and E. Knobel. "Structural models for intermediate care areas". In: Critical Care Medicine 27.10 (1999), pp. 2266–2271. doi: 10.1097/00003246-199910000-00034.
[10] D. López-Padilla, M. Corral Blanco, S. Ferrer Espinosa, et al. "Unidades de Cuidados Respiratorios Intermedios: preguntas y respuestas [Intermediate Respiratory Care Units: Questions and Answers]". Spanish. In: Open Respiratory Archives 4.4 (2022), p. 100220. doi: 10.1016/j. opresp.2022.100220.
[11] S. Nava, M. Confalonieri, and C. Rampulla. "Intermediate respiratory intensive care units in Europe: a European perspective". In: Thorax 53.9 (1998), pp. 798–802. doi: 10.1136/thx.53.9. 798.
[12] J. D. J. Plate, L. P. H. Leenen, M. Houwert, et al. "Utilisation of Intermediate Care Units: A Systematic Review". In: Critical Care Research and Practice 2017 (2017), p. 8038460. doi: 10.1155/2017/8038460.
[13] O. Mediano, M. Luján, and S. López-Monzoni. "An Intermediate Respiratory Care Unit in Every Hospital". In: Archivos de Bronconeumología (English Edition) 59.1 (2023), pp. 7–9. doi: 10.1016/j.arbres.2022.09.021.
[14] P. López-Jardón, M. C. Martínez-Fernández, R. García-Fernández, et al. "Utility of Intermediate Care Units: A Systematic Review Study". In: Healthcare 12.3 (2024), p. 296. doi: 10.3390/ healthcare12030296.
[15] B. P. Krieger, P. Ershowsky, and D. Spivack. "One year's experience with a noninvasively monitored intermediate care unit for pulmonary patients". In: JAMA 264.10 (1990), pp. 1277– 1280. doi: 10.1001/jama.1990.03450090079031.
[16] S. Nava and N. Hill. "Non-invasive ventilation in acute respiratory failure". In: Lancet 374.9685 (2009), pp. 250–259. doi: 10.1016/S0140-6736(09)60496-7.
[17] M. Antonelli, G. Conti, M. Rocco, et al. "A comparison of noninvasive positive-pressure ventilation and conventional mechanical ventilation in patients with acute respiratory failure". In: The New England Journal of Medicine 339.7 (1998), pp. 429–435. doi: 10.1056/NEJM199808133390703.
[18] A. M. Esquinas, ed. Noninvasive Mechanical Ventilation. Theory, Equipment, and Clinical Applications. 2nd ed. Springer Cham, 2016, p. 959. doi: 10.1007/978-3-319-21653-9.
[19] M. Nishimura. "High-Flow Nasal Cannula Oxygen Therapy in Adults: Physiological Benefits, Indication, Clinical Benefits, and Adverse Effects". In: Respiratory Care 61.4 (2016), pp. 529–541. doi: 10.4187/respcare.04577.
[20] M. Luján, J. J. Monclou, P. Pérez, et al. "Multidisciplinary Consensus on the Management of Non-Invasive Respiratory Support in the COVID-19 Patient". In: Archivos de Bronconeumología 60.5 (2024), pp. 285–295. doi: 10.1016/j.arbres.2024.02.017.
[21] C. Caballero-Eraso, N. Pascual Martínez, O. Mediano, et al. "Unidades de Cuidados Respiratorios Intermedios (UCRI) durante la pandemia COVID-19. La realidad [Intermediate Respiratory Care Units (RICUs) During the COVID-19 Pandemic. The Reality]". In: Archivos de Bronconeumología 58.4 (2022), pp. 284–286. doi: 10.1016/j.arbres.2021.10.004.
[22] O. Grosgurin, A. Soulier, P. Soulier, et al. "Role of Intermediate Care Unit Admission and Noninvasive Respiratory Support during the COVID-19 Pandemic: A Retrospective Cohort Study". In: Respiration 100.9 (2021), pp. 786–793. doi: 10.1159/000516329.
[23] M. Matute-Villacís, C. Díaz-Ángel, G. Burbano, et al. "Role of respiratory intermediate care units during the SARS-CoV-2 pandemic". In: BMC Pulmonary Medicine 21.1 (2021), p. 228. doi: 10.1186/s12890-021-01593-5.
[24] G. Suarez-Cuartin, J. Osorio, M. Salgueira, et al. "Clinical Outcomes of Severe COVID-19 Patients Admitted to an Intermediate Respiratory Care Unit". In: Frontiers in Medicine 8 (2021), p. 711027. doi: 10.3389/fmed.2021.711027.
[25] R. Scala and L. Pisani. "Non-invasive ventilation in acute respiratory failure: Which recipe?" In: European Respiratory Review 27.150 (2018), p. 180029. doi: 10.1183/16000617.0029-2018.
[26] P. K. Lindenauer, M. S. Stefan, M.-S. Shieh, et al. "Outcomes Associated With Invasive and Noninvasive Ventilation Among Patients Hospitalized With Exacerbations of Chronic Obstructive Pulmonary Disease". In: JAMA Internal Medicine 174.12 (2014), pp. 1982–1993. doi: 10.1001/ jamainternmed.2014.5430.
[27] G. Hilbert, D. Gruson, F. Vargas, et al. "Noninvasive ventilation in immunosuppressed patients with pulmonary infiltrates, fever, and acute respiratory failure". In: New England Journal of Medicine 344.7 (2001), pp. 481–487. doi: 10.1056/NEJM200102153440703.
[28] M. Ferrer, A. Esquinas, M. Leon, et al. "Noninvasive ventilation in severe hypoxemic respiratory failure: a randomized clinical trial". In: American Journal of Respiratory and Critical Care Medicine 168.12 (2003), pp. 1438–1444. doi: 10.1164/rccm.200301-072OC.
[29] S. Heili-Frades, R. Ferrer-Luca, A. Zapata-Fenor, et al. "Cost and Mortality Analysis of an Intermediate Respiratory Care Unit. Is It Really Efficient and Safe?" In: Archivos de Bronconeumología (English Edition) 55.12 (2019), pp. 634–641. doi: 10.1016/j.arbres.2019.06.008.
[30] M. Capuzzo, C. A. Volta, T. Tassinati, et al. "Hospital mortality of adults admitted to Intensive Care Units in hospitals with and without Intermediate Care Units: a multicentre European cohort study". In: Critical Care 18.5 (2014), p. 551. doi: 10.1186/s13054-014-0551-8.
[31] M. Prin and H. Wunsch. "The role of stepdown beds in hospital care". In: American Journal of Respiratory and Critical Care Medicine 190.11 (2014), pp. 1210–1216. doi: 10.1164/rccm. 201406-1117PP.
[32] E. Sala, C. Balaguer, M. Carrera, et al. "Actividad de una unidad de cuidados respiratorios intermedios dependiente de un servicio de neumología [Activity of an intermediate respiratory care unit attached to a respiratory medicine department]". Spanish. In: Archivos de Bronconeumología 45.4 (2009). Epub 2009 Mar 14, pp. 168–172. doi: 10.1016/j.arbres.2008.09.003.
[33] M. Galdeano Lozano, D. Catalan-Matamoros, and M. Trapero-Bertran. "Effectiveness of Intermediate Respiratory Care Units as an Alternative to Intensive Care Units during the COVID-19 Pandemic in Catalonia". In: International Journal of Environmental Research and Public Health 19.10 (2022), p. 6034. doi: 10.3390/ijerph19106034.
[34] R. Scala and L. Pisani. "Noninvasive ventilation in acute respiratory failure: which recipe for success?" In: European Respiratory Review 27.150 (2018), p. 180029. doi: 10.1183/16000617. 0029-2018.
[35] N. Ambrosino and G. Vagheggini. "Noninvasive positive pressure ventilation in the respiratory intermediate care unit: back to the future?" In: American Journal of Respiratory and Critical Care Medicine 184.1 (2011), pp. 10–11. doi: 10.1164/rccm.201103-0551ED.
[36] R. Cabestre García. "Respiratory Nursing: Soul, Brain, and Heart of Intermediate Respiratory Care Units". In: Archivos de Bronconeumología (English Edition) 59.12 (2023), pp. 789–790. doi: 10.1016/j.arbres.2023.07.022.
[37] L. Rose, S. Nelson, L. Johnston, et al. "Decisions made by critical care nurses during mechanical ventilation and weaning in an Australian intensive care unit". In: American Journal of Critical Care 16.5 (2007). Quiz 444, pp. 434–443. doi: 10.4037/ajcc2007.16.5.434.
[38] A. S. Multz, D. B. Chalfin, I. M. Samson, et al. "A "closed" medical intensive care unit (MICU) improves resource utilization when compared with an "open" MICU". In: American Journal of Respiratory and Critical Care Medicine 157.5 Pt 1 (1998), pp. 1468–1473. doi: 10.1164/ajrccm. 157.5.9708039.
[39] L. H. Aiken, D. M. Sloane, L. Bruyneel, et al. "Nurse staffing and education and hospital mortality in nine European countries: a retrospective observational study". In: The Lancet 383.9931 (2014), pp. 1824–1830. doi: 10.1016/S0140-6736(13)62631-8.
[40] J. Needleman, P. Buerhaus, V. S. Pankratz, et al. "Nurse staffing and inpatient hospital mortality". In: New England Journal of Medicine 364.11 (2011), pp. 1037–1045. doi: 10.1056/ NEJMsa1001025.
[41] P. Griffiths, J. Ball, K. Bloor, et al. "Nurse staffing levels, missed vital signs and mortality in hospitals: retrospective longitudinal observational study". In: (2018). doi: 10.3310/hsdr06380.
[42] R. L. Kane, T. A. Shamliyan, C. Mueller, et al. "The Association of Registered Nurse Staffing Levels and Patient Outcomes: Systematic Review and Meta-Analysis". In: Medical Care 45.12 (2007), pp. 1195–1204. doi: 10.1097/MLR.0b013e3181468ca3.
[43] P. J. Pronovost, D. C. Angus, T. Dorman, et al. "Physician Staffing Patterns and Clinical Outcomes in Critically Ill Patients: A Systematic Review". In: JAMA 288.17 (2002), pp. 2151– 2162. doi: 10.1001/jama.288.17.2151.
[44] N. S. Hill, J. Brennan, E. Garpestad, et al. "Noninvasive ventilation in acute respiratory failure". In: Critical Care Medicine 35.10 (2007), pp. 2402–2407. doi: 10.1097/01.CCM.0000284587. 36541.7F.
[45] M. Davies, T. Quinnell, N. Oscroft, et al. "British Thoracic Society Quality Standards for acute non-invasive ventilation in adults". In: BMJ Open Respiratory Research 5.1 (2018), e000283. doi: 10.1136/bmjresp-2018-000283.
[46] J. L. Vincent, S. Einav, R. Pearse, et al. "Improving detection of patient deterioration in the general hospital ward environment". In: European Journal of Anaesthesiology 35.5 (2018), pp. 325– 333. doi: 10.1097/EJA.0000000000000798.
[47] A. Carlucci, L. Pisani, P. Ceriana, et al. "Patient-ventilator asynchronies: may the respiratory mechanics play a role?" In: Critical Care 17.2 (2013), R54. doi: 10.1186/cc12580.
[48] A. Corrado, C. Roussos, N. Ambrosino, et al. "Respiratory intermediate care units: A European survey". In: European Respiratory Journal 20.5 (2002), pp. 1343–1350. doi: 10.1183/09031936. 02.00058202.
[49] E. Garpestad, J. Brennan, and N. S. Hill. "Noninvasive ventilation for critical care". In: Chest 132.2 (2007), pp. 711–720. doi: 10.1378/chest.06-2643.
[50] R. Hall, ed. Patient Flow: Reducing Delay in Healthcare Delivery. 2nd ed. International Series in Operations Research & Management Science. Springer New York, NY, 2013, pp. XII, 553. doi: 10.1007/978-1-4614-9512-3.
[51] E. Litvak and M. Bisognano. "More patients, less payment: increasing hospital efficiency in the aftermath of health reform". In: Health Affairs (Millwood) 30.1 (2011), pp. 76–80. doi: 10.1377/hlthaff.2010.1114.
[52] N. C. Proudlove, K. Gordon, and R. Boaden. "Can simulating queues help practitioners plan better health services?" In: Journal of Health Organization and Management 22.4 (2008), pp. 341– 355. doi: 10.1108/14777260810893900.
[53] J. A. Sánchez-Martínez, C. Morales-García, B. C. Gómez-Hernández, et al. "Development of a Non-Invasive Respiratory Support Pathway in COVID-19 patients with severe respiratory failure in an Intermediate Respiratory Care Unit". In: European Respiratory Journal 60.suppl 66 (2022), p. 2040. doi: 10.1183/13993003.congress-2022.2040.
[54] R. M. Kleinpell. "Acute care nurse practitioner practice: Results of a 5-year longitudinal study". In: American Journal of Critical Care 14.3 (2005), pp. 211–219. doi: 10.4037/ajcc2005.14.3.211.
[55] S. Nava, N. Ambrosino, E. Clini, et al. "Non-invasive mechanical ventilation in the weaning of patients with respiratory failure due to chronic obstructive pulmonary disease. A randomized, prospective trial". In: Annals of Internal Medicine 128.9 (1998), pp. 721–728. doi: 10.7326/0003- 4819-128-9-199805010-00004.
[56] O. Diekmann, H. Heesterbeek, and T. Britton. Mathematical Tools for Understanding Infectious Disease Dynamics. Princeton University Press, 2013.
[57] M. L. Brandeau, F. Sainfort, and W. P. Pierskalla, eds. Operations Research and Health Care. A Handbook of Methods and Applications. 1st ed. International Series in Operations Research & Management Science. New York, NY: Springer New York, 2004, pp. VIII, 874. doi: 10.1007/ b106574.
[58] A. Alban, S. Chick, D. Dongelmans, et al. "ICU capacity management during the COVID-19 pandemic using a process simulation". In: Intensive Care Medicine 46.8 (2020). Epub 2020 May 7, pp. 1624–1626. doi: 10.1007/s00134-020-06066-7.
[59] T. Wang, Y. Wu, J. Y.-N. Lau, et al. "A four-compartment model for the COVID-19 infection—implications on infection kinetics, control measures, and lockdown exit strategies". In: Precision Clinical Medicine 3.2 (2020), pp. 104–112. doi: 10.1093/pcmedi/pbaa018.
[60] W. O. Kermack and A. G. McKendrick. "A Contribution to the Mathematical Theory of Epidemics". In: Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 115.772 (1927), pp. 700–721. doi: 10.1098/rspa.1927.0118.
[61] S. Robinson. Simulation: The Practice of Model Development and Use. Bloomsbury Academic, 2014.
[62] S. H. Jacobson, S. N. Hall, and J. R. Swisher. "Discrete-Event Simulation of Health Care Systems". In: Patient Flow: Reducing Delay in Healthcare Delivery. Ed. by R. W. Hall. Boston, MA: Springer US, 2006, pp. 211–252. doi: 10.1007/978-0-387-33636-7\_8.
[63] S. Zhang, J. Ponce, Z. Zhang, et al. "An integrated framework for building trustworthy datadriven epidemiological models: Application to the COVID-19 outbreak in New York City". In: PLoS Computational Biology 17.9 (2021), e1009334. doi: 10.1371/journal.pcbi.1009334.
[64] A. Saltelli, M. Ratto, T. Andres, et al. Global Sensitivity Analysis: The Primer. John Wiley & Sons, 2008.
[65] H. T. Banks and K. Kunisch. Estimation Techniques for Distributed Parameter Systems. Birkhäuser, 1989. doi: 10.1007/978-1-4612-3700-6.
[66] B. Oner, F. D. Zengul, N. Oner, et al. "Nursing-sensitive indicators for nursing care: A systematic review (1997–2017)". In: Nursing Open 8.3 (2021), pp. 1005–1022. doi: https://doi.org/10. 1002/nop2.654.
[67] K. Gopal, L. S. Lee, and H.-V. Seow. "Parameter Estimation of Compartmental Epidemiological Model Using Harmony Search Algorithm and Its Variants". In: Applied Sciences 11.3 (2021), p. 1138. doi: 10.3390/app11031138.
[68] D. López-Padilla, F. González Martínez, J. Torres-Macho, et al. "Unidades de Cuidados Respiratorios Intermedios: preguntas y respuestas [Intermediate Respiratory Care Units: Questions and Answers]". In: Open Respiratory Archives 4.4 (2022), p. 100220. doi: 10.1016/j.opresp.2022. 100220.
[69] R. Mellado-Artigas, B. L. Ferreyro, F. Angriman, et al. "High-flow nasal oxygen in patients with COVID-19-associated acute respiratory failure". In: Critical Care 25.1 (2021), p. 58. doi: 10.1186/s13054-021-03469-w.
[70] J. Needleman. "Nurse staffing: the knowns and unknowns". In: Nursing Economic\$ 33.1 (2015), pp. 5–7.
[71] J. Needleman, P. Buerhaus, S. Mattke, et al. "Nurse-staffing levels and the quality of care in hospitals". In: New England Journal of Medicine 346.22 (2002), pp. 1715–1722. doi: 10.1056/ NEJMsa012247.
[72] L. H. Aiken, S. P. Clarke, R. B. Cheung, et al. "Educational levels of hospital nurses and surgical patient mortality". In: JAMA 290.12 (2003), pp. 1617–1623. doi: 10.1001/jama.290.12.1617.
[73] C. Caballero-Eraso, N. Pascual Martínez, O. Mediano, et al. "Unidades de Cuidados Respiratorios Intermedios (UCRI) durante la pandemia COVID-19. La realidad [Intermediate Respiratory Care Units (RICUs) During the COVID-19 Pandemic. The Reality]". In: Archivos de Bronconeumología 58.4 (2022), pp. 284–286. doi: 10.1016/j.arbres.2021.10.004.
[74] T. Al-Karkhi and K. Byatt. "A compartmental model to describe acute medical in-patient flow through a hospital". In: Heliyon 11.3 (2025), e42260. doi: 10.1016/j.heliyon.2025.e42260.
[75] C. Kahraman and Y. I. Topcu, eds. Operations Research Applications in Health Care Management. 1st ed. International Series in Operations Research & Management Science 255. Cham: Springer Cham, 2018, pp. XV, 604. doi: 10.1007/978-3-319-65455-3.
[76] J.-P. Frat, A. W. Thille, A. Mercat, et al. "High-Flow Oxygen through Nasal Cannula in Acute Hypoxemic Respiratory Failure". In: New England Journal of Medicine 372.23 (2015), pp. 2185– 2196. doi: 10.1056/NEJMoa1503326.
[77] C. E. Rasmussen and C. K. I. Williams. Gaussian Processes for Machine Learning. The MIT Press, 2005. doi: 10.7551/mitpress/3206.001.0001.
[78] B. M. Jiménez-Rodríguez, I. Bellido-Estévez, I. Hernández-García, et al. "On the single and multiple associations of COVID-19 post-acute sequelae: 6-month prospective cohort study". In: Scientific Reports 12.1 (2022), p. 3402. doi: 10.1038/s41598-022-07433-8.
[79] J. S. Whittle, I. Pavlov, A. D. Sacchetti, et al. "Respiratory support for adult patients with COVID-19". In: Journal of the American College of Emergency Physicians Open 1.2 (2020), pp. 95–101. doi: 10.1002/emp2.12071.
[80] World Health Organization. Clinical management of COVID-19: interim guidance, 27 May 2020. Tech. rep. WHO/2019-nCoV/clinical/2020.5. World Health Organization, 2020.
[81] Inkscape. https://inkscape.org/. Vector Graphics Editor, Version 1.3 used. Accessed March 30, 2025. 2024.
[82] F. Faul, E. Erdfelder, A.-G. Lang, et al. "G*Power 3: A flexible statistical power analysis program for the social, behavioral, and biomedical sciences". In: Behavior Research Methods 39.2 (2007), pp. 175–191. doi: 10.3758/bf03193146.
[83] J. Cohen. Statistical Power Analysis for the Behavioral Sciences. Taylor & Francis, 2013.
Ana Carmen Navas-Ortega!1 , José Antonio Sánchez-Martínez!1,2, Paula García-Flores1,2 , Concepción Morales-García1,2, and Rene Fabregas3
1Department of Pneumology, University Hospital Virgen de Las Nieves, Granada, Spain. 2Biosanitary Research Institute of Granada-Ibs, Granada, Spain.
Corresponding Author: Concepción Morales-García, concepcion.morales.sspa@juntadeandalucia.es and Rene Fabregas, rfabregas@ugr.es.
Thursday 12th June, 2025
Effective management of Intermediate Respiratory Care Units (IRCUs) requires a quantitative understanding of patient flow dynamics, particularly under varying admission loads and evolving clinical interventions. This section introduces a discrete-time compartmental model specifically designed to capture patient progression within an IRCU, explicitly accounting for the use of non-invasive ventilation (NIV) and the impact of targeted clinical management strategies. By subdividing the IRCU population based on NIV status (state X for non-NIV, state Y for active NIV), the model utilizes a system of difference equations to characterize patient trajectories including admission, NIV initiation, recovery pathways, transfer to the Intensive Care Unit (ICU), and mortality (exitus), focusing on flows originating within the IRCU environment. Crucially, the model integrates parameters representing key clinical interventions aimed at modulating patient outcomes, providing a framework for evaluating their effectiveness.
We define the state variables at discrete time t (typically days) as follows:
3Department of Applied Mathematics and Modeling Nature (MNat) Research Unit, Faculty of Sciences, University of Granada, Granada, Spain.
‡ These authors contributed equally.
Let A(t) denote the influx of new patient admissions directed to the IRCU system during the time interval [t, t + 1) [patients/time step]. The model dynamics are governed by baseline transition parameters (fractional rates or probabilities per time step) and specific intervention modulation parameters:
*(Note: Dynamics post-ICU transfer, such as ICU length of stay or outcomes from the ICU, are considered outside the scope of this IRCU-focused model).*
The structure of patient flow, incorporating the points of intervention, is illustrated in Figure S1. Admissions A(t) feed into the non-NIV state X(t). Patients in X(t) can recover directly (R), be transferred to ICU (Z), experience exitus (W), or initiate NIV therapy, moving to state Y (t). From the NIV state Y (t), patients can recover (R), be transferred to ICU (Z), or experience exitus (W). The intervention modulators (*µ0, ε0,* ∆¹, ¼) specifically act upon the transition rates originating from the active NIV state Y (t), reflecting targeted efforts to improve outcomes for this high-acuity patient group. Compartmental diagrams of this nature are standard tools for visualizing system dynamics in mathematical biology and healthcare modeling [7, 8, 9].
Figure S1| Integrated Compartmental Model of IRCU Patient Flow with Intervention Pathways. Nodes represent patient states: Admissions (A(t)), patients within the IRCU distinguished by NIV status (X(t): no NIV; Y (t): active NIV, highlighted within the IRCU region), and cumulative outcomes (R(t): recovered, Z(t): ICU transfer, W(t): exitus). Solid arrows depict patient flow rates. Parameters α, ρ, η, ν govern transitions from state X. Parameters γ, ε, θ0 represent baseline transition rates from the active NIV state Y . The 'Intervention Modulators' (*γ0, ε0,* ∆θ, λ, depicted in the pink box) quantify the impact of clinical strategies by modifying the effective rates of ICU transfer (γ(1 − *γ*0)), exitus (ε(1 − *ε*0)), and recovery (θ(t) = θ0 + ∆θe−λt) from state Y , with conceptual links indicated by faded arrows.
The mathematical formulation translates the flows depicted in Figure S1 into a system of firstorder difference equations. Each term corresponds to an inflow or outflow from a compartment, incorporating the relevant baseline rates and intervention modulators:
$$X(t+1) = X(t) + A(t) - (\alpha + \rho + \eta + \nu)X(t)$$ (S1)
$$Y(t+1) = Y(t) + \alpha X(t) - \left[ \gamma(1-\gamma_0) + \underbrace{\left(\theta0 + \Delta \theta e^{-\lambda t}\right)}{\theta(t)} + \varepsilon(1-\varepsilon_0) \right] Y(t)$$ (S2)
$$Z(t+1) = Z(t) + \gamma(1 - \gamma_0)Y(t) + \eta X(t)$$ (S3)
$$W(t+1) = W(t) + \varepsilon(1 - \varepsilon_0)Y(t) + \nu X(t)$$ (S4)
$$R(t+1) = R(t) + \left[\theta_0 + \Delta\theta e^{-\lambda t}\right] Y(t) + \rho X(t) \tag{S5}$$
Equation (S1) describes the change in the non-NIV population (X), balancing admissions A(t) against outflows due to NIV initiation (³), direct recovery (Ä), direct ICU transfer (¸), and direct exitus (¿). Equation (S2) tracks the active NIV population (Y ), with inflow from state X and outflows to ICU, exitus, and recovery, where these rates are potentially modified by interventions (*µ0, ε0,* ∆¹, ¼). Equations (S3), (S4), and (S5) accumulate patients transitioning into the ICU, exitus, and recovered states, respectively, incorporating contributions from both X and the modulated Y pathways.
This model structure aligns with observed clinical workflows in intermediate care settings [10, 11] and adheres to established principles of compartmental modeling in mathematical biology and epidemiology [12, 13]. Its explicit incorporation of intervention parameters enhances its utility for healthcare operational analysis. The framework allows for quantitative exploration of various clinical scenarios:
Such compartmental models are foundational for understanding population dynamics [7, 8, 14] and provide critical tools for healthcare planning, resource allocation, and system optimization [15, 16]. Furthermore, the deterministic core can be integrated with data-driven approaches (e.g., Bayesian inference, machine learning) for parameter estimation, uncertainty quantification, and real-time adaptation [9, 17], enhancing predictive accuracy. Consequently, this model serves as a robust computational platform—a form of "digital twin"—for the *in silico* evaluation of operational policies and clinical protocols, facilitating evidence-informed decision-making to improve patient outcomes and system efficiency within the IRCU [6].
While the discrete-time model provides valuable insights into patient flow dynamics at specific intervals, transitioning to a continuous-time framework offers complementary advantages. Continuous models, expressed as systems of ordinary differential equations (ODEs), allow for the application of calculus-based analytical techniques and can represent underlying processes that evolve smoothly over time. This section details the derivation of such a continuous model from our discrete formulation.
Recall the discrete model incorporates patient movements between states during each time interval ∆t. Considering small intervals, we approximate the change in patient numbers using instantaneous rates. Let A(t) be the rate of new admissions [patients/time] and ³, Ä, µ, etc., be instantaneous transition rates [1/time]. The net change in a compartment f over ∆t is approximately the sum of inflows minus outflows:
$$\begin{split} X(t+\Delta t) - X(t) &\approx \left[ A(t) - (\alpha + \rho + \eta + \nu) X(t) \right] \Delta t, \ Y(t+\Delta t) - Y(t) &\approx \left[ \alpha X(t) - \left( \gamma (1-\gamma_0) + \theta(t) + \varepsilon (1-\varepsilon_0) \right) Y(t) \right] \Delta t, \ Z(t+\Delta t) - Z(t) &\approx \left[ \gamma (1-\gamma_0) Y(t) + \eta X(t) \right] \Delta t, \ W(t+\Delta t) - W(t) &\approx \left[ \varepsilon (1-\varepsilon_0) Y(t) + \nu X(t) \right] \Delta t, \ R(t+\Delta t) - R(t) &\approx \left[ \theta(t) Y(t) + \rho X(t) \right] \Delta t. \end{split}$$
Using the Taylor expansion f(t + ∆t) ≈ f(t) + *df dt*∆t and taking the limit as ∆t → 0, we derive the continuous-time ODE system:
$$\frac{dX}{dt} = A(t) - (\alpha + \rho + \eta + \nu)X(t), \tag{S6}$$
$$\frac{dY}{dt} = \alpha X(t) - \left[\gamma(1 - \gamma_0) + \theta(t) + \varepsilon(1 - \varepsilon_0)\right] Y(t), \tag{S7}$$
$$\frac{dZ}{dt} = \gamma(1 - \gamma_0)Y(t) + \eta X(t),\tag{S8}$$
$$\frac{dW}{dt} = \varepsilon (1 - \varepsilon_0) Y(t) + \nu X(t), \tag{S9}$$
$$\frac{dR}{dt} = \theta(t)Y(t) + \rho X(t). \tag{S10}$$
This continuous system explicitly incorporates the clinical interventions: Mortality Reduction (*ε*0): Scaled rate ε(1−*ε*0) reflects mitigation from enhanced protocols [1]. ICU Transfer Mitigation (*µ*0): Reduced rate µ(1−*µ*0) models optimized NIV management preventing ICU transfers [4]. Transient Staffing Effects (∆*¹, ¼*): Time-dependent recovery ¹(t) = ¹0 + ∆*¹e*−¼t represents enhanced potential during staffing surges [5, 6]. The inclusion of ¹(t) renders the system non-autonomous.
To understand core dynamics, we employ dimensional analysis. We define characteristic scales: Ä = (³ + Ä + ¸ + ¿) −1 (characteristic time in state X) and N0 = *A*0Ä (characteristic population size, with A0 a baseline admission rate).
Dimensionless variables (*t,*˜ X, ˜ Y , ˜ Z, ˜ W, ˜ R, ˜ *A*˜) and parameters (Ã) are introduced:
$$\tilde{t} = t/\tau, \quad \tilde{X} = X/N_0, \quad \dots \quad \tilde{A}(\tilde{t}) = A(t)/A_0$$
$$\pi{\alpha} = \alpha \tau, \quad \pi{\rho} = \rho \tau, \quad \pi{\eta} = \eta \tau, \quad \pi{\nu} = \nu \tau, \quad \text{(Note: } \pi{\alpha} + \pi{\rho} + \pi{\eta} + \pi{\nu} = 1\text{)}$$
$$\pi{\gamma} = \gamma \tau, \quad \pi{\varepsilon} = \varepsilon \tau, \quad \pi{\theta}(\tilde{t}) = \theta(t)\tau, \quad \pi{\lambda} = \lambda \tau$$
Substituting these yields the dimensionless system:
$$\frac{d\tilde{X}}{d\tilde{t}} = \tilde{A}(\tilde{t}) - \tilde{X},\tag{S11}$$
$$\frac{d\tilde{Y}}{d\tilde{t}} = \pi{\alpha}\tilde{X} - \left[\pi{\gamma}(1 - \gamma0) + \pi{\theta}(\tilde{t}) + \pi_{\varepsilon}(1 - \varepsilon_0)\right]\tilde{Y},\tag{S12}$$
$$\frac{d\tilde{Z}}{d\tilde{t}} = \pi_{\gamma}(1 - \gamma0)\tilde{Y} + \pi{\eta}\tilde{X},\tag{S13}$$
$$\frac{d\tilde{W}}{d\tilde{t}} = \pi_{\varepsilon}(1 - \varepsilon0)\tilde{Y} + \pi{\nu}\tilde{X},\tag{S14}$$
$$\frac{d\tilde{R}}{d\tilde{t}} = \pi{\theta}(\tilde{t})\tilde{Y} + \pi{\rho}\tilde{X}. \tag{S15}$$
where ù(*t*˜) = ù0 + *Ã*∆¹e −*Ãλt*˜ .
While the general non-autonomous system ((S11)–(S15)) typically requires numerical integration due to the time-dependent coefficient ù(*t*˜), an exact analytical solution can be derived for the important special case where the system is autonomous. This occurs when the transient staffing effect is absent (*Ã∆¹* = 0, thus ù(*t*˜) = ù0 = constant) and the admission rate is constant (*A*˜(*t*˜) = *A*˜). Assuming zero initial conditions (*X*˜(0) = *Y*˜ (0) = *Z*˜(0) = *W*˜ (0) = *R*˜(0) = 0), the solution unfolds sequentially.
First, the equation for *X*˜ (S11) yields:
$$\tilde{X}(t) = \tilde{A}(1 - e^{-t}).$$
Next, defining the constant total outflow rate from state Y as K1 = õ(1 − *µ*0) + ù0 + Ãε(1 − *ε*0), the equation for *Y*˜ (S12) becomes a standard linear first-order ODE with constant coefficients forced by *X*˜(t). Assuming K1 ̸= 1 and K1 ̸= 0, its solution is:
$$\tilde{Y}(t) = \frac{\pi_{\alpha}\tilde{A}}{K_1(K_1 - 1)} \left[ (K_1 - 1) - K_1 e^{-t} + e^{-K_1 t} \right].$$
(In the specific case where K1 = 1*, the solution takes the form Y*˜ (t) = *óA*˜(1 − e −t − *te−t* )).
Finally, the cumulative variables *Z*˜(t), *W*˜ (t), and *R*˜(t) are obtained by direct integration of their respective differential equations (S13), (S14), and (S15), substituting the derived expressions for *X*˜(t) and *Y*˜ (t). The resulting solutions are linear combinations of a term linear in time t, a constant term, and exponential terms e −t and e −*K*1t . For instance, the solution for cumulative ICU transfers *Z*˜(t) takes the form:
$$\tilde{Z}(t) = A_Z t + B_Z'(1 - e^{-t}) + C_Z'(1 - e^{-K_1 t}),$$
where AZ = *A*˜(õ ′*ó/K*1+ø) is the asymptotic accumulation rate (with õ ′ = õ(1−*µ*0)), and *B*′ *Z , C*′ Z are constants derived from the integration involving BZ and CZ as defined in the detailed derivation. Analogous expressions hold for *W*˜ (t) and *R*˜(t). This analytical solution for the autonomous case provides a valuable baseline for understanding the system's fundamental dynamics and validating numerical solvers used for the more general non-autonomous system.
We analyze the system's behavior under two conditions: the autonomous case, representing baseline operations with constant parameters, and the non-autonomous case, incorporating time-dependent interventions like staffing surges.
First, consider the system under constant baseline conditions, where the admission rate is constant (*A*˜(*t*˜) = 1) and the recovery rate from NIV is also constant (ù(*t*˜) = ù0 , i.e., *Ã∆¹* = 0). In this autonomous scenario, we can determine the equilibrium points (steady states) by setting the time derivatives in Equations (S11) and (S12) to zero. This yields the steady-state populations for the IRCU compartments:
$$\tilde{X}^* = 1, \quad \tilde{Y}^* = \frac{\pi{\alpha}}{\pi{\gamma}(1 - \gamma0) + \pi{\theta0} + \pi{\varepsilon}(1 - \varepsilon_0)}.$$
The steady-state occupancy levels depend directly on the relative transition rates (Ã values) and the baseline intervention levels (*µ0, ε*0). A standard stability analysis of the Jacobian matrix for the (X, ˜ *Y*˜ ) subsystem confirms that this equilibrium (*X*˜ ∗ *, Y*˜ ∗ ) is locally asymptotically stable for all physically realistic (non-negative) parameter values.
For the cumulative compartments (Z, ˜ W, ˜ *R*˜), their governing equations involve only non-negative inflow terms when evaluated at the steady state (*X*˜ ∗ *, Y*˜ ∗ ). For instance, the rate of change for cumulative ICU transfers becomes:
$$\left(\frac{d\tilde{Z}}{d\tilde{t}}\right)^* = \pi_{\gamma}(1-\gamma0)\tilde{Y}^* + \pi{\eta}\tilde{X}^* = \frac{\pi_{\gamma}(1-\gamma0)\pi{\alpha}}{\pi_{\gamma}(1-\gamma0) + \pi{\theta0} + \pi{\varepsilon}(1-\varepsilon0)} + \pi{\eta}.$$
Since this rate is constant and positive, the cumulative states *Z*˜(t), *W*˜ (t), and *R*˜(t) do not reach a finite steady state but instead grow linearly over time at constant rates determined by the steady-state populations *X*˜ ∗ and *Y*˜ ∗ . The primary operational concern is therefore not instability, but whether these constant accumulation rates (particularly for *Z*˜) lead to exceeding downstream capacities (like ICU beds) over extended periods.
When interventions introduce time-dependency, such as the transient staffing surge modeled by ù(t) with *Ã∆¹ >* 0 and ü > 0, the system becomes non-autonomous. By definition, such systems do not possess fixed steady-state equilibrium points while the parameters are changing. The system's state (*X*˜(t)*, Y*˜ (t), etc.) evolves dynamically according to the time-varying recovery rate.
However, we can analyze the long-term (asymptotic) behavior. As *t*˜ → ∞, the transient term *Ã*∆¹e −*Ãλt*˜ decays to zero. Consequently, ù(t) → ù0 . This means that the non-autonomous system asymptotically approaches the behavior of the autonomous system discussed above. The state variables *X*˜(t) and *Y*˜ (t) will converge towards the steady-state values *X*˜ ∗ and *Y*˜ ∗ calculated using the baseline recovery rate ù0 . Similarly, the *rates* of accumulation for *Z*˜(t), *W*˜ (t), and *R*˜(t) will asymptotically approach the constant rates derived for the autonomous case.
Therefore, while the system exhibits transient dynamics during the period where the staffing surge effect is significant, its long-term equilibrium behavior is governed by the baseline parameters. The primary impact of the transient intervention (*Ã*∆¹, ü) lies in altering the *trajectory* towards this asymptotic state, potentially mitigating peak loads or accelerating recovery during critical periods, as explored in the numerical simulations (Section S3.2).
To understand which parameters most significantly influence the system's long-term behavior, we perform a parametric sensitivity analysis. This analysis focuses on the steady-state Non-Invasive Ventilation (NIV) occupancy (*Y*˜ ∗ ) derived from the **autonomous system** (or, equivalently, the asymptotic state approached by the non-autonomous system as t → ∞). We quantify the influence of each parameter p using dimensionless relative sensitivity indices, defined as Sp = *∂Y*˜ ∗ ∂p · *p Y*˜ ∗ . This index measures the proportional change in *Y*˜ ∗ resulting from a proportional change in the parameter p. A positive index indicates that increasing the parameter increases *Y*˜ ∗ , while a negative index signifies the opposite relationship. This analysis helps identify the most effective levers for managing NIV resource utilization within the IRCU.
Sensitivity to NIV Admission Rate (*ó*). The sensitivity of the equilibrium NIV occupancy to the relative rate at which non-NIV patients transition to requiring NIV is calculated as:
$$S{\pi{\alpha}} = \frac{\partial \tilde{Y}^}{\partial \pi{\alpha}} \cdot \frac{\pi{\alpha}}{\tilde{Y}^} = 1$$
This index value of exactly 1 signifies a direct, linear, and proportional relationship. A 10% increase in ó (representing, for example, an influx of patients with higher acuity requiring NIV sooner) leads directly to a 10% increase in the long-term equilibrium occupancy of NIV beds (*Y*˜ ∗ ). This highlights the system's direct vulnerability to changes in the clinical characteristics of incoming patients or factors that accelerate the need for NIV initiation. Managing this sensitivity relies primarily on external factors influencing admission acuity or internal processes affecting the ³ transition, rather than downstream interventions from state Y.
Sensitivity to Baseline NIV Recovery Rate (ù0 ). The influence of the baseline rate at which patients recover directly from NIV is given by:
$$S{\pi{\theta0}} = \frac{\partial \tilde{Y}^*}{\partial \pi{\theta0}} \cdot \frac{\pi{\theta0}}{\tilde{Y}^*} = \frac{-\pi{\theta0}}{\pi{\gamma}(1-\gamma0) + \pi{\theta0} + \pi{\varepsilon}(1-\varepsilon_0)}$$
This sensitivity index is inherently negative (SÃθ0 < 0), indicating that enhancing the baseline recovery rate from NIV (ù0 ) serves to decrease the long-term occupancy *Y*˜ ∗ . The magnitude of this effect is determined by the proportion that baseline recovery contributes to the total outflow rate from the NIV state (Y). Interventions aimed at improving general NIV care efficacy, shortening duration of NIV dependency, or accelerating recovery would leverage this sensitivity to reduce the equilibrium burden on NIV resources.
Sensitivity to Mitigated ICU Transfer Rate (õ(1 − *µ*0)). The impact of the rate at which NIV patients deteriorate and require transfer to the ICU, considering mitigation efforts (*µ*0), is:
$$S{\pi{\gamma}(1-\gamma0)} = \frac{\partial \tilde{Y}^*}{\partial (\pi{\gamma}(1-\gamma0))} \cdot \frac{\pi{\gamma}(1-\gamma0)}{\tilde{Y}^*} = \frac{-\pi{\gamma}(1-\gamma0)}{\pi{\gamma}(1-\gamma0) + \pi{\theta0} + \pi{\varepsilon}(1-\varepsilon_0)}$$
This index is also negative (SÃγ(1−*µ*0) < 0). Increasing the effectiveness of protocols designed to prevent avoidable ICU transfers (i.e., increasing *µ*0, which decreases the overall term õ(1 − *µ*0)) leads to a reduction in the equilibrium NIV occupancy *Y*˜ ∗ . The magnitude of this sensitivity depends on the relative importance of the ICU transfer pathway compared to recovery and exitus pathways from state Y. If ICU transfer is a dominant pathway (large õ), interventions targeting µ0 can be particularly impactful levers for managing *Y*˜ ∗ .
Sensitivity to Mitigated NIV Exitus Rate (Ãε(1 − *ε*0)). Finally, the sensitivity to the rate of mortality occurring directly from the NIV state, accounting for mortality reduction interventions (*ε*0), is:
$$S{\pi{\varepsilon}(1-\varepsilon0)} = \frac{\partial \tilde{Y}^*}{\partial (\pi{\varepsilon}(1-\varepsilon0))} \cdot \frac{\pi{\varepsilon}(1-\varepsilon0)}{\tilde{Y}^*} = \frac{-\pi{\varepsilon}(1-\varepsilon0)}{\pi{\gamma}(1-\gamma0) + \pi{\theta0} + \pi{\varepsilon}(1-\varepsilon_0)}$$
Similar to recovery and ICU transfer, this sensitivity is negative (SÃε(1−*ε*0) < 0). Implementing protocols or enhancing care to reduce mortality while on NIV (i.e., increasing *ε*0, thus decreasing the term Ãε(1 − *ε*0)) contributes to lowering the long-term NIV occupancy *Y*˜ ∗ , as fewer patients remain in state Y prior to death. The relative impact compared to modifying recovery or ICU transfer rates depends on the baseline mortality rate Ãε and the effectiveness of the reduction intervention *ε*0.
It must be emphasized that these sensitivity indices quantify the impact of parameter variations on the *equilibrium* state *Y*˜ ∗ of the autonomous system, or equivalently, the *asymptotic* state approached by the non-autonomous system long after transient effects have subsided. During periods where time-dependent interventions are active (specifically, when the staffing surge term *Ã*∆¹e −*Ãλt*˜ is non-negligible), the instantaneous sensitivity of the system state (e.g., *Y*˜ (t)) to parameter changes will differ from these calculated equilibrium sensitivities. The parameters governing the transient enhancement (*Ã*∆¹, ü) directly shape the system's trajectory during that phase but do not influence the final asymptotic state *Y*˜ ∗ upon which this sensitivity analysis is based.
Complementary dimensionless simulations (Eqs. (S11)–(S15)) were performed to provide generalized insights into the system's structural response to interventions and dynamic changes, independent of absolute scaling. These studies focus on the impact of Non-Invasive Ventilation (NIV) related parameters. Fixed dimensionless parameters for baseline conditions are noted in the figures unless specified otherwise.
We first analyzed the system's response to a temporary enhancement in the NIV recovery rate (ù(*t*˜) = ù0 + *Ã*∆¹e −*Ãλt*˜ ) under constant admission (ÃA = 1), mimicking interventions like shortterm staffing surges (Figure S2). The results quantitatively demonstrate the direct clinical benefits of improved NIV efficiency. Increasing the intervention magnitude (*Ã*∆¹) significantly lowers peak and transient normalized NIV occupancy (*Y*˜ , Panel B), thereby reducing the concurrent demand for this resource-intensive therapy. This upstream effect directly translates to improved patient outcomes: cumulative ICU transfers (*Z*˜, Panel C) and exitus (*W*˜ , Panel D) are substantially mitigated, while cumulative recovery (*R*˜, Panel E) is accelerated. The phase-space trajectories (Panel F) visually confirm that stronger interventions guide the system towards more favorable outcome states (lower Z, ˜ *W*˜ ; higher *R*˜). This underscores the clinical value proposition of the IRCU: targeted improvements in NIV care can dynamically reshape patient flow to optimize resource utilization and patient outcomes during the intervention period.
Figure S2| Impact of Transient NIV Recovery Enhancement (π∆*θ*) on Dimensionless IRCU Dynamics. (A-E) Normalized state variables (*X*˜-*R*˜) vs. time (*t*˜) comparing baseline (π∆θ = 0, dashed) with scenarios of increasing transient recovery boost (π∆θ > 0, solid lines). (F) Normalized phase-space trajectories (*Y*˜ -*Z*˜-*W*˜ ) colored by normalized recovery (*R*˜). Enhanced recovery reduces NIV load (B), mitigates adverse outcomes (C, D), and accelerates recovery (E). Fixed simulation parameters noted.
To understand performance under external pressure, we simulated the impact of varying baseline NIV effectiveness during a significant admission surge (ÃA(*t*˜) as a Gaussian peak; Figure S3). Three effectiveness levels (Less Effective - LE, Baseline - BE, Highly Effective - HE), defined by concurrently varying *ù0, õ, Ãϵ* (parameters in figure caption), were compared. While the admission surge dominates the temporal pattern of active populations (X, ˜ *Y*˜ , Panels B, C), the simulation clearly shows that higher NIV effectiveness (HE) acts as a crucial buffer. It markedly attenuates the negative consequences of the surge on cumulative adverse outcomes, leading to substantially lower ICU transfers (*Z*˜, Panel D) and exitus (*W*˜ , Panel E), while maximizing recoveries (*R*˜, Panel F) compared to the LE and BE scenarios [18, 19]. Such effectiveness is recognized to be multifactorial in clinical practice, depending on staff proficiency, protocols, patient selection, and timeliness [20, 21]. This simulation thus highlights that an IRCU's value lies not merely in capacity [22], but significantly in its capability for high-quality care delivery which actively improves outcomes and alleviates downstream ICU pressure during crises [23].
Figure S3| Influence of NIV Effectiveness Level during an Admission Surge. (A) Normalized Gaussian admission surge *A*˜(t). (B-F) Resulting normalized state dynamics (*X*˜-*R*˜) for Less Effective (LE, dotted), Baseline (BE, solid), and Highly Effective (HE, dashed) NIV parameter sets (definitions in Panel B legend). Higher effectiveness substantially buffers the surge's impact, reducing adverse outcomes (D, E) while maximizing recovery (F).
We also assessed system resilience to internal fluctuations by simulating a time-varying *NIV recovery rate* (ù(t), assuming this is Panel A's variable based on parameters) under constant admission (ÃA = 1) (Figure S4). Comparing LE, BE, and HE effectiveness levels reveals that higher effectiveness not only yields consistently better outcomes (lower W, ˜ *Z*˜; higher *R*˜) but also enhances system stability. The HE scenario shows dampened oscillations in NIV occupancy (*Y*˜ , Panel C) and mitigated outcome disruptions following the transient change in ù(t). This suggests that high-quality NIV care, influenced by clinical expertise, technology, and management strategies [24, 25, 26], is critical not only for managing external load but also for maintaining robust performance despite internal system variability or policy shifts (e.g., evolving ICU admission criteria).
Figure S4| System Response to a Time-Varying NIV Recovery Rate (πθ(t)) across Effectiveness Levels. (A) Profile of the simulated transient change in normalized parameter πθ(t) under constant admission *A*˜ = 1. (B-F) Comparison of resulting normalized state dynamics (*X*˜-*Z*˜) for Less Effective (LE, dotted), Baseline (BE, solid), and Highly Effective (HE, dashed) NIV effectiveness levels. Higher effectiveness improves overall outcomes and enhances system stability against internal perturbations.
Finally, recognizing that interventions during surges often involve resource tradeoffs, we simulated the combined effect of simultaneously varying ICU transfer mitigation (*µ*0, 0 to 0.5) and NIV mortality reduction (*ε*0, 0 to 0.5) effectiveness under admission surge conditions (Figure S5). Plotting the resulting cumulative ICU transfers averted versus cumulative exitus averted reveals the achievable efficiency landscape. The distribution forms a Pareto-like frontier [27], visually representing the tradeoff: prioritizing µ0 (larger points) maximally averts ICU transfers, while prioritizing ε0 (brighter colors) maximally averts exitus events. This quantitative tradeoff analysis provides a valuable tool for supporting evidence-informed strategic decision-making [28]. Depending on the primary constraint or objective during a crisis (e.g., preserving ICU capacity vs. minimizing mortality), IRCU leadership can use such insights to target quality improvement efforts or allocate resources more effectively towards protocols impacting either µ0 or ε0 [29].
Figure S5| Tradeoffs Between Averting ICU Transfers and Exitus via Combined Interventions During a Surge. Efficiency landscape plotting cumulative ICU transfers averted versus cumulative exitus averted (relative to baseline surge 'X') when simultaneously varying ICU mitigation effectiveness (*γ*0, point size, 0-0.5) and mortality reduction effectiveness (*ε*0, color, 0-0.5). The frontier illustrates achievable outcomes and the inherent tradeoffs for resource prioritization.
Accurate forecasting of patient occupancy is essential for effective resource management within Intermediate Respiratory Care Units (IRCUs), particularly given the inherent variability in patient demand and discharge patterns commonly observed in healthcare settings [30, 31]. While traditional queuing models offer analytical tractability, their reliance on strong parametric assumptions (e.g., exponential LOS, Poisson arrivals) often limits their ability to capture the complexities of real-world systems, such as non-standard LOS distributions or time-varying admission rates [31, 32, 33, 34, 35]. To overcome these limitations, we developed a data-driven, probabilistic modeling framework for IRCU occupancy grounded in empirical observations. This approach integrates a non-parametric representation of the Length-of-Stay (LOS) distribution with flexible, stochastic modeling of the admission process using Gaussian Processes (GPs), allowing for robust occupancy prediction and uncertainty quantification.
The model fundamentally links daily occupancy to past admissions and the duration patients stay. Let d be the LOS in days. Rather than assuming a specific parametric form, we directly characterize the LOS distribution using the empirical Probability Mass Function (PMF) derived from observed patient data:
$$f(d) = \frac{\text{Number of patients with observed LOS exactly } d}{\text{Total number of observed patients}},$$ (S16)
calculated for all observed LOS durations d up to a maximum *d*max. This non-parametric method preserves the potentially complex shape (e.g., multimodality, heavy tails) of the observed LOS distribution [7, 36].
From the PMF, the corresponding empirical survival function S(d), representing the probability that a patient admitted on day s is still present in the unit d days later (i.e., P(LOS g d)), is calculated as:
$$S(d) = 1 - \sum_{k=0}^{d-1} f(k). \tag{S17}$$
Note that S(0) = 1 by definition.
Let A(t) be the number of new admissions to the IRCU on day t. The total IRCU occupancy X(t) on day t represents the accumulation of patients admitted on or before day t who have not yet been discharged. This relationship is formalized through the discrete convolution of the admission history with the survival function:
$$X(t) = \sum_{s=0}^{t} A(s) S(t-s).$$ (S18)
Under steady-state conditions where the admission rate A(t) is constant (A(t) = A), this convolution leads to the expected equilibrium occupancy X = A × E[LOS], where E[LOS] is the mean LOS calculated from f(d), linking the model to fundamental queuing principles [16, 9].
Alternatively, the occupancy dynamics can be expressed recursively via a daily flow balance:
$$X(t+1) = X(t) + A(t+1) - L(t+1), (S19)$$
$$L(t+1) = \sum{d=1}^{d{\text{max}}} f(d) A(t+1-d) \mathcal{L}{t+1-d \ge 0},$$ (S20)
where A(t + 1) are the new admissions during day t + 1, L(t + 1) are the discharges during day t + 1 (calculated as the sum of patients admitted d days prior who are discharged on day t + 1, determined by the PMF f(d)), and X(t) is the occupancy at the start of day t. The indicator function ⊮{·} ensures causality.
Recognizing that daily hospital admissions A(t) are often complex and time-varying, we model them explicitly as draws from a Gaussian Process (GP), a flexible, non-parametric statistical tool suitable for time series [37]. Formally, we represent the admission on day t as:
$$A(t) \sim \mathcal{GP}(m(t), k(t, t')),$$ (S21)
where m(t) is the mean function (often set to zero or a simple trend a priori) and k(*t, t*′ ) is the covariance or kernel function, which defines the structure and smoothness of the admission patterns. The GP framework inherently treats the underlying admission function as a stochastic process, providing not only a predictive mean µA(t) after conditioning on observed data, but also principled uncertainty quantification via a predictive variance à 2 A(t).
The kernel function k(*t, t*′ ) specifies the correlation between admissions at different times t and t ′ . A common initial choice for capturing smooth trends is the Radial Basis Function (RBF) kernel combined with a White Noise kernel to account for independent observation uncertainty:
$$k(t,t') = C \cdot \exp\left(-\frac{(t-t')^2}{2\ell^2}\right) + \sigman^2 \,\delta{t,t'},\tag{S22}$$
where C is the amplitude variance, ℓ is the characteristic length-scale governing how quickly correlations decay with time separation, and à 2 n is the noise variance. These hyperparameters (collectively ¹ = {C, ℓ, Ã2 n}) are typically learned from observed admissions data D = {(ti , Ai)} by maximizing the marginal likelihood.
To better capture specific temporal patterns often present in healthcare admissions, such as weekly seasonality, composite kernels can be constructed by summing basic kernels [38]. For example, incorporating a periodic kernel allows explicit modeling of recurring weekly patterns (T = 7 days):
$$k_{\text{comp}}(t, t') = C_1 \cdot \text{RBF}(t, t') + C_2 \cdot \text{Periodic}(t, t'; T = 7) + \sigman^2 \delta{t, t'}, \tag{S23}$$
where *C1, C2,* and parameters within the Periodic kernel (like its own length-scale) are also learned from data. Further enhancements for capturing complex dependencies can involve alternative kernel choices or non-zero structured mean functions m(t), potentially informed by traditional time-series models [38].
By fitting the chosen GP model (i.e., selecting a kernel and optimizing hyperparameters) to historical admission data, we obtain a posterior distribution over the admission function A(t). Sampling trajectories from this posterior and simulating the occupancy dynamics using Eq. (S18) or (S19) for each sample allows the propagation of admission uncertainty through to the occupancy predictions, yielding probabilistic forecasts X(t) with associated credible intervals. This integrated approach leverages the empirical LOS distribution and flexible GP admission modeling to generate occupancy forecasts reflecting system nonlinearities and uncertainties, crucial for robust resource planning [39]. Furthermore, the framework allows for potential integration with real-time data streams and adaptive methods (e.g., Bayesian updating, AI techniques) to dynamically refine f(d) or GP predictions as new information becomes available [12], enhancing its practical utility in dynamic healthcare environments.
The baseline parameters for the ODE model (Eqs. (S6)-(S10)), representing average transition rates [days−1 ], were calibrated to reflect the observed dynamics of the UHVN IRCU cohort (n=249).
Methodology: The estimation combined observed outcome proportions with assumed average lengths of stay (Ä ) in the origin states, using the relationship: **Rate ≈ Proportion / Ä**.
(Avg. time in State Y before transition). This is a necessary simplifying assumption for this modeling approach.
4. Rate Calculation: ³ ≈ *PX→Y /ÄX* ≈ 0*.3092/7 ≈ 0.*04418 days−1 , Ä ≈ *PX*→R/ÄX ≈ 0*.6908/7 ≈ 0.*09868 days−1 , µ ≈ PY →Z/ÄY ≈ 0*.2338/10 = 0.*02338 days−1 , ε ≈ PY →W /ÄY ≈ 0*.0909/10 ≈ 0.*00909 days−1 , ¹0 ≈ PY →R/ÄY ≈ 0*.6753/10 ≈ 0.*06753 days−1 .
The resulting baseline parameters, consistent with total outflow rates (³ + Ä) ≈ 1/ÄX and (µ +ε + *¹0) ≈ 1/ÄY* , are summarized in Table S1.
Table S1| Baseline ODE model parameters estimated from cohort data and analysis assumptions.
| Parameter | Description | Estimated Value [days−1 ] |
|---|---|---|
| ³ | Rate X → Y (NIV Initiation) | 0.04418 |
| Ä | Rate X → R (Direct Recovery) | 0.09868 |
| ¸ | Rate X → Z (Direct ICU) | 0 (Fixed by observation) |
| ¿ | Rate X → W (Direct Exitus) | 0 (Fixed by observation) |
| µ | Rate Y → Z (ICU from NIV) | 0.02338 |
| ε | Rate Y → W (Exitus from NIV) | 0.00909 |
| ¹0 | Rate Y → R (Recovery from NIV) | 0.06753 |
| Parameters assumed for estimation: | ||
| ÄX | Assumed Avg. LOS in State X | 7 days |
| ÄY | Assumed Avg. LOS in State Y | 10 days |
Note: These are baseline estimates for the autonomous model. For simulations involving time-varying recovery (e.g., Main Paper Fig. 3G-I, J-L), the rate θ(t) was modified as described. Intervention modulators *γ0, ε*0 were set to 0 for baseline.
Table S2| Summary of the model's dynamic state variables.
| Variable | Description [Unit] |
|---|---|
| X(t) | Patients in IRCU, not receiving NIV at time t [patients] |
| Y (t) | Patients in IRCU, actively receiving NIV at time t [patients] |
| Z(t) | Patients transferred from IRCU pathway to ICU by time t [patients] |
| W(t) | Exitus from IRCU pathway (X or Y) by time t [patients] |
| R(t) | Recovered (discharged/to ward) from IRCU pathway (X or Y) by time t [patients] |
Table S3| Summary of the model parameters including intervention modulators.
| Parameter | Description [Unit] | |
|---------------------------------------------------|------------------------------------------------------------------|--|
| Admission | | |
| A(t) | Influx rate of new patient admissions to IRCU [patients/time] | |
| | Transitions from Non-NIV State (X) | |
| ³ | Rate of NIV initiation (X → Y) [time−1
] | |
| Ä | Direct recovery/discharge rate (from X) [time−1
] | |
| ¸ | Direct ICU transfer rate (from X) [time−1
] | |
| ¿ | Direct exitus rate (from X) [time−1
] | |
| | Baseline Transitions from NIV State (Y) | |
| µ | Baseline rate of ICU transfer (from Y) [time−1
] | |
| ε | Baseline rate of exitus (from Y) [time−1
] | |
| ¹0 | Baseline rate of recovery (from Y) [time−1
] | |
| Intervention Modulators (acting on Y transitions) | | |
| µ0 | ICU transfer mitigation factor (effectiveness) [dimensionless] | |
| ε0 | Exitus reduction factor (effectiveness) [dimensionless] | |
| ∆¹ | Initial magnitude of transient recovery enhancement [time−1
] | |
| ¼ | Decay rate of transient recovery enhancement [time−1
] | |
[8] M. J. Keeling and P. Rohani. Modeling Infectious Diseases in Humans and Animals. Princeton University Press, 2008. doi: 10.1515/9781400841035.
[9] O. Diekmann, H. Heesterbeek, and T. Britton. Mathematical Tools for Understanding Infectious Disease Dynamics. Princeton University Press, 2013.
[10] T. Al-Karkhi and K. Byatt. "A compartmental model to describe acute medical in-patient flow through a hospital". In: Heliyon 11.3 (2025), e42260. doi: 10.1016/j.heliyon.2025.e42260.
[11] J. Garrido, D. Martínez-Rodríguez, F. Rodríguez-Serrano, et al. "Mathematical model optimized for prediction and health care planning for COVID-19". In: Medicina Intensiva (English Edition) 46.5 (2022), pp. 248–258. doi: https://doi.org/10.1016/j.medine.2022.02.020.
[12] H. T. Banks and K. Kunisch. Estimation Techniques for Distributed Parameter Systems. Birkhäuser, 1989. doi: 10.1007/978-1-4612-3700-6.
[13] O. Diekmann and J. A. P. Heesterbeek. Mathematical Epidemiology of Infectious Diseases: Model Building, Analysis and Interpretation. John Wiley & Sons, 2000.
[14] F. Brauer. "Mathematical epidemiology: Past, present, and future". In: Infectious Disease Modelling 2.2 (2017), pp. 113–127. doi: 10.1016/j.idm.2017.02.001.
[15] L. J. S. Allen. An Introduction to Stochastic Processes with Applications to Biology. Pearson Education, 2003.
[16] S. M. Ross. Introduction to Probability Models. 11th. Academic Press, 2014.
[17] N. Bacaër. A Short History of Mathematical Population Dynamics. Springer, 2011. doi: 10. 1007/978-0-85729-115-8.
[18] B. Rochwerg, L. Brochard, M. W. Elliott, et al. "Official ERS/ATS clinical practice guidelines: noninvasive ventilation for acute respiratory failure". In: European Respiratory Journal 50.2 (2017), p. 1602426. doi: 10.1183/13993003.02426-2016.
[19] R. Scala and L. Pisani. "Non-invasive ventilation in acute respiratory failure: Which recipe?" In: European Respiratory Review 27.150 (2018), p. 180029. doi: 10.1183/16000617.0029-2018.
[20] E. Ozyilmaz, A. O. Ugurlu, and S. Nava. "Timing of noninvasive ventilation failure: causes, risk factors, and potential remedies". In: BMC Pulmonary Medicine 14.1 (2014), p. 19. doi: 10.1186/1471-2466-14-19.
[21] I. Ruzsics, P. Matrai, P. Hegyi, et al. "Noninvasive ventilation improves the outcome in patients with pneumonia-associated respiratory failure: Systematic review and meta-analysis". In: Journal of Infection and Public Health 15.3 (2022), pp. 349–359. doi: 10.1016/j.jiph.2022.02.004.
[22] M. Prin and H. Wunsch. "The role of stepdown beds in hospital care". In: American Journal of Respiratory and Critical Care Medicine 190.11 (2014), pp. 1210–1216. doi: 10.1164/rccm. 201406-1117PP.
[23] S. A. Nasraway, I. L. Cohen, R. C. Dennis, et al. "Guidelines on admission and discharge for adult intermediate care units". In: Critical Care Medicine 26.3 (1998), pp. 607–610. doi: 10.1097/00003246-199803000-00039.
[24] A. M. Esquinas, ed. Noninvasive Mechanical Ventilation. Theory, Equipment, and Clinical Applications. Publicado el 1 de enero de 2010. Springer Berlin Heidelberg, 2010, p. 402. doi: 10.1007/978-3-642-11365-9.
[25] A. M. Esquinas, ed. Noninvasive Mechanical Ventilation. Theory, Equipment, and Clinical Applications. 2nd ed. Springer Cham, 2016, p. 959. doi: 10.1007/978-3-319-21653-9.
[26] A. M. Esquinas, ed. Noninvasive Mechanical Ventilation. Theory, Equipment, and Clinical Applications. 3rd ed. Publicado el 28 de agosto de 2023. Springer Cham, 2023, pp. XXXVII, 847. doi: 10.1007/978-3-031-28963-7.
[27] K. Deb. "Multi-Objective Optimization". In: Search Methodologies: Introductory Tutorials in Optimization and Decision Support Techniques. Ed. by E. K. Burke and G. Kendall. Boston, MA: Springer US, 2005, pp. 273–316. doi: 10.1007/0-387-28356-0\_10.
[28] M. L. Brandeau, F. Sainfort, and W. P. Pierskalla, eds. Operations Research and Health Care. A Handbook of Methods and Applications. 1st ed. International Series in Operations Research & Management Science. New York, NY: Springer New York, 2004, pp. VIII, 874. doi: 10.1007/ b106574.
[29] C. Kahraman and Y. I. Topcu, eds. Operations Research Applications in Health Care Management. 1st ed. International Series in Operations Research & Management Science 255. Cham: Springer Cham, 2018, pp. XV, 604. doi: 10.1007/978-3-319-65455-3.
[30] E. Litvak and M. Bisognano. "More patients, less payment: increasing hospital efficiency in the aftermath of health reform". In: Health Affairs (Millwood) 30.1 (2011), pp. 76–80. doi: 10.1377/hlthaff.2010.1114.
[31] R. Hall, ed. Patient Flow: Reducing Delay in Healthcare Delivery. 2nd ed. International Series in Operations Research & Management Science. Springer New York, NY, 2013, pp. XII, 553. doi: 10.1007/978-1-4614-9512-3.
[32] J. B. Jun, S. H. Jacobson, and J. R. Swisher. "Application of Discrete-Event Simulation in Health Care Clinics: A Survey". In: The Journal of the Operational Research Society 50.2 (1999), pp. 109–123. doi: 10.2307/3010560.
[33] D. Fone, S. Hollinghurst, M. Temple, et al. "Systematic review of the use and value of computer simulation modelling in population health and health care delivery". In: Journal of Public Health 25.4 (2003), pp. 325–335. doi: 10.1093/pubmed/fdg075.
[34] S. H. Jacobson, S. N. Hall, and J. R. Swisher. "Discrete-Event Simulation of Health Care Systems". In: Patient Flow: Reducing Delay in Healthcare Delivery. Ed. by R. W. Hall. Boston, MA: Springer US, 2006, pp. 211–252. doi: 10.1007/978-0-387-33636-7\_8.
[35] J. D. Griffiths, N. Price-Lloyd, M. Smithies, et al. "A queueing model of activities in an intensive care unit". In: IMA Journal of Management Mathematics 17.3 (2006), pp. 277–288. doi: 10. 1093/imaman/dpi042.
[36] J. M. Heffernan, R. J. Smith, and L. M. Wahl. "Perspectives on the basic reproductive ratio". In: Journal of the Royal Society Interface 2.4 (2005), pp. 281–293. doi: 10.1098/rsif.2005.0042.
[37] C. E. Rasmussen and C. K. I. Williams. Gaussian Processes for Machine Learning. The MIT Press, 2005. doi: 10.7551/mitpress/3206.001.0001.
[38] S. Roberts, M. Osborne, M. Ebden, et al. "Gaussian processes for time-series modelling". In: Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 371.1984 (2013), p. 20110550. doi: 10.1098/rsta.2011.0550.
[39] W. Whitt and X. Zhang. "Forecasting arrivals and occupancy levels in an emergency department". In: Operations Research for Health Care 21 (2019), pp. 1–18. doi: 10.1016/j.orhc. 2019.01.002.