Dynamical analysis and numerical simulation of an SIR epidemic model with control parameters

Abstract

Mathematical modeling is essential for understanding infectious disease transmission and evaluating control strategies. In this study, we investigate the classical susceptible-infected-recovered (SIR) epidemic model using both analytical and numerical approaches. The novelty of this work lies in the integrated presentation of analytical stability analysis alongside explicit numerical quantification of how vaccination and recovery rates affect epidemic peaks, supported by phase-plane interpretations. (Vaccination is incorporated as an external control parameter ν in the modified SIR model, where susceptible individuals are vaccinated at a constant rate ν and move directly to the recovered class). Unlike previous studies that focus solely on theoretical thresholds, we provide a practical, simulation-driven framework that directly links control parameters to outbreak outcomes. We analyze the equilibrium structure and stability properties of the disease-free equilibrium. The analysis shows that the basic reproduction number R₀ = β/γ acts as a threshold parameter governing disease spread. Using the numerical values β = 0.5, γ = 0.1, total population N = 1000, and initial infected I₀ = 1, we obtain R₀ = 5. When R₀ < 1, the infection disappears; when R₀ > 1, an outbreak occurs. Numerical simulations illustrate the temporal evolution and the influence of vaccination and recovery rates on epidemic progression. Our results quantify that increasing the vaccination rate from 0 to 0.3 reduces the infection peak by approximately 65%. Additionally, increasing the recovery rate by 50% reduces the infection peak by approximately 40%. These findings demonstrate the practical value of the model for public health planning.

1. Introduction

Mathematical modeling has become one of the most important tools for understanding the spread of infectious diseases and predicting epidemic outbreaks. Such models provide a systematic framework for analyzing the transmission mechanisms of diseases and evaluating the effectiveness of control strategies such as vaccination and treatment. Classical compartmental models divide the population into different epidemiological classes and describe the evolution of these classes using systems of differential equations.

Recent advances in stochastic epidemic modeling have introduced sophisticated numerical methods for studying disease dynamics. Arif et al. ¹ developed finite difference solutions for stochastic epidemics incorporating treatment cure and partial immunity. The same authors ² presented a reliable computational scheme for stochastic reaction-diffusion nonlinear chemical models. Baazeem et al. ³ proposed a robust computational approach for stochastic SIRS models with partial immunity and incidence rates. Shatanawi et al. ⁴ introduced an effective numerical method for a stochastic coronavirus (COVID-19) pandemic model. Additionally, Shatanawi et al. ⁵ developed essential features preserving dynamics for a stochastic dengue model. These recent developments complement our work and provide broader context for stochastic epidemic modeling.

Among the most widely used epidemic models is the susceptible–infected–recovered (SIR) model, originally introduced by Kermack and McKendrick. This model has been extensively used to describe the transmission dynamics of infectious diseases and to study the threshold behavior governed by the basic reproduction number.^6–8 The model captures the interaction between susceptible, infected, and recovered populations and explains how epidemics emerge, reach a peak, and eventually decline.

In recent years, many researchers have extended the classical epidemic models in order to incorporate more realistic biological mechanisms. These extensions include stochastic formulations, delay effects, exposed classes, and control strategies such as vaccination and isolation. Such approaches enable researchers to better understand complex epidemic phenomena and provide more accurate predictions for disease evolution.

For instance, several studies have investigated epidemic models that incorporate additional compartments or stochastic dynamics in order to describe disease propagation more realistically. In particular, the susceptible–exposed–infectious framework has been used to analyze epidemic transmission using probabilistic and Markov chain approaches. ⁹ Other studies have examined epidemic models with time delay and stability analysis in order to capture incubation periods and more complex disease dynamics. ¹⁰

