SIRV fractional epidemic model of influenza with vaccine game theory and stability analysis

1.   Introduction

Influenza is a global infectious disease that affects humans. It has taken away millions of lives in the 20th century ^[1]. These viruses become seasonal influenza after widespread transmission. Influenza typically spreads through respiratory droplets, such as coughing and sneezing. For most individuals, influenza is mild, and people usually recover fully within 3 to 5 days. However, in some cases, the severity of the disease can lead to hospitalization and death. Symptoms include fever, fatigue, cough, and in severe cases, viral pneumonia. Furthermore, the emergence of novel influenza strains, including pandemic strains like H1N1 in 2009, underscores the necessity for ongoing monitoring, accurate prediction, and effective control strategies. Generally, the compartment model is used to simulate disease progression and control. This model is widely used in the epidemiological study of infectious diseases, a concept first proposed by Kermack and McKendrick ^[2]. Many researchers have previously attempted to model influenza using the compartment model ^[3,4]. The complexity of influenza transmission dynamics, which involves intricate interactions, necessitates the development of sophisticated mathematical models to gain insights into its spread and to inform public health interventions.

While medications are used to treat influenza, the most effective way to prevent the disease is vaccination ^[1]. Mass vaccination campaigns ^[5] have successfully controlled the spread of smallpox and polio, ultimately leading to their eradication. Some individuals who require vaccination can develop severe complications from influenza, such as viral pneumonia, and may even die from flu-related complications ^[6]. However, the public's understanding of influenza vaccines is limited, and many people fail to appreciate the seriousness of the disease fully. Influenza vaccination is primarily conducted voluntarily. Therefore, it is essential to develop a more effective vaccination program to establish a comprehensive strategy for vaccine distribution. MacDonald was the first to propose the infectious disease compartment model that incorporates vaccination ^[7]. Since then, numerous simple compartmental mathematical models have been widely utilized to assess the impact or potential impact of vaccination on the transmission of certain human diseases. Shim ^[8] studied the SIR and SIS models of vaccination. Gumel et al. ^[9] used the SEVIR model to study the spread of SARS. It is important to recognize that vaccination is not a one-time solution. Yan et al. ^[10] studied a new class of SIRV models for influenza vaccination. To effectively control outbreaks, we must adapt our vaccination strategies as needed, which includes promoting booster shots and developing vaccines to combat newly mutated strains.

Game theory constitutes a mathematical framework that examines strategic decision-making among rational individuals or groups in circumstances entailing conflict or cooperation. Decision-making behaviors of individuals can be effectively simulated using game theory. Its application in epidemiology has gained significant traction, especially in comprehending how individual behaviors transform during an epidemic, encompassing vaccination decisions, social distancing, and adherence to public health measures. We can consider the strategic interactions among individuals as they make decisions that impact their own and others' health. For example, individuals might assess the costs and benefits of vaccination, considering factors such as the perceived risk of infection, the vaccine's efficacy, and potential side effects. By modeling human decision-making processes within game-theoretic frameworks, researchers can investigate how individual actions affect the overall dissemination of influenza and how these actions can be influenced by incentives. Game theory has been utilized to model the behavioral choices of individuals during the transmission of infectious diseases. Researchers' utilization of game theory to predict the spread of epidemic infectious diseases has emerged as a significant branch of modern bio-mathematics. We need to contemplate how individuals respond to vaccination. Hence, it is indispensable to employ game theory to investigate vaccination strategies. We assume that all players are rational and that the actual strategies of each individual interact within the game-theoretic context. Christ Bauch established a dynamic SIR model using game theory ^[11]. Bauch and Bhattacharyya proposed a model for vaccine delay strategies that account for age structure ^[12]. Liu et al. ^[13] used game theory to study chickenpox vaccination with age structure. In some models, it is assumed that individuals will imitate the behavior of their neighbors in a certain grid ^[14,15,16]. Reluga developed game theory models that use stochastic differential equations for modeling ^[17,18]. Game-theoretic models that consider individual differences are typically represented using individual networks, which consider individual heterogeneity ^[19,20]. These studies have demonstrated that gaming strategies are efficacious in preventing and controlling the spread of infectious diseases. They can be optimized to minimize costs, budgets, and the number of infections.

Traditional epidemiological models, often based on integer-order differential equations, have been instrumental in understanding the fundamental mechanisms of disease transmission. However, these models may not fully capture the long-memory effects, non-local interactions, and heterogeneous responses frequently observed in real-world epidemic systems. With the development of fractional calculus theory, an increasing number of researchers are beginning to apply it in biology to more accurately describe the behavior of model systems. In contrast to integer-order models, fractional derivative epidemic models have better systematic memory characteristics. The Caputo fractional derivative, a specific type of fractional derivative, has attracted considerable attention in recent years due to its capacity to model systems that show memory effects and non-Markovian behavior—traits intrinsic to many biological and epidemiological processes. Additionally, fractional models possess one more degree of freedom than integer models when fitting data. Furthermore, fractional differential equations are more suitable than integer differential equations for describing non-local phenomena and can provide more realistic explanations. The fractional-order infectious disease model can not only describe the transmission process of infectious diseases in the population more accurately but also predict the development trend of infectious diseases better. Fractional derivatives are introduced into the model in a specific way, such as replacing the integer derivative with the Caputo fractional derivative. In recent years, studying of models involving fractional calculus has proliferated. Many researchers have studied diseases such as dengue ^[21], measles ^[22], AIDS ^[23], and influenza ^[24] using fractional-order derivative. Li et al. ^[25] used new parameter estimation techniques to study multi-term fractional dynamical epidemic models. Defterli ^[26] studied temperature-dependent fractional order dengue models. Kozioł et al. ^[27] studied a fractional order SIR model for COVID-19. Miskovic-Stankovic et al. ^[28] studied a fractional order derivative model for different drugs. These studies demonstrate the utility of the fractional-order system by illustrating its ability to provide more accurate and nuanced insights into complex dynamical processes, such as those observed in infectious disease modeling.

