A Nine-Compartment Nonlinear Epidemic Model with Spline-Based Identification of Time-Varying Transmission and Vaccination Dynamics: Application to the COVID-19 Third Wave in Italy

Abstract: We develop a nine-compartment nonlinear epidemic model incorporating two co-circulating viral strains (ancestral I1 and the Alpha variant B.1.1.7 I2, which is 43-90% more transmissible, c2=1.5), a super-spreader subpopulation, partial vaccine-induced immunity with waning, and explicit hospitalization dynamics with differentiated mortality. Transmission and vaccination rates are treated as time-varying control inputs and identified from Italian COVID-19 data (January-May 2021) via a Piecewise Cubic Hermite Interpolating Polynomial (PCHIP) control-node parameterization, reducing calibration to a fourteen-variable Sequential Quadratic Programming (SQP) problem with monotonicity and box constraints. A parametric bootstrap (n=1000) quantifies parameter uncertainty. The calibrated model achieves R^2=0.966 for active hospitalizations, R^2=0.987 for cumulative fatalities, and R^2=0.999 for cumulative vaccinations. Well-posedness, the basic reproduction number in closed form, and local and global stability of the disease-free equilibrium are established analytically. An L-infinity approximation error bound shows that the PCHIP control-node parameterization converges to the true time-varying parameters at rate O(h^2) as the node spacing vanishes. Local identifiability and a noise stability bound are established via the Fisher information matrix. A sufficient threshold condition proves epidemic decay under time-varying suppression whenever the effective reproduction number remains persistently below one. Sensitivity analyses consistently rank hospital throughput parameters above the transmission rate, providing a mathematical basis for the observation that reactive containment measures cannot prevent a hospitalization peak already driven by the pre-existing latent viral load.

• The paper presents a nine-compartment ODE model that differentiates SARS-CoV-2 variants, super-spreader events, and vaccination dynamics using a spline-based identification approach.
• The paper employs a low-dimensional PCHIP spline method to calibrate time-dependent transmission and vaccination rates, achieving high fidelity with R² values over 0.96 in Italian COVID-19 data.
• The paper demonstrates that hospital throughput parameters, rather than immediate transmission controls, predominantly drive peak hospital burden, guiding effective pandemic interventions.

Spline-Based Identification of Time-Varying Transmission and Vaccination in a Nine-Compartment COVID-19 Model

Model Architecture and Innovations

This paper presents a nine-compartment ODE epidemic model structured to discriminate between co-circulating SARS-CoV-2 strains, super-spreader events, explicit hospitalization/recovery/fatality pathways, partial (waning) immunity after vaccination, and time-dependent intervention effects. The compartmental topology is depicted in Figure 1, capturing direct and indirect flows, with two explicit control inputs governing force of infection and vaccination rollout.

Figure 1: Transition diagram of the nine-compartment epidemic model; solid arrows denote transitions, dashed arrows denote force-of-infection terms, and $\beta(t), w_1(t)$ are time-varying controls.

The system extends classical SEIR-type formulations by: (i) parametrizing the competitive displacement of an ancestral and a more transmissible variant (B.1.1.7/Alpha), (ii) modeling a small, explicitly identified super-spreader compartment with an elevated force-of-infection multiplier, and (iii) accommodating partial vaccine protection, vaccine waning, and explicit hospitalization/fatality stratification. All transmission pathways employ frequency-dependent transmission. By scaling variant and super-spreader infection rates by clinical literature parameters ($c_2=c_P=1.5$), structural parameter identifiability is achieved even in the presence of unobserved compartments.

Spline-Based Parameter Identification Framework

The identification of time-dependent parameters $\beta(t)$ (transmission rate) and $w_1(t)$ (vaccination rate) is formalized as an ODE-constrained nonlinear inverse problem. Rather than fitting daily parameters (leading to severe overfitting and loss of identifiability), the trajectories are parameterized using a low-dimensional PCHIP spline at five control nodes. The model calibration then reduces to optimizing 14 parameters (10 PCHIP nodes, 4 scalars for initial conditions and hospital fatality rate), subject to monotonicity and box constraints, via SQP with multi-start and parametric bootstrap. This design leverages the local monotonicity and non-negativity of PCHIP and ensures fast and robust convergence compared to MCMC or high-dimensional unsmoothed approaches.

The authors provide a formal $L^\infty$ approximation error bound, showing $O(h^2)$ convergence of the PCHIP parameterization to the underlying smooth time-series as node density increases, demonstrating that the representational error is dominated by measurement uncertainty at their chosen granularity.

Analytical Properties: Stability, Reproduction Numbers, Identifiability