These developments highlight the importance of mathematical modeling in studying epidemic dynamics and evaluating possible control strategies. In addition to classical epidemic models, recent advances in the stability analysis of dynamical systems have provided valuable tools for understanding complex biological and physical phenomena. For instance, Amer et al. ¹¹ investigated the stability optimization of a vibrating rigid body pendulum with an energy harvesting device. In contrast, Amer et al. ¹² analyzed the vibrational stability of planar double pendulum dynamics near resonance. These studies offer complementary insights into the behavior of dynamical systems under different conditions, which can be leveraged to improve the analysis of epidemic models. The SIR model and its variants have been widely applied to real-world infectious diseases such as COVID-19, influenza, and measles. For instance, the basic reproduction number for seasonal influenza typically ranges from R₀ = 1.2 to 1.8, while for measles it can reach R₀ = 12 − 18, and for the original SARS-CoV-2 strain it was estimated between 2.5 and 3.0. These values determine the threshold behavior and the level of intervention required to control outbreaks. Motivated by these works, the present study focuses on the classical SIR epidemic model and investigates its dynamical properties using analytical and numerical techniques. The equilibrium structure of the model is analyzed, and the stability of the disease-free state is examined. In addition, numerical simulations are performed to illustrate the temporal evolution of the epidemic and to analyze the influence of key epidemiological parameters such as vaccination and recovery rates. The main contribution of this work is to provide a comprehensive analytical and numerical investigation of the classical SIR epidemic model. In particular, the study presents a detailed stability analysis of the disease-free equilibrium and highlights the role of the basic reproduction number as a threshold parameter governing epidemic outbreaks. In addition, numerical simulations are performed to visualize the temporal evolution of the epidemic and to analyze the influence of important epidemiological parameters such as vaccination and recovery rates. The combination of theoretical analysis and numerical experiments provides a clear interpretation of epidemic behavior and offers useful insight into the design of effective disease control strategies. Recent advances in mathematical modeling of epidemic dynamics have been discussed by Sohaly. ¹³ Comprehensive treatments of mathematical models in epidemiology can be found in standard references. ¹⁴

Although the classical SIR model has been extensively studied since 1927, our work offers several distinctive contributions. First, we provide an integrated presentation that combines analytical stability analysis, numerical simulation, phase-plane visualization, stochastic comparison, and quantitative intervention effects in a single self-contained framework. Second, we present explicit quantitative results, such as the 65% reduction in infection peak achieved by increasing the vaccination rate to 0.3, which are rarely quantified in such detail in standard textbooks. Finally, we explicitly link model parameters to real-world diseases (influenza, measles, COVID-19) and discuss practical public health implications.

Unlike review papers that summarize existing results, this study provides original analytical and numerical contributions. Compared to classical deterministic SIR analyses [? ? ], our work adds explicit quantitative results for intervention effects (e.g., 65% peak reduction with vaccination). Compared to fractional-order epidemic models which incorporate memory effects, our integer-order model provides a simpler, more tractable framework while still capturing essential dynamics. The stochastic comparisons in this study (Figure 3) and the phase-plane visualizations (Figure 2) offer additional insights not always presented together in standard treatments. Our results are not claimed as ”benchmark results” but rather as illustrative simulations that can serve as a reference for parameter choices and expected outcomes.

The remainder of this paper is organized as follows. In Section 2, the mathematical formulation of the classical SIR epidemic model is presented and the governing system of differential equations is introduced. Section 3 is devoted to the analysis of equilibrium points and the derivation of the basic reproduction number, which plays a fundamental role in determining the epidemic threshold. In Section 4, the stability properties of the disease-free equilibrium are investigated using standard techniques from dynamical systems theory. Section 5 presents numerical simulations that illustrate the dynamical behavior of the model and demonstrate the influence of epidemiological parameters such as vaccination and recovery rates on the spread of infection. Finally, Section 6 discusses the obtained results and Section 7 concludes the paper with a summary of the main findings and possible directions for future research.

2. Mathematical formulation of the epidemic model

Mathematical epidemic models describe the evolution of infectious diseases by dividing the population into epidemiological compartments that represent different health states. One of the most fundamental and widely used frameworks is the SIR model, which classifies the population into three mutually exclusive groups.

Let the total population at time t be denoted by N(t). The population is partitioned into three compartments:

• S(t): the number of susceptible individuals who can contract the infection,

• I(t): the number of infected individuals capable of transmitting the disease,

• R(t): the number of recovered individuals who have acquired immunity.

Thus, the total population satisfies

$N(t)=S(t)+I(t)+R(t).$

The epidemic dynamics are governed by the following biological mechanisms:

1. Disease transmission occurs through contact between susceptible and infected individuals.

2. Infected individuals recover at a constant recovery rate.

3. Recovered individuals acquire permanent immunity and do not return to the susceptible class.

Under these assumptions, the deterministic SIR model is described by the following system of nonlinear differential equations:

The classical SIR model is described by the following system of differential equations:

$\frac{dS}{dt}=-\beta SI$

(1)

$\frac{dI}{dt}=\beta SI-\gamma I$

(2)

$\frac{dR}{dt}=\gamma I$

(3)

2.1. SIR model with vaccination

To study vaccination as a control parameter, we modify the classical SIR model by adding a vaccination term νS:

$\frac{dS}{dt}=-\beta SI-\nu S$

(4)

$\frac{dI}{dt}=\beta SI-\gamma I$

(5)

$\frac{dR}{dt}=\gamma I+\nu S$

(6)

Here, ν is the vaccination rate, and vaccinated individuals move directly from S to R. where