In this paper, we investigate the effects of vaccination information and infection information on individuals' behaviors regarding vaccination, based on the transmission dynamics between individuals. Drawing on the objective laws of infectious diseases, we consider a fractional-order model that incorporates the vaccination game and conduct a mathematical analysis of this model. Our innovation primarily involves the introduction of a game-theoretic framework for fractional-order systems, as well as the incorporation of a reality-based individual behavioral response function into the model. We develop a simulation based on the aforementioned factors. Furthermore, the model we establish can serve as a valuable tool for policymakers and public health officials. By simulating various scenarios and assessing the impact of different interventions, policymakers can identify the most effective strategies for mitigating the spread of influenza. This can lead to more informed and targeted public health policies that are better equipped to address the challenges posed by infectious diseases. Consider that integrating a fractional-order model with game theory can yield a model that more accurately represents real-world situations. This is an aspect that is lacking in this study. Not only that, but the novelty of this paper lies in our development of a nonlinear response function for a more realistic simulation of vaccination. Additionally, we analyze the Ulam-Hyers stability in a fractional-order system.

The structure of this article is as follows: In the next section, we introduce the establishment of the model. In Section 3, we extend the model to fractional order and analyze it. In Section 4, we analyze the Lyapunov stability and the Ulam-Hyers global stability of the model. In Section 5, we analyze the model's behavior, investigate the impact of different parameters on the spread of influenza, and discuss the implications of our findings for public health policy and practice. We conclude the paper in Section 6.

2.   SIRV Mathematical model based on ODE

In this section, we present the novel model. Considering the transmission characteristics of influenza and the vaccination status, we first consider the classical SIR model and add the V compartment based on the SIR model. We assume that the spread of influenza is uniform and that the probability of transmission is proportional to the frequency of contact between susceptible and infected individuals. We also assume that deaths due to the disease are negligible and that the total population remains constant. People are divided into four compartments: $S$: susceptible, $I$: infected, $R$: recovered, and $V$: immunized. Susceptible individuals will enter the infection compartment after effective contact with infected individuals, and infected individuals recover from a certain proportion and then enter the recovery compartment. Because the continuous mutation of the influenza virus requires annual updates to the vaccine to prevent the disease effectively, we only consider only infection for a short period (less than a year). We assume that the recovered population is permanently immune. Our influenza epidemic model is shown in Figure 1. From the figure, we can deduce the following system of differential equations

$S$ represents susceptible individuals, $I$ represents infected individuals, $R$ represents recovered individuals, and $V$ represents immunized individuals. $\sigma$ represents the effectiveness of vaccination, with $\sigma \in [0, 1]$. $\phi$ represents the rate at which vaccine effectiveness diminishes. $\beta$ is the effective transmission rate, $\gamma$ denotes the recovery rate, and $1/\mu$ represents the average life expectancy of this model population. We ignore the mortality rate from influenza. Furthermore, $x\in[0, 1]$ stands for the vaccination rate, which is simulated by means of game theory. In this model, we assume that the total population is $N$. Let $N = S+V+I+R$, and $s = S/N$, $i = I/N$, $r = R/N$, and $v = V/N$, such that $s+v+i+r = 1$. The model can be simplified to

It was noted that $x$ followed vaccination behavior. The vaccination situation is closely related to many factors. It is assumed here that voluntary vaccination is available to individuals across the population, and that individuals are influenced by and respond nonlinearly to information about the costs of strategies in their current state. The probability that the considered strategy will be adopted can generally be assumed to be proportional to the difference between the expected returns of the two strategies. Therefore, we first quantify the human behavioral tendencies by dividing the costs into the costs of infection and the costs of vaccination. We first assume that the cost of inoculation is $c_v$ and the cost of infection as $c_i$. From the meaning of $\beta$, we set the rate of individual infection as $\beta i$. We assume that vaccination is effective. Thus, the cost of infection for an unvaccinated individual as $\beta ic_i$. We define information function $g_1(t)$ as follows

If $g_1(x) > 0$, people are more likely to get vaccinated, and vice versa if $g_1(x) < 0$.

The impact of the severity of the outbreak on vaccination also requires consideration here. Hence, we also consider the new infection to be another information function. When people realize that the number of new cases is relatively large, they consider the current outbreak more serious. This means that the need for vaccination has become even greater. Hence, the next thing we consider is the number of new cases per day. We set it to $g_2(t)$.