Rigorous analysis is given for model well-posedness (existence, uniqueness, and non-negativity of solutions), and the full next-generation matrix for the disease subsystem yields a closed-form $\mathcal{R}_0$ decomposing infection into original, variant, and super-spreader contributions. The disease-free equilibrium is proven locally and globally stable whenever $\mathcal{R}_0<1$, with explicit expressions for the effective reproduction number under time-varying input policies and vaccination coverage.

The identifiability analysis is two-pronged: a structural algebraic argument for the constant-parameter restriction, and, crucially, a local identifiability and noise-stability result for the full nonlinear time-varying parameterization, established via positive-definiteness of the Fisher information matrix at the solution and further confirmed by bootstrap. The main theoretical result here is a bound on the solution perturbation induced by observational noise.

Empirical Results: Calibration, Goodness-of-Fit, and Pandemic Dynamics

The model is calibrated to national Italian data for hospitalizations, fatalities, and vaccinations during January–May 2021 (the third/Alpha-driven wave and initial vaccination campaign), with data preprocessing to mitigate reporting artifacts and calibration restricted to avoid known protocol discontinuities. Smoothed data is shown in Figure 2, and the fit to these targets alongside reconstructed $\beta(t),w_1(t)$ is shown in Figure 3.

Figure 2: Time-series of daily hospitalizations, cumulative fatalities, and vaccinations for Italy; shaded window denotes calibration interval (third wave and vaccine ramp-up).

Figure 3: Model fits (lines) to observed daily hospitalizations, cumulative fatalities, and cumulative vaccinations (points), with lower right showing PCHIP-estimated $\beta(t)$ and $c_2=c_P=1.5$0 trajectories with bootstrap $c_2=c_P=1.5$1 bands (vertical dashed lines: node positions).

The estimator achieves high fidelity ($c_2=c_P=1.5$2 for $c_2=c_P=1.5$3, $c_2=c_P=1.5$4 for $c_2=c_P=1.5$5, $c_2=c_P=1.5$6 for $c_2=c_P=1.5$7) and credible intervals for estimated parameters are narrow under repeated bootstrap. Numerical results demonstrate that the Alpha variant's per-contact transmission advantage and a highly nonlinear increase in daily vaccination capacity (from $c_2=c_P=1.5$8 to $c_2=c_P=1.5$9 doses/day) are both robustly identified and quantitatively consistent with external epidemiological data.

Importantly, when predicting held-out data (May out-of-sample), the predicted active hospitalizations retain strong fidelity ($\beta(t)$0, $\beta(t)$1 relative error), validating the spline-based approach as neither underfitting nor overfitting the underlying signals.

Sensitivity Analysis and Policy Implications

Sensitivity analysis, via both OAT and Morris EE methods, yields the robust and somewhat counter-intuitive finding that hospital throughput parameters dominate the effect on peak hospital burden relative to contemporaneous transmission rates. Specifically, the hospital discharge ($\beta(t)$2) and admission ($\beta(t)$3) rates exert the strongest influence on the hospitalization peak, while transmission-rate perturbations (even when temporally aligned) are strongly attenuated by the intrinsic epidemic lag structure. This quantifies that reactive intervention at a pandemic peak cannot significantly alter the hospital overload dynamic, as it is largely set by latent infections acquired in prior weeks.

Model Selection and Limitations

Model selection using AIC/BIC confirms that the 14-parameter, five-node PCHIP parameterization is optimally parsimonious under both in-sample and out-of-sample error criteria; increasing node count yields negligible further fit improvement but is penalized by BIC. The model's abstraction as a well-mixed national system limits applicability to regional/heterogeneous contexts, but internal analysis indicates that critical features (identifiability, forecast skill, policy conclusions) are not compromised by this simplification. Theoretical input-to-state stability for arbitrary time-varying policy inputs is not fully closed but a sufficient condition for epidemic extinction is proven.

Conclusions and Future Directions

This work establishes a rigorously analyzed, spline-based approach for time-dependent parameter identification in nonlinear epidemic models with explicit biological realism. Contrary to standard intuition, hospital throughput rates supersede NPIs/transmission rate controls in determining peak medical system risk, a finding with direct relevance to pandemic response strategy. The practical parameter identification technique, grounded in projection-based regularization and robust uncertainty quantification, is shown to be computationally efficient and empirically accurate. The framework naturally enables extension to multi-region or model-predictive control (MPC) architectures, where time-dependent $\beta(t)$4 and $\beta(t)$5 serve as data-driven proxies for optimal intervention planning.

The methodology and findings are broadly instructive for future applications in epidemic forecasting, intervention evaluation, and real-time policy adaptation using ODE-constrained optimization and nonlinear inverse problem theory (2606.07413).