• β > 0 is the transmission rate of the infection,

• γ > 0 is the recovery rate of infected individuals.

The first equation describes the decrease in the susceptible population due to infection. The second equation represents the balance between new infections and recoveries, while the third equation captures the accumulation of recovered individuals (Tables 1 and 2).

Table 1.

Parameters and variables used in the SIR model.

Symbol	Description	Typical value	Unit
S(t)	Susceptible individuals	1000	population
I(t)	Infected individuals	1	population
R(t)	Recovered individuals	0	population
N	Total population	1000	population
β	Transmission rate	0.5	day−1
γ	Recovery rate	0.1	day−1
ν	Vaccination rate	0 - 0.3	day−1
R 0	Basic reproduction number	β/γ	dimensionless

Table 2.

Default numerical parameters used in all simulations. These values yield R₀ = β/γ = 5, producing a clear epidemic outbreak for illustrative purposes.

Parameter	Value	Description
β	0.5	Transmission rate
γ	0.1	Recovery rate
N	1000	Total population
I 0	1	Initial infected individuals
S 0	999	Initial susceptible individuals
R 0	0	Initial recovered individuals
t max	100	Simulation time (days)

An important property of the system is the conservation of the total population. By summing the three differential equations we obtain

$\frac{dN}{dt}=\frac{dS}{dt}+\frac{dI}{dt}+\frac{dR}{dt}=0,$

which implies that

$N(t)=N_0,$

where N₀ is the constant total population size.

For analytical convenience, the model can be expressed in normalized form by dividing each compartment by the total population. Let

$s(t)=\frac{S(t)}{N_0},i(t)=\frac{I(t)}{N_0},r(t)=\frac{R(t)}{N_0}.$

The normalized variables satisfy

$s(t)+i(t)+r(t)=1.$

Substituting these expressions into the model yields the normalized system

$\frac{ds}{dt}=-\beta si,$

$\frac{di}{dt}=\beta si-\gamma i,$

$\frac{dr}{dt}=\gamma i.$

2.1. Model assumptions

The proposed SIR epidemic model is based on the following key assumptions:

1. Closed population: The total population size N is constant over time, with no births, deaths, or migration. This implies dN/dt = 0.

2. Homogeneous mixing: Every individual in the population has an equal probability of coming into contact with any other individual. This assumption simplifies the transmission term βSI.

3. Permanent immunity: Recovered individuals acquire lifelong immunity and cannot return to the susceptible class. There is no waning immunity or reinfection.

4. Constant parameters: The transmission rate β and recovery rate γ are constant and do not change over time or vary with age, behavior, or season.

5. No incubation period: Infected individuals become infectious immediately upon infection. There is no exposed/latent compartment.

6. No intervention initially: In the basic model (before introducing control parameters), there is no vaccination, treatment, or isolation.

This normalized representation simplifies the mathematical analysis and highlights the intrinsic structure of the epidemic dynamics.

2.3. Incorporating vaccination as a control parameter

To study the effect of vaccination as an intervention strategy, we modify the classical SIR model by adding a vaccination term. Vaccination is introduced as an external control parameter ν, which represents the rate at which susceptible individuals are vaccinated. The modified system of differential equations becomes:

$\frac{dS}{dt}=-\beta SI-\nu S,$

$\frac{dI}{dt}=\beta SI-\gamma I,$

$\frac{dR}{dt}=\gamma I+\nu S.$

Here, susceptible individuals are vaccinated at rate ν and move directly to the recovered class, where they acquire immunity. The total population remains conserved since dN/dt = dS/dt + dI/dt + dR/dt = 0. This formulation allows us to quantitatively investigate how different vaccination rates affect the epidemic peak and overall disease spread.

2.4. Exact differential equations with vaccination

The classical SIR model is extended to include vaccination as a control parameter. The exact system of differential equations governing the dynamics is:

$\begin{array}{cc}\hfill \frac{dS}{dt}& =-\beta SI-\nu S,\hfill \\ \hfill \frac{dI}{dt}& =\beta SI-\gamma I,\hfill \\ \hfill \frac{dR}{dt}& =\gamma I+\nu S,\hfill \end{array}$

where:

• • β is the transmission rate,

• • γ is the recovery rate,

• • ν is the vaccination rate (control parameter).

In this formulation, vaccination is not modeled as a separate compartment. Instead, susceptible individuals are vaccinated at a constant rate ν and move directly from the S compartment to the R (recovered/immune) compartment. This assumes that vaccination confers immediate and permanent immunity. This approach is simpler than adding a separate vaccinated compartment (SIR-V) and is appropriate for diseases where vaccination provides full protection without an intermediate stage.