These two information functions are not equivalent, so let the weights be $m$, and we set the net information function by

Here, $0 < m < 1$. If $m = 0$, vaccination behavior is determined by $g_2(t)$. When $m = 1$, it is completely determined by $g_1(t)$. This models the effects of different weights.

In fact, depending on the situation, the response function is not always linear. Moreover, we consider the impact of information transmission on individuals' vaccination behavior. This is so the individuals will imitate their neighbors. In real life, individuals also have limited social circles. We assume that individuals imitate the behavior of others, but not in a completely arbitrary manner. Instead, we assume that individuals can only learn from their neighbors within a specified network. We assume that the response function for the behavior of inoculation conforms to the Fermi rule. Using the Fermi function^[29], we define the behavioral response as

here $b > 0$ represents the strength of selection. When $b$ is larger, the individual is more sensitive to information, suggesting a greater level of rationality.

Finally, we derive the equation for behavioral changes with $x$ by classical imitation dynamics method ^[11]. Similar to this, we assume that individuals randomly select other members of the population at a fixed, constant rate. If the strategy of the sampled individual yields a higher payoff, it is adopted with a probability proportional to the anticipated increase in payoff. In this model, $f_R(t)$ represents the payoff gain associated with switching to the vaccinator strategy. When $f_R(t) > 0$, nonvaccinators probably become vaccinators, but not vice versa. $x$ denote the relative frequency of vaccinators. If individuals sample at rate of $\sigma$, then a non-vaccinator samples vaccinators at rate $\sigma x$. Consequently, a non-vaccinator adopts the vaccinator strategy with probability $af_R(t)$, where $a$ is the constant. Hence, the dynamic of $x$ when $f_R(t) > 0$ is governed by:

$(1-x)$ represents the frequency of individuals who are not vaccinated. Conversely, if $f_R(t) < 0$, it follows from a similar step that

Taking $k = \sigma a$ is a coefficient for imitation intensity, we obtain

Afterwards, it is concluded that the system we want to study as follows:

In summary, we developed a new ODE model by incorporating an additional compartment into the traditional SIR model and integrating it with vaccination game theory. This extended model captures the dynamics of vaccination strategies within the population and examines how these strategies influence disease transmission. To simplify the analysis, we streamlined the model by eliminating one of the equations, thereby maintaining its essential features while enhancing clarity.

3.   Fractional-order model and its analysis

In this section, given the advantages of fractional derivatives, we introduce the definition of fractional calculus and generalize our proposed model to fractional-order systems.

3.1.   Fractional-order derivatives: concepts and definitions

Fractional derivatives are considered a more effective approach to modeling disease epidemics because they accurately capture the memory effect that often influences people's responses to diseases. The non-negativity and boundedness of the solutions of this fractional-order system are then proved. Finally, the existence and uniqueness of the solution under certain conditions are proved. We will begin by providing the definition of the Caputo derivative.

Definition 1. ^[30] For an integrable function f, set $\alpha\in \mathbb{R}^+$, $m-1 < \alpha < m$, the $\alpha$ order Caputo derivative of function $f$ is given by

The Riemann-Liouville fractional integral of $\alpha$ is given by

Definition 2. ^[30] The Laplace transform of order $\alpha$ of the Caputo differential operator is given by

or

Definition 3. ^[30] The Mittag-Leffler function of order $n > 0$ is defined by

The two-parameter Mittag-Leffler function with $n_1, n_2 > 0$ is defined by

And we give the Laplace transform of two-parameter Mittag-Leffler function

Finally, we give the generalized mean value theorem.

Lemma 1. Generalized mean value theorem ^[31]. For $0 < \alpha \leq 1$, let $g(t) \in C([a, b])$ and $^C \! D_t^{\alpha} g(t) \in(a, b]$. Then it holds

We generalize the model to the following system in terms of the Caputo derivative.

3.2.   Nonnegative solution of model

In this section, we demonstrate that the fractional system (3.9) is biologically well-posed. Next, we demonstrate that the solution of system (3.9) is nonnegative and bounded at all times, or that $\Lambda$ is an invariant set,

Lemma 2. For $S(0)\geq0$, $I(0)\geq0$, $R(0)\geq0$ and $V(0)\geq0$, the solution of the system (3.9) is nonnegative.

Proof. First, we consider the trajectory of the solution of (3.9) along the $S$-axis, i.e., $I(0) = R(0) = V(0) = 0$, $S(0)\geq 0$. From the first equation in (3.9), we have

Taking the Laplace transform of both sides of the Eq (3.10),

Then,

By inverse Laplace transform, we have

Using similar argumentation, we respectively get trajectories of the solution of (3.9) along the $I$-axis, $R$-axis, and $V$-axis,

Finally, we state that the solution of the model (3.9) is positive in the $\Lambda$ as follows: We set $\tilde{t} > 0$ such that $S(\tilde{t}) = 0$, $I(\tilde{t}) > 0$, $R(\tilde{t}) > 0$, $V(\tilde{t}) > 0$, and $S(t) < S(\tilde{t})$. On this plane, $D_t^\alpha S(t)|_{t = \tilde{t}} = \mu N > 0.$ By the generalized mean value theorem, we obtain

