This paper investigates a fractional-order mathematical model for the co-infection dynamics of pneumonia and typhoid fever using the Liouville–Caputo derivative. We establish the existence, uniqueness, non-negativity, and boundedness of solutions using Banach’s fixed point theorem and fractional comparison principles. The Hyers–Ulam and generalized Ulam–Hyers–Rassias stability of the system are rigorously proved; this stability analysis is epidemiologically significant because it guarantees that small perturbations in initial conditions or model parameters—inevitable in real-world data collection—do not lead to unbounded deviations in disease trajectory predictions. To approximate solutions numerically, we develop a Laplace-Based Optimized Decomposition Method (LODM) and validate its convergence against a modified predictor–corrector scheme. The LODM provides a semi-analytical series solution, while the predictor–corrector method serves as a numerical benchmark; this dual approach ensures reliability of simulations. Numerical simulations illustrate the influence of the fractional order ξ on system dynamics. Quantitative comparison between ξ=1 (integer order) and ξ<1 (fractional order) demonstrates that fractional modeling reduces peak infection by 12–18% and delays epidemic peaks by 15–30 days, confirming that memory effects capture long-term epidemiological dependencies that integer-order models fail to reproduce. A biological interpretation links the fractional order to immune memory, pathogen persistence, and intervention latency. This study provides both theoretical and numerical evidence supporting the use of fractional calculus in epidemiological modeling.
Alhamzi et al. (2026) studied this question.