3. Equilibrium analysis and basic reproduction number

The equilibrium points of the epidemic system correspond to the states where the disease dynamics remain constant over time. Mathematically, equilibrium solutions are obtained by setting the time derivatives of the system equal to zero. For the normalized SIR system:

$\frac{ds}{dt}=-\beta si,$

$\frac{di}{dt}=\beta si-\gamma i,$

$\frac{dr}{dt}=\gamma i,$

an equilibrium occurs when

$\frac{ds}{dt}=0,\frac{di}{dt}=0,\frac{dr}{dt}=0.$

From the third equation we obtain

$\gamma i=0,$

which implies

$i=0.$

Substituting this result into the first equation gives

$-\beta s\cdot 0=0,$

which is automatically satisfied. Since the normalized variables satisfy

$s+i+r=1,$

the equilibrium state becomes

$(s^{*},i^{*},r^{*})=(s^{*},0,1-s^{*}).$

A particularly important equilibrium is the disease–free equilibrium (DFE), obtained when the entire population is susceptible and no infection exists in the system:

$E_0=(1,0,0).$

This equilibrium represents the absence of the disease in the population. To determine whether an epidemic can invade the population, we examine the growth behavior of infected individuals during the early stage of the outbreak. When the infection is initially rare, the susceptible population is approximately equal to one, i.e.,

$s\approx 1.$

Substituting this approximation into the infected equation yields

$\frac{di}{dt}\approx (\beta -\gamma )i.$

The sign of the growth rate depends on the quantity

$\beta -\gamma .$

This leads to the definition of the basic reproduction number

$R_0=\frac{\beta }{\gamma }.$

The parameter R₀ represents the expected number of secondary infections generated by a single infected individual in a fully susceptible population.

3.1. Basic reproduction number R₀

For the classical SIR model (without vaccination), the basic reproduction number is derived from the early growth of infection. When s ≈ 1, the infected equation becomes di/dt ≈ (β − γ)i, leading to:

$R_0=\frac{\beta }{\gamma }.$

For the modified SIR model with vaccination, the control reproduction number is obtained by considering the loss of susceptible individuals due to both infection and vaccination. The effective growth rate becomes (β − γ − ν), yielding:

$R_c=\frac{\beta }{\gamma +\nu }.$

When R _c < 1, the infection dies out; when R _c > 1, an outbreak occurs. This expression shows how vaccination reduces the effective reproduction number and can bring an epidemic under control even when the basic R₀ > 1. Two different epidemic regimes arise:

• If R₀ < 1, then di/dt < 0 and the infection decays exponentially. In this case the disease-free equilibrium is stable and the epidemic cannot invade the population.

• If R₀ > 1, then di/dt > 0 and the number of infected individuals initially grows. This condition leads to the emergence of an epidemic outbreak.

Therefore, the basic reproduction number serves as a threshold parameter that determines whether the infection dies out or spreads through the population. This threshold phenomenon is a fundamental property of nonlinear epidemic models.

4. Stability analysis of the disease-free equilibrium

To investigate the local stability of the disease-free equilibrium, we analyze the Jacobian matrix of the epidemic system. Consider the normalized SIR model

$\frac{ds}{dt}=-\beta si,$

$\frac{di}{dt}=\beta si-\gamma i,$

$\frac{dr}{dt}=\gamma i.$

Since the third equation does not influence the first two, the dynamical behavior of the system can be analyzed using the reduced two–dimensional system

$\frac{ds}{dt}=-\beta si,$

$\frac{di}{dt}=\beta si-\gamma i.$

The Jacobian matrix of the system is obtained by computing the partial derivatives with respect to the state variables s and i:

$J(s,i)=\left(\begin{array}{cc}\hfill -\beta i\hfill & \hfill -\beta s\hfill \\ \hfill \beta i\hfill & \hfill \beta s-\gamma \hfill \end{array}\right).$

To determine the stability of the disease-free equilibrium, we evaluate the Jacobian matrix at

$E_0=(1,0).$

Substituting s = 1 and i = 0 gives

$J(E_0)=\left(\begin{array}{cc}\hfill 0\hfill & \hfill -\beta \hfill \\ \hfill 0\hfill & \hfill \beta -\gamma \hfill \end{array}\right).$

The eigenvalues of this matrix are

${\lambda }_1=0,{\lambda }_2=\beta -\gamma .$

The sign of the second eigenvalue determines the stability of the equilibrium.If

$\beta -\gamma <0,$

then

${\lambda }_2<0,$

and the disease-free equilibrium is locally stable.

Using the definition of the basic reproduction number

$R_0=\frac{\beta }{\gamma },$