We conclude that $S(t) > S(\tilde{t})$. This contradicts our previous assumption $\tilde{t}$. This analysis is similar for $I(t), R(t)$ and $V(t).$ □

Lemma 3. The solutions $S(t), I(t), R(t)$ and $V(t)$ of the model (3.9) are bounded when they are positive.

Proof. According to the conclusion of Lemma 2 we can get that when $t > 0$, $^C \! D_t^\alpha S(t)$, $^C \! D_t^\alpha I(t)$, $^C \! D_t^\alpha R(t)$, $^C \! D_t^\alpha V(t) > 0$. If we set $N(t) = S(t)+I(t)+R(t)+V(t)$, with the simple calculation, obtain that

Clearly, when $t > 0$, $N(t) = N$. It shows that $S(t), I(t), R(t)$ and $V(t)$ are bounded and the set $\Lambda$ is a positively invariant region of $\mathbb{R}^4$. □

3.3.   Existence and uniqueness of solution

Next, we introduce the existence and uniqueness of a solution of the system (3.9) by imposing certain conditions.

We notice that the system (3.9) can be rewritten as follows

Let $X(t) = (S, I, R, V, x)^T$, and we define $f(t, X(t)) = (f_i)^T, i = 1, 2, \dots, 5$ and $T\in \mathbb{R}^+$. Then, the system (3.9) is equivalent to

We turn it to integral representation:

For our analysis, we consider $E = C([0, T); \mathbb{R}^5)$, the Banach space of continuous functions with the norm

Moreover, we define a operator $P:E\to E$ by

we observe that $P$ is well-defined due to the continuity of $f$.

Theorem 1. The function $f$ defined in (3.18) is a Lipschitz function on $E$.

Proof. Set $\bar{X} = (\bar{S}, \bar{I}, \bar{R}, \bar{V}, \bar{x})$ represents another point in $E$. According to the Lemma 2 and Lemma 3, we find the $b_1$, $b_2$, $b_3$, and $b_4$ satisfy

First, discuss the $f_1$ of $f$, we can observe that

Here, we set $L_1 = \max\{\frac{\beta b_2}{N}+\sigma+\mu, \frac{\beta b_1}{N}, \sigma b_1 \}$. Similarly, we obtain $L_2, L_3, L_4, L_5$. $L_2 = \max\{\frac{\beta b_1}{N}+\gamma+\mu, \frac{\beta b_2}{N}\},$ $L_3 = \max\{\mu, \gamma\},$ $L_4 = \phi+\mu,$ $L_5 = k.$ Finally, claim that

where $L = L_1+L_2+L_3+L_4+L_5$. □

Theorem 2. Set $M = \frac{N}{\Gamma(\alpha+1)}$.If $ML < 1$, then there exists a unique solution of (3.19) on $[0, T)$.

Proof. We utilize Banach contraction mapping principle on $P$. By definition,

We define $\kappa > \frac{\|X_0\|+MT^\alpha}{1-ML}$ and a closed convex set $B_\kappa = \{X\in E: \|X\|_E\leq\kappa\}$. We show that $PB_\kappa\subseteq B_\kappa$. Note that $f$ satisfies the Lipschitz condition by inequality (3.26). Let a point $X\in B_\kappa$, we obtain

Therefore, $PB_\kappa\subseteq B_\kappa$ and $P$ is a self-map. Now, we prove that $P$ is a contraction. We set $X, \tilde{X}\in E$ satisfying (3.19), and observe that

According the condition $ML < 1$, $P$ is a contraction mapping. By the Banach contraction mapping principle, $P$ has a unique fixed point on $[0, T)$. □

We prove that the solution of the model (3.19) can be extended to $[0, +\infty)$.

Theorem 3. Set $M = \frac{N}{\Gamma(\alpha+1)}$. If $ML < 1$, then there exists a unique solution of Eq (3.19) on $[0, +\infty)$.

Proof. Note that $t\in[0, T)$, $X(t)$ is bounded. We show that there exists $\lim_{t\to T^-}X(t) = X^*$. It means $\exists \eta\in[0, T)$ such that for $t\in(\eta, T)$, $\|X(t)-X^*(t)\| < \varepsilon$. Set $t_k$ and $\tau_k$, and for $\forall\varepsilon > 0$,

Then, for a large $k\geq k_0$, we conclude that

This makes a contradiction and, thus, the limit $\lim_{t\to T^-}X(t)$ exists. To make the solution $X$ continuous on $[0, T)$, we define $X(T) = \lim_{t\to T^-}X(t)$. In this way, we prove that this solution is continuous on $[0, T)$.

Set $h < +\infty$. According to the integral representation (3.21) and Theorem 2, we obtain

Set $\kappa_h > \frac{\|X_1\|+MT^\alpha}{1-ML}$ and a closed convex set $B_h = \{X\in E: \|X\|_E\leq\kappa\}$. Similar to Theorem 2, we conclude that according to the Banach contraction mapping principle, $P$ has a unique fixed point $X^*$ on $[T, T+h]$. Hence, we get the new solution