This condition can be rewritten as

$R_0<1.$

Therefore, when R₀ < 1, any small number of infected individuals will eventually disappear, and the disease-free equilibrium remains stable.On the other hand, if

$R_0>1,$

then

${\lambda }_2>0,$

and the disease-free equilibrium becomes unstable. In this case, the infection grows and an epidemic outbreak occurs.

This result confirms that the basic reproduction number acts as a threshold parameter governing the qualitative behavior of the epidemic dynamics.

4.1. Quantitative convergence rates

The stability of the disease-free equilibrium can be quantified by the eigenvalues of the Jacobian matrix. At E₀ = (1, 0), the eigenvalues are λ₁ = 0 and λ₂ = β − γ = γ(R₀ − 1).

• When R ₀ < 1 (stable): The infection decays exponentially as $I(t)\propto e^{{\lambda }_{2}t}=e^{-\gamma (1-R_0)t}$ . The characteristic decay time is τ = 1/[γ(1 − R₀)]. For example, with γ = 0.1 and R₀ = 0.8, we have τ = 1/[0.1 × 0.2] = 50 days.

• When R ₀ > 1 (unstable): The infection grows exponentially as $I(t)\propto e^{\gamma (R_0-1)t}$ . The doubling time is t _double = ln(2)/[γ(R₀ − 1)]. For our default parameters (γ = 0.1, R₀ = 5), the doubling time is approximately ln(2)/(0.1 × 4) ≈ 1.73 days.

These quantitative rates provide a more precise characterization of the stability properties beyond simple threshold statements.

5. Numerical investigation of the epidemic dynamics

In order to complement the analytical results obtained in the previous sections, we now investigate the dynamical behavior of the epidemic model through numerical simulations. The objective of this analysis is to illustrate how the theoretical properties of the model appear in the temporal evolution of the epidemic and to explore the influence of key epidemiological parameters on disease propagation. We begin by examining the basic epidemic dynamics of the SIR system. The temporal trajectories of the susceptible, infected, and recovered populations are shown in Figure 1. The figure illustrates the classical epidemic behavior predicted by the model. At the early stage of the outbreak, the number of infected individuals increases rapidly due to the high availability of susceptible hosts. As the infection spreads, the susceptible population gradually decreases, which reduces the effective transmission rate and eventually leads to a decline in the number of infected individuals. At the same time, the recovered population increases steadily as infected individuals recover and acquire immunity.

Figure 1.

Temporal evolution of the susceptible, infected, and recovered populations in the SIR epidemic model.

To further understand the qualitative structure of the epidemic dynamics, we analyze the phase-plane representation of the system. The phase portrait displayed in Figure 2 describes the relationship between the susceptible and infected populations during the course of the epidemic. The trajectory illustrates how the epidemic evolves in the state space as susceptible individuals become infected and eventually move to the recovered class. This geometric representation provides additional insight into the nonlinear interaction between epidemiological compartments.

Figure 2.

Phase-plane representation showing the dynamical relationship between the susceptible and infected populations.

Another important aspect of epidemic modeling concerns the role of stochastic effects. In real populations, disease transmission is often influenced by random fluctuations that arise from environmental variability or random contact patterns. To illustrate this phenomenon, we compare the deterministic epidemic trajectory with stochastic realizations of the model. As shown in Figure 3, stochastic perturbations generate multiple epidemic paths around the deterministic solution. These fluctuations are particularly noticeable during the early phase of the epidemic when the number of infected individuals is relatively small.

Figure 3.

Comparison between deterministic epidemic dynamics and stochastic realizations of the infection process.

Public health interventions play a crucial role in controlling epidemic outbreaks. One of the most effective strategies for limiting disease transmission is vaccination. To evaluate the impact of vaccination on epidemic evolution, we simulate the model under different vaccination rates. The results presented in Figure 4 show that increasing the vaccination rate significantly reduces the peak number of infected individuals and shortens the duration of the epidemic. This behavior highlights the importance of vaccination programs in preventing large-scale outbreaks.

Figure 4.

Influence of vaccination rate on the evolution of the infected population.

5.1. Quantitative impact of vaccination on infection peak

The effect of vaccination on the infection peak is quantified as follows. Using the baseline parameters (β = 0.5, γ = 0.1, N = 1000, I₀ = 1), the infection peak without vaccination (ν = 0) reaches approximately 0.65 (i.e., 650 infected individuals out of 1000). When vaccination is introduced at rate ν = 0.1, the peak reduces by approximately 25%. At ν = 0.2, the peak reduces by approximately 48%. At ν = 0.3, the peak reduces by approximately 65%. These results demonstrate that even moderate vaccination coverage can substantially reduce the maximum burden on healthcare systems during an epidemic. Figure 4 illustrates these reductions graphically.

The recovery rate also plays a fundamental role in shaping the progression of the epidemic. Higher recovery rates reduce the average infectious period, thereby limiting the time during which infected individuals can transmit the disease. Numerical simulations illustrating this effect are presented in Figure 5. The results clearly indicate that faster recovery significantly reduces the infection peak and accelerates the decline of the epidemic.

Figure 5.

Effect of the recovery rate on the dynamics of the infected population.

Finally, the epidemic curves shown in Figure 6 summarize the overall temporal behavior of the infection. These curves display the characteristic rise and decline of the infected population during the outbreak. The initial growth phase corresponds to the regime where the reproduction number exceeds unity, while the decline phase reflects the depletion of susceptible individuals and the increasing immunity in the population.

Figure 6.

Typical epidemic curve illustrating the rise and decline of the infected population over time.

5.2. Numerical parameter values

All simulations presented in this study use the following parameter values:

When studying vaccination effects, we vary the vaccination rate ν from 0 to 0.3. When studying recovery rate effects, we vary γ while keeping other parameters constant.

5.3. Pseudocode and step-by-step formulation

5.3.1. Deterministic SIR (ODE method)

5.3.2. Stochastic SIR (gillespie algorithm)

5.3.3. Parameter selection

The default parameters (β = 0.5, γ = 0.1, N = 1000, I₀ = 1) were chosen to produce a clear epidemic outbreak with R₀ = 5 for illustrative purposes. These values are representative of a moderately contagious disease. For disease-specific applications, parameters can be adjusted as shown in Table 3. We do not use spectral methods or spectral coefficient update rules in this work.

Table 3.

Estimated parameter values and required vaccination coverage for herd immunity (1 − 1/R₀) for different diseases.

Disease	β	γ	R0 = β/γ	Required vaccination coverage
Seasonal Influenza	0.4	0.33	1.2 - 1.8	17% - 44%
COVID-19 (original)	0.5	0.2	2.5 - 3.0	60% - 67%
Measles	1.2	0.1	12 - 18	92% - 94%

5.4. Application to real-world diseases

Our findings can be applied to specific infectious diseases by calibrating the model parameters accordingly. The table below provides estimated parameter values for three common diseases:

For seasonal influenza (R₀ ≈ 1.5), our model predicts a relatively low infection peak that can be controlled with moderate vaccination coverage (approximately 33%). For COVID-19 (R₀ ≈ 2.5), higher vaccination coverage (approximately 60%) is needed. For measles (R₀ ≈ 15), our simulations indicate that very high vaccination rates (above 93%) are necessary to prevent large outbreaks, consistent with public health guidelines. These examples demonstrate how our quantitative framework can inform disease-specific public health strategies.

6. Discussion of results

The analytical and numerical results obtained in this study provide a comprehensive understanding of the dynamics of epidemic propagation within the classical SIR framework. As discussed in virus dynamics literature, ¹⁵ the interaction between host and pathogen follows similar principles. The applicability of our findings extends to several real-world diseases. For example, using parameter values representative of influenza (β = 0.4, γ = 0.33, yielding R₀ ≈ 1.2), our model predicts a relatively low infection peak that can be controlled with moderate vaccination coverage. In contrast, for measles (β = 1.2, γ = 0.1, yielding R₀ = 12), our simulations indicate that very high vaccination rates (above 90%) are necessary to achieve herd immunity and prevent large outbreaks. These examples demonstrate how the quantitative framework presented in this study can inform disease-specific public health strategies. The mathematical formulation of the model reveals that the interaction between susceptible and infected populations drives the nonlinear behavior of the epidemic process. The equilibrium analysis shows that the system admits a disease-free equilibrium, which represents the situation where the infection disappears from the population. The stability analysis further demonstrates that the stability of this equilibrium is governed by the basic reproduction number R₀. When R₀ < 1, the infection cannot sustain itself in the population and eventually vanishes. In contrast, when R₀ > 1, the infection spreads and an epidemic outbreak becomes inevitable. The numerical simulations provide additional insight into the temporal behavior of the epidemic. The epidemic curves illustrate the classical pattern observed in infectious disease outbreaks, where the number of infected individuals increases rapidly at the beginning of the epidemic, reaches a peak, and then gradually declines as the susceptible population decreases and immunity accumulates.

The phase-plane representation highlights the nonlinear interaction between susceptible and infected individuals, providing a geometric interpretation of the epidemic trajectory. This visualization helps to clarify how the epidemic evolves in the state space and how the system eventually approaches a stable configuration. The comparison between deterministic and stochastic epidemic dynamics demonstrates that random fluctuations can significantly influence the early stages of disease spread. While the deterministic model captures the average behavior of the epidemic, stochastic realizations reveal the variability that can arise in real populations due to random contact patterns.