Hence, $[0, T]$ is not the maximal interval of existence for the solution. Finally, we prove the uniqueness of the model similar to Theorem 2. Note that $f$ satisfies the Lipschitz condition by inequality (3.26). Set 2 solutions $X$ and $\tilde{X}$ satisfy (3.19) observe that

According the condition $ML < 1$, $P$ is a contraction mapping. By the Banach contraction mapping principle, $P$ has a unique fixed point on $[0, +\infty)$.

□

4.   Stability analysis

In this section, we provide an analysis of the equilibria in this model. We first set the left hand side of (3.19) to zero to solve each equilibrium point. It is equivalent to finding equilibria for an integer-order model. Then, we analyze the local stability of each equilibrium point. The uniformly Lyapunov stability of solution in the fractional model follows from ^[32].

4.1.   Boundary equilibria and their existence

We find the equilibrium $E_n = (s_n, i_n, v_n, x_n)$ of the system (2.2). First, we set the left hand side of (2.9) to zero to solve each equilibrium point. Based on the reality, we know that $x\in [0, 1]$, and firstly consider the condition $x = 0$ and $x = 1$.

Set $x = 0$, and we obtain the following equations

When $i = 0$, there is a disease-free equilibrium

If $i \neq 0$, the endemic equilibrium is

According to the given method of ^[33], we obtain the basic reproduction number $\mathscr{R}_0$ of the model. Here, we set $\mathscr{F}$ and $\mathscr{V}$ as

We obtain that $F = (\beta S_0), V = (\mu+\gamma)$. Finally, we set $S_0 = 1$, and obtain

Due to the significance of the model, we consider that the equilibria should, in the nonnegative cone, be

We verify the conditions under which equilibria exist within the region. $E_1$ belongs to $\Lambda$. When $i_2 > 0$, we conclude that $E_2 \in \Lambda$. By computation, $i_2 > 0$ is equal to $\mathscr{R}_0 > 1$.

Next, we consider that the condition of $x = 1$. We obtain

When $i = 0$, there is a disease-free equilibrium

There is no doubt that $E_3 \in \Lambda$ by simple calculation.

From $i \neq 0$, we obtain the equilibrium $E_4$

We indirectly represent $i_4$ and state the existence of $E_4$. We set

called effective reproduction number. Then we obtain

By calculating this, we can obtain the existence condition $E_4\in\Lambda$ is

4.2.   Internal equilibrium and its existence

Finally, we set the right sides of all equations to 0, and obtain

Solve the equation, we obtain $E_5 = (s_5, i_5, v_5, x_5)$. Moreover $s_5 = \frac{1}{ \mathscr{R}_0 }.$ $x'(t) = 0$, such that $g(t) = 0$, by calculating, we obtain

Using $s_5, i_5$, we shall express $v_5, x_5$.

Next, we will determine the conditions for existence of $E_5$, i.e.,

When $\mathscr{R}_0 > 1$, we find that $0 < s_5 < 1$.

By $0\leq m \leq 1$, such that $i_5 > 0$. Hence, by (4.6), we obtain $v_5 < 1$. For simplicity, we define $R$ as

Assume that $\mathscr{R}_v < R < \mathscr{R}_0$, (4.6) and (4.7) mean $0 < x_5 < 1$. From $R > 1$ we can verify that $E_5 \in \Gamma_5$.

4.3.   Lyapunov stability

In this section, we prove the Lyapunov stability of the equilibria points in system (3.19). For the Caputo fractional-order linear model, the French mathematician Denis Matignon ^[34] obtained the following conclusion:

Theorem 4. For a fractional-order system (3.19), its eigenvalue $\lambda_i$ satisfies

if and only if the fractional-order system (3.19) is locally Lyapunov asymptotically stable.

As shown in Figure 2, when $\alpha = 1$, it is the normal Lyapunov stability problem for integer-order systems, and the Lyapunov stability condition for integer-order models is necessary for the stability of fractional-order systems with order $\alpha$. Hence, it is equivalent to analyzing the local stability of an integer-order system. The integer-order Lyapunov stability of the model is discussed next.

First, we linearize the system using Taylor series expansions. For example, we linearize $s'(t)$.

We write the right hand side of $s'(t)$ as $f(s, i, v, x)$. After expanding the function around $E_n$, we have:

$o(\rho)$ represents infinitesimals of higher order. We know that $f(E_n) = 0$. We set $y_1 = s-s_n$, $y_2 = i-i_n$, $y_3 = v-v_n$, $y_4 = x-x_n$, and $\frac{ds}{dt} = \frac{dy_1}{dt}$, $\frac{di}{dt} = \frac{dy_2}{dt} \dots$, we obtain:

Similarly,

where

where $\overline{g(t)} = m(\beta i_nc_i-c_v)+(1-m)\beta N s_n i_n$.

Next, we analyze the Lyapunov stability of $E_n$.

According to the theory of linear ODE systems, we set the solution of the system is $(s, i, v, x) = e^{\lambda t}(\bar{c_1}, \bar{c_2}, \bar{c_3}, \bar{c_4})$. $\lambda$ and $\bar{c_i}$ are constants. We obtain the characteristic equations of $E_n$.