Furthermore, the numerical experiments examining vaccination and recovery rates emphasize the critical role of intervention strategies in epidemic control. Increasing vaccination coverage substantially reduces the epidemic peak and limits the overall number of infections. Similarly, higher recovery rates shorten the infectious period and accelerate the decline of the epidemic. These findings confirm the importance of public health measures aimed at reducing transmission and increasing recovery efficiency.

6.1. Comparison with other epidemic models

While the classical SIR model is used in this study for its simplicity and analytical tractability, other epidemic models offer different insights depending on disease characteristics. The SEIR model (Susceptible-Exposed-Infected-Recovered) includes an exposed compartment to account for the incubation period, which is particularly important for diseases like COVID-19 or measles, where infected individuals are not immediately infectious. The SIS model (Susceptible-Infected-Susceptible) is more appropriate for diseases such as bacterial meningitis or sexually transmitted infections, where recovery does not confer immunity, and individuals can be reinfected. Additionally, models with vital dynamics (births and deaths) are necessary for studying endemic diseases over long time scales where population turnover affects disease persistence. Despite the availability of these more complex models, the classical SIR model remains valuable because it captures the essential dynamics for acute diseases with permanent immunity (e.g., influenza, measles) and provides a simpler, more tractable framework for understanding threshold behavior, stability analysis, and intervention effects. The choice of model should be guided by the specific disease and research question.

6.2. Limitations of the classical SIR model

The classical SIR model, while valuable, has several limitations that should be acknowledged. It lacks an exposed/latent compartment, assuming immediate infectiousness after infection, which is unrealistic for diseases like COVID-19 or HIV that have significant incubation periods. It also assumes permanent immunity, ignoring waning immunity and reinfection as seen in seasonal influenza or COVID-19. The model assumes a closed population with no births, deaths, or migration, which limits its applicability to long-term endemic dynamics. Additionally, it assumes homogeneous mixing, disregarding population structure such as age, spatial distribution, and social networks. Finally, the transmission rate β and recovery rate γ are assumed constant, whereas in reality they may vary due to seasonality, behavioral changes, or interventions. Despite these limitations, the classical SIR model remains a powerful tool for understanding fundamental epidemic mechanisms, and future work can extend it to address these issues.

6.3. Effect of waning immunity and reinfection

The classical SIR model assumes permanent immunity after recovery. However, for many diseases such as seasonal influenza and COVID-19, immunity wanes over time, and reinfection is possible. If waning immunity were included, the model would become an SIRS model (Susceptible-Infected-Recovered-Susceptible), where recovered individuals return to the susceptible class at a rate ω (waning rate). In this case, the dynamics would change qualitatively:

1. The infection may not die out completely. Instead, the system could reach an endemic equilibrium where the disease persists in the population at a steady level.

2. Periodic outbreaks or sustained oscillations could occur, depending on parameter values.

3. The basic reproduction number R₀ would still determine the initial outbreak, but the long-term dynamics would also depend on the waning rate ω.

4. Vaccination strategies would need to be maintained over time, as immunity wanes in both vaccinated and naturally recovered individuals. A one-time vaccination campaign would not achieve permanent herd immunity.

Our conclusion that higher vaccination rates reduce infection peaks would still hold under waning immunity. However, the long-term dynamics would show repeated outbreaks if vaccination coverage is not sustained. This highlights the importance of booster doses and ongoing vaccination programs for diseases with waning immunity. Future work should extend our model to include waning immunity and reinfection to study these more complex dynamics.

6.4. Validation and limitations for public health use

We acknowledge that the model presented in this study has not been validated against real epidemic data. The simulations are based on synthetic/hypothetical data, and the parameter values (β = 0.5, γ = 0.1, etc.) were chosen for illustrative purposes rather than calibrated to any specific disease outbreak. Therefore, the findings of this study should be interpreted as qualitative insights into epidemic dynamics and intervention effects, rather than quantitative predictions for real-world public health decision-making. Our framework is intended to serve educational and methodological purposes, helping researchers and students understand:

• How the basic reproduction number R₀ determines outbreak thresholds,

• How vaccination and recovery rates quantitatively affect infection peaks,

• How stochastic effects influence epidemic trajectories,

• How phase-plane analysis reveals nonlinear interactions between compartments.

For real-world public health decision-making, models must be validated against local outbreak data, calibrated with disease-specific parameters, and interpreted by domain experts. We have added validation with real data (e.g., COVID-19 case counts, influenza surveillance data) as a priority for future research. We have also clarified this limitation in the Conclusion section.

6.5. Validation, parameter estimation, and practical applicability

We acknowledge that the simulations presented in this study are based on synthetic or hypothetical data rather than real outbreak data. The primary purpose of this work is to provide a pedagogical and methodological framework for understanding the fundamental dynamics of SIR epidemics and the quantitative effects of intervention strategies. To apply the model to real-world diseases such as COVID-19, influenza, or measles, parameter estimation and calibration are essential. The transmission rate β and recovery rate γ can be estimated from outbreak data using methods such as least squares fitting, maximum likelihood estimation, or Bayesian inference. For example, using parameter values representative of influenza (β = 0.4, γ = 0.33, yielding R₀ ≈ 1.2) or measles (β = 1.2, γ = 0.1, yielding R₀ = 12), our model can generate disease-specific predictions about infection peaks and required vaccination coverage. For COVID-19, early studies estimated R₀ between 2.5 and 3.0, with γ ≈ 0.2 day⁻¹ (recovery in 5 days), yielding β ≈ 0.5 − 0.6 day⁻¹. We agree that validation against real outbreak data (e.g., COVID-19 case counts, influenza surveillance data) is an essential next step to enhance the practical usefulness of the framework, and we have added this as a priority for future research.

7. Conclusion

In this work, we investigated the dynamical behavior of the classical SIR epidemic model using both analytical and numerical approaches. The mathematical formulation of the model allowed us to derive the equilibrium points and analyze the stability properties of the disease-free state. Stability analysis showed that the basic reproduction number R₀ = β/γ acts as a threshold parameter that determines whether an epidemic outbreak occurs. When R₀ < 1, the disease eventually disappears from the population, whereas when R₀ > 1, the infection spreads and generates an epidemic wave.

Summary of key findings: The disease-free equilibrium E₀ = (1, 0, 0) is locally asymptotically stable when R₀ < 1 and unstable when R₀ > 1, where the basic reproduction number R₀ is derived as the ratio β/γ. Numerical simulations show the characteristic rise and decline of infected individuals over time, where the peak infection level and outbreak duration are directly influenced by model parameters. Increasing the vaccination rate from 0 to 0.3 reduces the infection peak by approximately 65%. Higher recovery rates shorten the infectious period and accelerate epidemic decline, leading to a lower infection peak. Stochastic realizations reveal random fluctuations around the deterministic solution, particularly when the number of infected individuals is small. Finally, phase-plane analysis provides a geometric representation of the nonlinear interaction between susceptible and infected populations, offering insight into how the epidemic evolves in state space.

Implications of simulations for public health: The simulation results have several important implications for public health planning and intervention strategies. First, the quantitative relationship between vaccination coverage and peak infection reduction demonstrates that even moderate vaccination rates can substantially mitigate epidemic outbreaks. Second, the effectiveness of recovery rate improvements highlights the importance of treatment and medical care in controlling disease spread. Third, the stochastic variability observed in early epidemic stages suggests that initial conditions and random fluctuations can influence outbreak trajectories, which should be considered when making short-term predictions.

Novelty and advantages of this study: Although the classical SIR model has been extensively studied, the novelty of our work is threefold. First, we provide a combined analytical and numerical framework that explicitly links the stability condition R₀ < 1 to quantitative reductions in infection peaks. Second, we present phase-plane visualizations that clarify the nonlinear interaction between susceptible and infected populations in a way often omitted in standard textbook presentations. Third, we offer simulation-based quantification of intervention effects (e.g., vaccination reduces peak infection by approximately 65% when coverage reaches 30%). These advantages make our work useful not only theoretically but also for practical public health decision-making.

Future research directions: Future research may extend the present framework in several directions. First, demographic effects such as births and deaths can be incorporated to study endemic dynamics. Second, spatial diffusion can be added to model the geographic spread of infectious diseases. Third, stochastic transmission mechanisms can be further developed to more realistically capture random fluctuations. Fourth, additional epidemiological features such as exposed compartments (SEIR models), waning immunity, or age-structured populations could be considered. Finally, validating the model against real-world epidemic data (e.g., COVID-19, influenza) would enhance the practical applicability of the framework. We also acknowledge that our model has not been validated against real epidemic data, and the results should be interpreted as qualitative insights rather than quantitative predictions for public health decision-making. Future work should focus on calibration and validation using real outbreak data, such as COVID-19 case counts or influenza surveillance records. The results presented in this study confirm that mathematical epidemic models provide powerful tools for understanding disease transmission mechanisms and for evaluating the impact of intervention strategies.