First, we analyze equilibria $E_1(1, 0, 0, 0)$

The equation has four solutions of $\lambda$ given by

$\lambda_1 = -\mu < 0$,

$\lambda_2 = \beta-(\gamma+\mu)$, when $\mathscr{R}_0 < 1$, $\lambda_2 < 0$,

$\lambda_3 = -(\phi+\mu) < 0$,

$\lambda_4 = A_x(E_1) = k \left(\frac{1}{1+ e ^{-b\overline{g(t)}}}-\frac{1}{2}\right) < 0$.

When all solutions $\lambda_i$ have negative real parts, we can judge that this point is locally asymptotically stable in this system. Hence, we have:

Proposition 1. $E_1$ is locally asymptotically stable if $\mathscr{R}_0 < 1$.

We consider $E_2(\frac{\gamma+\mu}{\beta}, \frac{\mu(\beta-\gamma-\mu)}{\beta(\gamma+\mu)}, 0, 0)$, we obtain

A solution $\lambda_1 = k \left(\frac{1}{1+ e ^{-b\overline{g(t)}}}-\frac{1}{2}\right)$, when $\mathscr{R}_0 < R$, $\lambda_1 < 0$. Actually, $\mathscr{R}_0 < R$ is equal to:

This means that: $g(t) > 0$.

$\lambda_2 = -(\phi+\mu) < 0$.

To determine if the real part of these two roots, $\lambda_3$ and $\lambda_4$, is less than 0, we use a simpler method. We solve this equation

or equivalently

where $a_1 = \beta i_2+\mu, a_2 = \beta^2s_2i_2$. We conclude that $a_1 > 0, a_2 > 0$ if $\mathscr{R}_0 > 1$. According to Routh-Hurwitz criterion, we have the following results:

Proposition 2. $E_2$ is locally asymptotically stable if $1 < \mathscr{R}_0 < R$.

Subsequently, we consider $E_3(\frac{\mu + \phi}{\mu + \phi + \sigma}, 0, \frac{\sigma}{ \mu + \phi + \sigma}, 1)$, we obtain:

Note that $\lambda_1 = A_x(E_3) = k \left(\frac{1}{1+ e ^{-b\overline{g(t)}}}-\frac{1}{2}\right) > 0$ owing to $\overline{g(t)} > 0$. Hence, we have:

Proposition 3. $E_3$ is always not stable.

Consider $E_4(\frac{1}{\mathscr{R}_0}, i_4, \frac{\sigma}{\mathscr{R}_0 (\mu+ \phi)}, 1)$:

A solution $\lambda_1 = A_x(E_4) > 0$ if $\mathscr{R}_v > R$. We discuss the other three solutions.

This is equal to $\lambda^3+b_1\lambda^2+b_2\lambda+b_3 = 0$ where $b_1 = (\mu+\phi)+(\beta i_4+\mu+\sigma)$, $b_2 = (\phi \gamma +(\mu+\phi)(\beta i_4+\mu+\sigma)+\frac{\beta^2i_4}{\mathscr{R}_0})$, $b_3 = (\mu+\phi)\frac{\beta^2i_4}{\mathscr{R}_0}$. We observe that: $b_1 > 0, b_2 > 0, b_3 > 0$.

By calculating, we obtain

then, we conclude that:

Proposition 4. $E_4$ is locally asymptotically stable if $\mathscr{R}_v > 1$.

Finally, we consider $E_5(\frac{1}{\mathscr{R}_0}, i_5, v_5, x_5)$:

This is equal to $\lambda^4+c_1\lambda^3+c_2\lambda^2+c_3\lambda+c_4 = 0$. Here,

By calculating, we obtain:

and we set

Hence, we obtain:

Proposition 5. $E_5$ is locally asymptotically stable when $\mathscr{R}_v < R < \mathscr{R}_0$ and $F_1 > 0, F_2 > 0$.

4.4.   Ulam-Hyers global stability

Next, we turn to the fractional model. The study of Ulam-Hyers stability has its origins in a question posed by S. Ulam in 1940 as to whether an exact solution must necessarily exist in the neighborhood of an approximate solution to a given equation ^[35]. We set $\bar{X}$ is a solution of (3.19), and define the function $\bar{X}\in E$ is a solution of (3.19) if and only if

The definition of Ulam-Hyers stability of the system (3.19) is

Definition 4. If there exists $C_K > 0$ such that for any $\varepsilon > 0$, there exists a function $X\in E$ for every solution $\bar{X}\in E$ that satisfies the system (3.19), and

Then the fractional-order model (3.19) is Ulam-Hyers stable.

We also define the generalized Ulam-Hyers stability by

Definition 5. If there exists $\psi_K:\mathbb{R}^+\to \mathbb{R}^+$ with $\psi_K(0) = 0$ such that for any $\varepsilon > 0$, there exists a function $X\in E$ for every solution $\bar{X}\in E$ that satisfies the system (3.19), and there are

Then, the fractional-order model (3.19) is generalized Ulam-Hyers stable.

We prove that the model (3.19) is Ulam-Hyers stable and generalized Ulam-Hyers stable.

Theorem 5. Let the hypothesis of Theorem 2 hold, then the fractional-order model (3.19) is Ulam-Hyers stable and generalized Ulam-Hyers stable.

Proof. Let $X$ be the only solution of model (3.19), and $\bar{X}$ is the solution that satisfies (4.32).

By the simple calculation and condition (4.32), we obtain

consequently,

We obtain $\|\bar{X}-X\|_{E} \leq C \varepsilon$, $C = \frac{M}{1-M L}$. Hence, the model is Ulam-Hyers stable. Considering the continuity of $M$ and $L$, we conclude that it is generalized Ulam-Hyers stable. The condition set here holds on domain $\Lambda$, so the solution is globally stable. □

5.   Simulation and application

In this section, in order to provide evidence for the study discussed in this paper, we use examples and models of influenza transmission, along with data from the paper, to demonstrate the impact of vaccination on influenza and its influence on model dynamics. Moreover, we verify the conclusions that we reached earlier and subsequently derive a strategy for vaccination. By conducting these simulations, we gain insights into how fractional derivatives and parameters fine-tuning affect the model's forecasts, thereby fostering a more profound comprehension of the model's dynamic properties.

5.1.   Fractional-order model simulation and discussion

We use Grünwald-Letnikov method for simulation ^[36,37]. To introduce our numerical simulation, we utilize the parameters of Example 1.

Example 1. We set $N = 1000, \mu = 0.05, \beta = 3.5, \gamma = 0.6, c_v = 1, b = 0.8, \phi = 0.05, c_i = 30, m = 0.9, k = 0.5, \sigma = 0.8,$ in the fractional model.

It can be seen that the disease erupts rapidly at $t < 10$ and then tends to a more stable state. In addition, we can see how the values change with different values of $\alpha$ in Figure 3. Here, the value of $\alpha$ clearly influences the state of model change. By varying the fractional order $\alpha$, we can identify a curve that best fits the actual data. This is also one of the advantages of the fractional system. That is fractional-order models have one more degree of freedom than integer order models. In fact, the model has some other properties which we will analyse below. We can verify that it satisfies the Lyapunov stability condition, i.e.,

Next, we take $x(0)$ as the boundary value, that is $x(0) = 0$, $x(0) = 1$ and $x(0) = 0.01$. Here $\alpha = 0.95$. Figure 4 verifies the stability of the equilibrium point. Additionally, it is possible to visualize the effect of changes in vaccination rates on the overall model. The effectiveness of vaccination for disease control is verified from the change in $I(t)$. Figure 5 illustrates the impact of vaccination on the prevalence of the disease.

Through a comprehensive analysis and comparison of these figures, we can clearly observe the significant changes that occurred before and after the implementation of the vaccination program. This leads us to the strong conclusion that a scientific, standardized, and widespread vaccination strategy is essential for effectively controlling the spread of the epidemic and significantly interrupting the chain of transmission of the disease. This reveals that as the vaccinated population expands, the number of new infections declines markedly, and the scale of outbreaks is effectively curtailed. When the vaccination strategy is implemented, even in the presence of an infectious individual, the disease can be maintained at a low prevalence, thereby alleviating pressure on the healthcare system.

5.2.   Effect of parameters on vaccination rate and disease control

We simulate $x$ to see the effect of different factors on vaccination rates. First, we simulate the Example 1 above for $\alpha = 1$ in Figure 6. We find vaccination rates oscillating.

These fluctuations not only reflect the dynamic changes in the public's willingness to receive vaccinations but also highlight the complexities and challenges associated with vaccination efforts. This observation aligns with our objective understanding of social behavior patterns and the effectiveness of public health interventions. As vaccination rates reach a certain high level, public enthusiasm for vaccination tends to diminish. Following the initial control of the epidemic, some individuals may develop an optimistic belief that the crisis has passed, leading to a decreased sense of urgency regarding continued vaccination. Furthermore, those who have been vaccinated may, after experiencing the short-term protection offered by the vaccine, underestimate the importance of its long-term benefits, thereby increasing the risk of vaccine expiration. This phenomenon is consistent with the observations derived from the model. Next, we observe the effect of fractional order $\alpha$ on the model. In continuously monitoring vaccination rates, it is evident that vaccination rates do not exhibit a simple linear increase or a stable state; instead, they display a pronounced oscillating pattern over time. We find that as the fractional order $\alpha$ decreases, the oscillations also decrease and converge to the equilibrium more quickly in Figure 7.

The fact that different initial values of $x_0$ also give the same convergence results for the model verifies the existence of internal equilibrium points in Figure 8.

The cost of vaccines are undoubtedly crucial considerations. This is not only directly related to the accessibility of vaccines but also have a profound impact on vaccination rates. Therefore, how vaccine costs jointly influence vaccination rates is essential for formulating scientific and rational vaccination strategies, optimizing resource allocation, and enhancing public health. Specifically, vaccine cost is a critical factor affecting the public's willingness and ability to vaccinate. When vaccine costs are low, more families and individuals can afford vaccinations, leading to a substantial increase in vaccination rates. The Figure 9 clearly indicates that as vaccine costs decrease, vaccination rates can reach relatively high levels. This phenomenon suggests that reducing vaccine costs is an effective strategy for expanding vaccination coverage and improving herd immunity.

For epidemic control, the most important thing is the internal equilibrium point. Next, we study the effect of parameters on the internal equilibrium point. We use influenza transmission data and models to illustrate the implementation of vaccination plans, establishment of immune barriers, and attainment of vaccine control effects. We put the set parameters into Table 1.

We examine the effect of individual parameters on the scale of the disease. First, we need to verify the conditions under which the equilibrium point exists.

If the cost of vaccination is higher, it will lead to a greater reluctance for people to get vaccinated in Figure 10(a). Similarly, if the cost of infection is low, it can also lead to the situation depicted in Figure 10(b). The graph also illustrates that the proportion of infections $i_5$ is an increasing function of $m$. This proves our assumption. When the cost of vaccination is very low (let $c_v = 0$), people are more likely to choose to get vaccinated. At this time, the infection rate is $i_5 = 0$. For influenza control, it is crucial to reduce the cost of vaccines.

Finally, we can conclusively state that reducing the cost of vaccines and increasing public awareness are not only essential but also pivotal strategies for effectively controlling outbreaks. By making vaccines more affordable, we can ensure broader access and higher vaccination rates among vulnerable populations, thereby significantly mitigating the spread of infectious diseases. Additionally, enhancing public awareness through educational campaigns and informative media outlets empowers individuals to adopt preventive measures and receive vaccinations, further reinforcing our collective defenses against influenza pandemics.

5.3.   Discussion

In this section, we discuss the results obtained in the previous sections. Through a comprehensive analysis and comparison of these figures, we can clearly observe the significant changes that occurred before and after the implementation of the vaccination program. Our findings underscore the essential role of a scientific, standardized, and widespread vaccination strategy in effectively controlling the spread of the epidemic and substantially interrupting the chain of disease transmission. The simulation reveals a compelling trend: as the vaccinated population increases, the number of new infections declines markedly, and the scale of outbreaks is effectively curtailed. This demonstrates that, even in the presence of infectious individuals, a robust vaccination strategy can maintain the disease at a low prevalence, thereby alleviating considerable pressure on healthcare systems.

The cost of vaccines is a crucial factor in this analysis. It is not only directly related to the accessibility of vaccines but also significantly impacts vaccination rates. Specifically, the cost of vaccines plays a vital role in determining the public's willingness and ability to vaccinate. When vaccine costs are low, more families and individuals can afford vaccinations, resulting in a substantial increase in vaccination rates. This indicates that reducing vaccine costs is an effective strategy for expanding vaccination coverage and improving herd immunity. Conversely, when the cost of vaccination is elevated, it results in increased reluctance among individuals to receive vaccinations. Lower costs promote wider access and higher vaccination rates among vulnerable populations, which significantly reduces the spread of infectious diseases. However, merely reducing costs is not enough; enhancing public awareness is equally crucial. Educational campaigns and informative media can empower individuals to adopt preventive measures and obtain vaccinations, thereby strengthening our collective defenses against influenza pandemics. By integrating affordability with awareness, we can ensure that vaccination programs are both accessible and effective.

In conclusion, reducing the cost of vaccines and increasing public awareness are not only essential but also pivotal strategies for effectively controlling outbreaks. By making vaccines more affordable, we can broaden access and elevate vaccination rates among those most susceptible to infection, significantly mitigating the spread of disease. Additionally, enhancing public awareness through targeted educational initiatives reinforces preventive behaviors and encourages vaccination uptake, further strengthening our collective resilience against influenza and other infectious diseases. Together, these strategies create a robust framework for public health intervention and epidemic control.

6.   Conclusions

We propose a novel Caputo fractional model that incorporates vaccination strategies. We have proven the existence and uniqueness of the solution by applying theories in Banach Space to this model and have analyzed the stability of Lyapunov and Ulam-Hyers. Additionally, we conducted numerical simulations with various parameter values. Due to the freedom of the $\alpha$, the findings indicate that the fractional model significantly enhances the accuracy of disease modeling. We analyze the oscillatory patterns of vaccination rates and determine their relevance. The effectiveness of this model in predicting disease outbreaks has been demonstrated. We analyzed the internal equilibrium to explore the relationship between disease prevalence and various parameters. This study is highly significant for improving the accuracy of disease models and developing more effective strategies to control the spread of influenza.

In the future, we will also explore the incorporation of additional real-world factors into the model to more accurately simulate the dynamics of influenza transmission, for example, the effect of social distance. We could consider incorporating stochastic dynamics into the model to enhance its realism in simulations. Additionally, parameter estimation techniques can be employed to determine specific parameters, thereby facilitating the creation of scenarios for disease control. Furthermore, optimal control methods can be utilized in this study to devise a more effective strategy.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

The authors would like to thank the School of Mathematics and Statistics at Changchun University of Science and Technology for the support of this study. This work was supported by Science and Technology Development Plan Project of Jilin Province, China (No. 20230101291JC). China Scholarship Council (No. 201807585008).

Conflict of interest

The authors declare there is no conflicts of interest.