Advanced analysis in epidemiological modeling: detection of waves

Abdon Atangana; Seda İğret Araz; Abdon Atangana; Seda İğret Araz

doi:10.3934/math.2022992

AIMS Mathematics

2022, Volume 7, Issue 10: 18010-18030. doi: 10.3934/math.2022992

Previous Article Next Article

Research article Special Issues

Advanced analysis in epidemiological modeling: detection of waves

Abdon Atangana ^{1,2
,
,},
Seda İğret Araz ^1,3

1.
Faculty of Natural and Agricultural Sciences, University of the Free State, South Africa
2.
Department of Medical Research, China Medical University Hospital, China Medical University, Taichung, Taiwan
3.
Faculty of Education, Siirt University, Siirt 56100, Turkey

Received: 16 February 2022 Revised: 23 June 2022 Accepted: 18 July 2022 Published: 08 August 2022
MSC : 26A33, 65C40, 92D30

Mathematical concepts have been used in the last decades to predict the behavior of the spread of infectious diseases. Among them, the reproductive number concept has been used in several published papers to study the stability of the mathematical model used to predict the spread patterns. Some conditions were suggested to conclude if there would be either stability or instability. An analysis was also meant to determine conditions under which infectious classes will increase or die out. Some authors pointed out limitations of the reproductive number, as they presented its inability to help predict the spread patterns. The concept of strength number and analysis of second derivatives of the mathematical models were suggested as additional tools to help detect waves. This paper aims to apply these additional analyses in a simple model to predict the future.

Keywords:

Citation: Abdon Atangana, Seda İğret Araz. Advanced analysis in epidemiological modeling: detection of waves[J]. AIMS Mathematics, 2022, 7(10): 18010-18030. doi: 10.3934/math.2022992

Related Papers:

[1]	Said Gounane, Jamal Bakkas, Mohamed Hanine, Gyu Sang Choi, Imran Ashraf . Generalized logistic model with time-varying parameters to analyze COVID-19 outbreak data. AIMS Mathematics, 2024, 9(7): 18589-18607. doi: 10.3934/math.2024905
[2]	I. H. K. Premarathna, H. M. Srivastava, Z. A. M. S. Juman, Ali AlArjani, Md Sharif Uddin, Shib Sankar Sana . Mathematical modeling approach to predict COVID-19 infected people in Sri Lanka. AIMS Mathematics, 2022, 7(3): 4672-4699. doi: 10.3934/math.2022260
[3]	Khalid Khan, Amir Ali, Muhammad Irfan, Zareen A. Khan . Solitary wave solutions in time-fractional Korteweg-de Vries equations with power law kernel. AIMS Mathematics, 2023, 8(1): 792-814. doi: 10.3934/math.2023039
[4]	Maha Alammari, Muhammad Abuzar, Solomon Manukure . Bifurcation analysis and the modulation instability in a nonlinear silica optical fibers. AIMS Mathematics, 2025, 10(7): 16692-16719. doi: 10.3934/math.2025749
[5]	Faeza Lafta Hasan, Mohamed A. Abdoon, Rania Saadeh, Ahmad Qazza, Dalal Khalid Almutairi . Exploring analytical results for (2+1) dimensional breaking soliton equation and stochastic fractional Broer-Kaup system. AIMS Mathematics, 2024, 9(5): 11622-11643. doi: 10.3934/math.2024570
[6]	Aly R. Seadawy, Bayan Alsaedi . Contraction of variational principle and optical soliton solutions for two models of nonlinear Schrödinger equation with polynomial law nonlinearity. AIMS Mathematics, 2024, 9(3): 6336-6367. doi: 10.3934/math.2024309
[7]	Ramsha Shafqat, Ateq Alsaadi . Artificial neural networks for stability analysis and simulation of delayed rabies spread models. AIMS Mathematics, 2024, 9(12): 33495-33531. doi: 10.3934/math.20241599
[8]	Ajmal Ali, Tayyaba Akram, Azhar Iqbal, Poom Kumam, Thana Sutthibutpong . A numerical approach for 2D time-fractional diffusion damped wave model. AIMS Mathematics, 2023, 8(4): 8249-8273. doi: 10.3934/math.2023416
[9]	Noha M. Kamel, Hamdy M. Ahmed, Wafaa B. Rabie . Solitons unveilings and modulation instability analysis for sixth-order coupled nonlinear Schrödinger equations in fiber bragg gratings. AIMS Mathematics, 2025, 10(3): 6952-6980. doi: 10.3934/math.2025318
[10]	Naseem Abbas, Amjad Hussain, Mohsen Bakouri, Thoraya N. Alharthi, Ilyas Khan . A study of dynamical features and novel soliton structures of complex-coupled Maccari's system. AIMS Mathematics, 2025, 10(2): 3025-3040. doi: 10.3934/math.2025141

Abstract

1. Introduction

Mathematical models of the spread of disease have been used successfully with limitations to predict future behaviors of the spread. The main aim is to have an asymptotic idea of what the spread will look like early enough that measures can be taken to help control and eradicate the spread of the virus. Several mathematical concepts have been suggested to help better analyze these models. During the modeling process, mathematicians or modelers first translate the observed facts into mathematical equations using either ordinary or partial differential equations ^[1,2,3,4]. The obtained system of equations, either linear or non-linear, is further analyzed. The first analysis consists of obtaining equilibrium points (disease-free and endemic). This analysis is often followed by an analysis of stability that includes global and asymptotic global stability, conditions under which the infectious classes will decrease, increase, or stay constant. The reproductive number, which can be obtained differently, for instance, using the next generation of the matrix, is obtained to help predict whether the spread will be stable or not. Then, finally, the study of the sign of the first derivative of the Lyapunov function is associated with the system to give a clear idea of the stability. In thousands of papers published on this topic, similar analyses can be found; so far, no significant contribution has been made after the idea of reproductive number was suggested. However, while checking with great care predictions made using this analysis, one will quickly see that this analysis cannot help humans predict waves in a spread. For example, several mathematical models have been suggested to predict the spread of Covid-19, but no analysis helps to indicate the number of waves that humans will face by Covid-19 spread ^{[5,6,7,8,9,10,11,12]}. Some researchers showed that the reproductive number was insufficient to provide an apparent behavior of the spread. However, no attention was paid to such analysis; one will agree that suggested mathematical models and research done in the last decades have not helped to predict waves. Very recently, additional sets of analyses were presented: first, an evaluation of the equilibrium point of the second derivative of the model; an evaluation of conditions under which the second derivatives of infections classes are positive, negative or zero, and a calculation of strength number using the next-generation matrix and evaluation of the sign of the second derivative of the Lyapunov function associated with the model. Additionally, it was also suggested that the model should be derived with piecewise derivatives; in this paper, we will apply that analysis to a simple epidemiological model.

This study is organized as follows. In Section 2, we present a deep analysis, including positiveness and boundedness of the solutions of an SEIR(Suscepted-Exposed-Infected-Recovered) model, reproduction and strength number and the first and second derivatives of the Lyapunov function. Moreover, using the second derivative, we aim to determine under what conditions waves are possible. In Section 3, we present a numerical scheme based on Newton polynomials for an SEIR model with a piecewise derivative, where the first part is a classical derivative, and the second part is a stochastic Caputo-Fabrizio derivative, and the numerical simulation is depicted for such model. In Section 4, we present the numerical solution of an SEIR model with piecewise derivative where the first part is a classical derivative, and the second part is a stochastic Atangana-Baleanu derivative, and the numerical simulation is performed for such model. In Section 5, the numerical solution of an SEIR model with piecewise derivative is presented. Here, the model is constructed such that the first part is a classical derivative and the second part is a stochastic Caputo derivative. The numerical simulation is also depicted for such a model.

2. A SEIR model

The SEIR model has been used and analyzed in many research works to predict the spread of some infectious diseases. This model has been used to calculate the so-called reproductive number; also, researchers have calculated its respective equilibrium points, including disease-free and endemic points. However, the model has not been used to check whether the model can predict waves or if the equilibrium points are local maximums or local minimums. In this section, an SEIR model will be subjected to a further analysis that may help determine if such a model is suitable for describing spread with waves.

$\begin{eqnarray} \frac{dS\left( t\right) }{dt} & = &\mu N-\mu S-\frac{\beta SI}{N} \\ \frac{dE\left( t\right) }{dt} & = &\frac{\beta SI}{N}-\left( \tau +\mu \right) E \\ \frac{dI\left( t\right) }{dt} & = &\tau E-\left( \gamma +\mu \right) I \\ \frac{dR\left( t\right) }{dt} & = &\gamma I-\mu R \end{eqnarray}$

(2.1)

where $N = S+E+I+R$ . The SEIR model will be subjected to the following initial conditions:

$\begin{equation} S\left( 0\right) = S^{0}, E\left( 0\right) = E^{0}, I\left( 0\right) = I^{0}, R\left( 0\right) = R^{0}. \end{equation}$

(2.2)

Here $S\left(t\right)$ is the class of the susceptible individuals. $E\left(t\right)$ is the class of the exposed individuals. $I\left(t\right)$ is the class of the infected individuals. $R\left(t\right)$ is the class of the recovered individuals. $\mu N$ is the death rate, $\beta$ is the contact rate, $\tau$ is the rate of progression from exposed to infectious, and $\gamma$ is the recovery rate of infectious individuals ^[13]. Note that

$\begin{equation} \frac{dN}{dt} = \mu N-\mu N = 0. \end{equation}$

(2.3)

We should point out that this model does not consider the number of deaths, and it is assumed that the total population is constant as a function of time.

2.1. Positiveness and boundedness of solutions

We start with a primary analysis, at least to show that the solutions are positive, since they depict real-world problems with positive values. In this subsection, we examine the conditions under which the positivity of the solutions of the considered model is satisfied. Let us start with the class $E\left(t\right)$ that verifies

$\begin{equation} E\left( t\right) \geq E^{0}e^{-\left( \tau +\mu \right) t}, \ \ \forall t\geq 0. \end{equation}$

(2.4)

For the function $I\left(t\right),$ we have the following inequalities:

$\begin{equation} I\left( t\right) {\rm{ }}\geq I^{0}e^{-\left( \gamma +\mu \right) t}, \ \ \forall t\geq 0, \end{equation}$

(2.5)

and

$\begin{equation} R\left( t\right) \geq R^{0}e^{-\mu t}, \ \ \forall t\geq 0. \end{equation}$

(2.6)

We shall define the norm

$\begin{equation} \left\Vert \lambda \right\Vert _{\infty } = \sup\limits_{t\in D_{\lambda }}\left\vert \lambda \left( t\right) \right\vert \end{equation}$

(2.7)

where $D_{\lambda }$ is the domain of $\lambda$ . Using the above norm, we get the following inequality for the function $S\left(t\right),$

$\begin{eqnarray} \overset{\cdot }{S}\left( t\right) {\rm{ }} & = &\mu N-\mu S-\frac{\beta SI}{N} , \ \ \forall t\geq 0, \\ &\geq &-\left( \mu +\frac{\beta \left\vert I\right\vert }{\left\vert N\right\vert }\right) S, \ \ \forall t\geq 0, \\ &\geq &-\left( \mu +\frac{\beta \sup\limits_{t\in D_{I}}\left\vert I\right\vert }{ \sup\limits_{t\in D_{N}}\left\vert N\right\vert }\right) S, \ \ \forall t\geq 0, \\ &\geq &-\left( \mu +\frac{\beta \left\Vert I\right\Vert _{\infty }}{ \left\Vert N\right\Vert _{\infty }}\right) S, \ \ \forall t\geq 0. \end{eqnarray}$

(2.8)

This yields

$\begin{equation} S\left( t\right) \geq S^{0}e^{-\left( \mu +\frac{\beta \left\Vert I\right\Vert _{\infty }}{\left\Vert N\right\Vert _{\infty }}\right) t}, \ \ \forall t\geq 0. \end{equation}$

(2.9)

2.2. Analysis of equilibrium points

In this subsection, we give a detailed analysis about equilibrium points. The disease-free equilibrium for this model is

$\begin{equation} \left( N, 0, 0, 0, 0\right) . \end{equation}$

(2.10)

To get the endemic equilibrium points, we need to solve the following system:

$\begin{equation} \left\{ \begin{array}{c} \mu N-\mu S-\frac{\beta SI}{N} = 0 \\ \frac{\beta SI}{N}-\left( \tau +\mu \right) E = 0 \\ \tau E-\left( \gamma +\mu \right) I = 0 \\ \gamma I-\mu R = 0 \end{array} \right. . \end{equation}$

(2.11)

Then, the endemic equilibrium points are

$\begin{eqnarray} S^{\ast } & = &\frac{N\left( \mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau \right) }{\tau \beta } \\ E^{\ast } & = &-\frac{\left( \mu ^{2}+\mu \tau +\mu \gamma -\tau \beta +\gamma \tau \right) \mu N}{\tau \beta \left( \tau +\mu \right) } \\ I^{\ast } & = &-\frac{N\mu \left( \mu ^{2}+\mu \tau +\mu \gamma -\tau \beta +\gamma \tau \right) }{\beta \left( \mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau \right) } \\ R^{\ast } & = &-\frac{N\gamma \left( \mu ^{2}+\mu \tau +\mu \gamma -\tau \beta +\gamma \tau \right) }{\beta \left( \mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau \right) }. \end{eqnarray}$

(2.12)

The above endemic equilibrium points are valid if the following inequality holds:

$\begin{equation} -\mu ^{2}-\mu \tau -\mu \gamma -\gamma \tau +\tau \beta > 0, \end{equation}$

(2.13)

namely,

$\begin{equation} 1 > \frac{\mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau }{\tau \beta }. \end{equation}$

(2.14)

2.3. Reproduction number

To proceed with our preliminary analysis, we shall present here the derivation of the so-called reproductive number. This is an essential value in the field of epidemiological modeling, as it helps to give an understanding of the stability conditions. It was documented that if this number is less than 1, we can expect stability; however, if this number is greater than 1, we could expect instability. Nonetheless, it was reported that the value is obtained in different ways. However, we will utilize the next generation matrix approach to determine the reproductive value in this case.

$\begin{eqnarray} \overset{\cdot }{E} & = &\frac{\beta SI}{N}-\left( \tau +\mu \right) E \\ \overset{\cdot }{I} & = &\tau E-\left( \gamma +\mu \right) I. \end{eqnarray}$

(2.15)

Using the next generation matrix approach ^[14], we calculate the matrices $F$ and $V^{-1}$ as

$\begin{equation} F = \left[ \begin{array}{cc} 0 & 0 \\ \frac{\beta S}{N} & 0 \end{array} \right] , V^{-1} = \left[ \begin{array}{cc} \frac{1}{\tau +\mu } & 0 \\ \frac{\tau }{\left( \tau +\mu \right) \left( \gamma +\mu \right) } & \frac{1 }{\gamma +\mu } \end{array} \right] . \end{equation}$

(2.16)

Then, the reproduction number is obtained as

$\begin{equation} R_{0} = \frac{\beta }{\left( \gamma +\mu \right) }. \end{equation}$

(2.17)

2.4. Strength number

The reproductive number has been used in thousands of research papers in epidemiological modeling with some success and significant limitations, but this concept has been employed to determine whether the spread will be severe, some researchers have pointed out some significant weaknesses. For example, they sadly realized that such value is not unique, as it can be obtained via different methods. Another issue raised was that this number should have been a function of time, not a constant value. It was also noticed that a reproductive number could not be used to indicate whether a model will predict waves or not. Very recently, an alternative number was suggested and was called the strength number; of course, this number will be subjected to several tests to see whether it can be used to help detect some complexities in the spread, or at least if this number can help detect waves in a spread ^[4,13]. The value was derived using the next-generation matrix by taking the second derivative of the infectious classes. This section will present the strength number associated with the SEIR model.

$\begin{eqnarray} \frac{\partial ^{2}}{\partial I^{2}}\left( \frac{\beta SI}{N}\right) & = &\beta S\frac{\partial }{\partial I}\left( \frac{\left( N-\overset{.}{N} I\right) }{N^{2}}\right) \\ & = &-\frac{\beta S}{N^{2}}. \end{eqnarray}$

(2.18)

In this case, we can have the following:

$\begin{equation} F_{A} = \left[ \begin{array}{cc} 0 & 0 \\ \frac{-\beta S}{N^{2}} & 0 \end{array} \right] . \end{equation}$

(2.19)

Then,

$\begin{equation} \det \left( F_{A}V^{-1}-\lambda I\right) = 0 \end{equation}$

(2.20)

yields

$\begin{equation} A_{0} = -\frac{S\beta \tau }{N^{2}\left( \tau +\mu \right) \left( \gamma +\mu \right) } < 0. \end{equation}$

(2.21)

Having the strength number be less than zero will lead us to a conclusion that will further be confirmed using a different analysis in finding the local minimum, local maximum, and inflection points. Thus, having a negative strength number is an indication that the SEIR model will have a single magnitude, either a maximum point with two inflection points, indicating a single wave, or an infection that will decrease rapidly from the disease-free equilibrium. With the renewal process, the infection will rise after a minimum point and then be stabilized or stopped later on. This will be confirmed while studying the sign of the second derivative of infectious classes.

2.5. First derivative of Lyapunov

For the endemic equilibrium point $\left(S^{\ast }, E^{\ast }, I^{\ast }, R^{\ast }\right),$ the associated Lyapunov function is analyzed below.

Theorem 1. When the reproductive number $R_{0} > 1$ , the endemic equilibrium points $E^{\ast }$ of the SEIR model are globally asymptotically stable.

Proof. For proof, the Lyapunov function can be written as

$\begin{eqnarray} L\left( S^{\ast }, E^{\ast }, I^{\ast }, R^{\ast }\right) & = &\left( S-S^{\ast }-S^{\ast }\log \frac{S^{\ast }}{S}\right) +\left( E-E^{\ast }-E^{\ast }\log \frac{E^{\ast }}{E}\right) \\ &&+\left( I-I^{\ast }-I^{\ast }\log \frac{I^{\ast }}{I}\right) +\left( R-R^{\ast }-R^{\ast }\log \frac{R^{\ast }}{R}\right) . \end{eqnarray}$

(2.22)

Therefore, applying the derivative with respect to $t$ on both sides yields

$\begin{eqnarray} \frac{dL}{dt} & = &\overset{\cdot }{L} = \left( \frac{S-S^{\ast }}{S}\right) \overset{\cdot }{S}+\left( \frac{E-E^{\ast }}{E}\right) \overset{\cdot }{E} \\ &&+\left( \frac{I-I^{\ast }}{I}\right) \overset{\cdot }{I}+\left( \frac{ R-R^{\ast }}{R}\right) \overset{\cdot }{R}. \end{eqnarray}$

(2.23)

Now, we can write their values for derivatives as follows:

$\begin{eqnarray} \frac{dL}{dt} & = &\left( \frac{S-S^{\ast }}{S}\right) \left( \mu N-\mu S- \frac{\beta SI}{N}\right) \\ &&+\left( \frac{E-E^{\ast }}{E}\right) \left( \frac{\beta SI}{N}-\left( \tau +\mu \right) E\right) \\ &&+\left( \frac{I-I^{\ast }}{I}\right) \left( \tau E-\left( \gamma +\mu \right) I\right) \\ &&+\left( \frac{R-R^{\ast }}{R}\right) \left( \gamma I-\mu R\right) . \end{eqnarray}$

(2.24)

Putting $S = S-S^{\ast },$ $E = E-E^{\ast },$ $I = I-I^{\ast },$ $R = R-R^{\ast }$ leads to

$\begin{eqnarray} \frac{dL}{dt} & = &\left( \frac{S-S^{\ast }}{S}\right) \left( \mu N-\mu \left( S-S^{\ast }\right) -\frac{\beta \left( S-S^{\ast }\right) \left( I-I^{\ast }\right) }{N}\right) \\ &&+\left( \frac{E-E^{\ast }}{E}\right) \left( \frac{\beta \left( S-S^{\ast }\right) \left( I-I^{\ast }\right) }{N}-\left( \tau +\mu \right) \left( E-E^{\ast }\right) \right) \\ &&+\left( \frac{I-I^{\ast }}{I}\right) \left( \tau \left( E-E^{\ast }\right) -\left( \gamma +\mu \right) \left( I-I^{\ast }\right) \right) \\ &&+\left( \frac{R-R^{\ast }}{R}\right) \left( \gamma \left( I-I^{\ast }\right) -\mu \left( R-R^{\ast }\right) \right) . \end{eqnarray}$

(2.25)

We can organize the above as follows:

$\begin{gather} \frac{dL}{dt} = \mu N-\mu N\frac{S^{\ast }}{S}-\frac{\left( S-S^{\ast }\right) ^{2}}{S}\mu -\frac{\left( S-S^{\ast }\right) ^{2}}{S}\frac{\beta I}{N}+\frac{ \left( S-S^{\ast }\right) ^{2}}{S}\frac{\beta I^{\ast }}{N} \\ -\frac{E^{\ast }}{E}\frac{\beta }{N}SI+\frac{E^{\ast }}{E}\frac{\beta }{N} SI^{\ast }+\frac{E^{\ast }}{E}\frac{\beta }{N}S^{\ast }I-\frac{E^{\ast }}{E} \frac{\beta }{N}S^{\ast }I^{\ast } \\ -\frac{\left( E-E^{\ast }\right) ^{2}}{E}\left( \tau +\mu \right) -\frac{ \left( I-I^{\ast }\right) ^{2}}{I}\left( \gamma +\mu \right) +\tau E-\tau E^{\ast } \\ -\frac{I^{\ast }}{I}\tau E+\frac{I^{\ast }}{I}\tau E^{\ast }+\gamma I-\gamma I^{\ast }-\frac{R^{\ast }}{R}\gamma I+\frac{R^{\ast }}{R}\gamma I^{\ast }- \frac{\left( R-R^{\ast }\right) ^{2}}{R}\mu. \end{gather}$

(2.26)

To avoid the complexity, the above can be written as

$\begin{equation} \frac{dL}{dt} = \Sigma -\Omega \end{equation}$

(2.27)

where

$\begin{gather} \Sigma = \mu N+\frac{\left( S-S^{\ast }\right) ^{2}}{S}\frac{\beta I^{\ast }}{ N}+\frac{E^{\ast }}{E}\frac{\beta }{N}SI^{\ast }+\frac{E^{\ast }}{E}\frac{ \beta }{N}S^{\ast }I \\ +\tau E+\frac{I^{\ast }}{I}\tau E^{\ast }+\gamma I+\frac{R^{\ast }}{R}\gamma I^{\ast } \end{gather}$

(2.28)

and

$\begin{gather} \Omega = \mu N\frac{S^{\ast }}{S}+\frac{\left( S-S^{\ast }\right) ^{2}}{S}\mu +\frac{\left( S-S^{\ast }\right) ^{2}}{S}\frac{\beta I}{N}+\frac{E^{\ast }}{E }\frac{\beta }{N}SI+\frac{E^{\ast }}{E}\frac{\beta }{N}S^{\ast }I^{\ast } \\ +\frac{\left( E-E^{\ast }\right) ^{2}}{E}\left( \tau +\mu \right) +\frac{ \left( I-I^{\ast }\right) ^{2}}{I}\left( \gamma +\mu \right) +\tau E^{\ast } \\ +\frac{I^{\ast }}{I}\tau E+\gamma I^{\ast }+\frac{R^{\ast }}{R}\gamma I+ \frac{\left( R-R^{\ast }\right) ^{2}}{R}\mu . \end{gather}$

(2.29)

It is concluded that if $\Sigma < \Omega,$ this yields $\frac{dL}{dt} < 0;$ however, when $S = S^{\ast }, E = E^{\ast }, I = I^{\ast }, R = R^{\ast }$

$\begin{equation} 0 = \Sigma -\Omega { \ \ \ \ }\Rightarrow \ \ \frac{dL}{dt} = 0. \end{equation}$

(2.30)

We can see that the largest compact invariant space( $\Gamma$ ) for the suggested model in

$\begin{equation} \left\{ \left( S^{\ast }, E^{\ast }, I^{\ast }, R^{\ast }\right) \in \Gamma : \frac{dL}{dt} = 0\right\} \end{equation}$

(2.31)

is the point $\left\{ E_{\ast }\right\}$ the endemic equilibrium of the considered model. By the help of LaSalle's invariance principle^[15], it follows that $E_{\ast }$ is globally asymptotically stable in $\Gamma$ if $\Sigma < \Omega.$

2.6. Second derivative of Lyapunov

While the sign of the first derivative of a Lyapunov function helps to describe the stability, one should note that the analysis of the first derivative of a given function could not fully help one to understand the variabilities of the function under study; for example, one will need the second derivative analysis to find inflection points, local maximums, and minimums. Further research is required on all the particularities of the variations. We thus present an analysis of the second derivative of the associated Lyapunov function of our model.

$\begin{eqnarray} \frac{d\overset{\cdot }{L}}{dt} & = &\frac{d}{dt}\left\{ \begin{array}{c} \left( 1-\frac{S^{\ast }}{S}\right) \overset{\cdot }{S}+\left( 1-\frac{ E^{\ast }}{E}\right) \overset{\cdot }{E}+\left( 1-\frac{I^{\ast }}{I}\right) \overset{\cdot }{I} \\ +\left( 1-\frac{R^{\ast }}{R}\right) \overset{\cdot }{R} \end{array} \right\} \\ & = &\left( \frac{\overset{\cdot }{S}}{S}\right) ^{2}S^{\ast }+\left( \frac{ \overset{\cdot }{E}}{E}\right) ^{2}E^{\ast }+\left( \frac{\overset{\cdot }{I} }{I}\right) ^{2}I^{\ast }+\left( \frac{\overset{\cdot }{R}}{R}\right) ^{2}R^{\ast }+\left( 1-\frac{S^{\ast }}{S}\right) \overset{\cdot \cdot }{S} \\ &&+\left( 1-\frac{E^{\ast }}{E}\right) \overset{\cdot \cdot }{E}+\left( 1- \frac{I^{\ast }}{I}\right) \overset{\cdot \cdot }{I}+\left( 1-\frac{R^{\ast } }{R}\right) \overset{\cdot \cdot }{R}. \end{eqnarray}$

(2.32)

Here, considering $\overset{\cdot }{N} = 0,$ we have

$\begin{eqnarray} \overset{\cdot \cdot }{S} & = &\mu \overset{\cdot }{N}-\mu \overset{\cdot }{S} -\beta \left( \frac{\left( \overset{\cdot }{S}I+\overset{\cdot }{I}S\right) }{N}\right) \\ \overset{\cdot \cdot }{E} & = &\beta \left( \frac{\left( \overset{\cdot }{S}I+ \overset{\cdot }{I}S\right) }{N}\right) -\left( \tau +\mu \right) \overset{ \cdot }{E} \\ \overset{\cdot \cdot }{I} & = &\tau \overset{\cdot }{E}-\left( \gamma +\mu \right) \overset{\cdot }{I} \\ \overset{\cdot \cdot }{R} & = &\gamma \overset{\cdot }{I}-\mu \overset{\cdot }{ R}. \end{eqnarray}$

(2.33)

Then, we have

$\begin{eqnarray} \frac{d\overset{\cdot }{L}}{dt} & = &\left( \frac{\overset{\cdot }{S}}{S} \right) ^{2}S^{\ast }+\left( \frac{\overset{\cdot }{E}}{E}\right) ^{2}E^{\ast }+\left( \frac{\overset{\cdot }{I}}{I}\right) ^{2}I^{\ast }+\left( \frac{\overset{\cdot }{R}}{R}\right) ^{2}R^{\ast } \\ &&+\left( 1-\frac{S^{\ast }}{S}\right) \left( -\mu \overset{\cdot }{S}-\beta \left( \frac{\left( \overset{\cdot }{S}I+\overset{\cdot }{I}S\right) }{N} \right) \right) \\ &&+\left( 1-\frac{E^{\ast }}{E}\right) \left( \beta \left( \frac{\left( \overset{\cdot }{S}I+\overset{\cdot }{I}S\right) }{N}\right) -\left( \tau +\mu \right) \overset{\cdot }{E}\right) \\ &&+\left( 1-\frac{I^{\ast }}{I}\right) \left( \tau \overset{\cdot }{E} -\left( \gamma +\mu \right) \overset{\cdot }{I}\right) \\ &&+\left( 1-\frac{R^{\ast }}{R}\right) \left( \gamma \overset{\cdot }{I}-\mu \overset{\cdot }{R}\right) \end{eqnarray}$

(2.34)

and

$\begin{eqnarray} \frac{d^{2}L}{dt^{2}} & = &\overset{\cdot }{\Pi }\left( S, E, I, R\right) -\mu \left( 1-\frac{S^{\ast }}{S}\right) \overset{\cdot }{S} \\ &&-\beta \left( 1-\frac{S^{\ast }}{S}\right) \left( \frac{\left( \overset{ \cdot }{S}I+\overset{\cdot }{I}S\right) }{N}\right) \\ &&+\left( 1-\frac{E^{\ast }}{E}\right) \beta \left( \frac{\left( \overset{ \cdot }{S}I+\overset{\cdot }{I}S\right) }{N}\right) -\left( 1-\frac{E^{\ast } }{E}\right) \left( \tau +\mu \right) \overset{\cdot }{E} \\ &&+\left( 1-\frac{I^{\ast }}{I}\right) \tau \overset{\cdot }{E}-\left( \gamma +\mu \right) \left( 1-\frac{I^{\ast }}{I}\right) \overset{\cdot }{I} \\ &&+\left( 1-\frac{R^{\ast }}{R}\right) \gamma \overset{\cdot }{I}-\left( 1- \frac{R^{\ast }}{R}\right) \mu \overset{\cdot }{R}. \end{eqnarray}$

(2.35)

After replacing $\overset{\cdot }{S}\left(t\right), \overset{\cdot }{E} \left(t\right), \overset{\cdot }{I}\left(t\right)$ and $\overset{\cdot }{R }\left(t\right)$ by their formulas from the considered model and putting it all together, we can get

$\begin{eqnarray} \frac{d^{2}L}{dt^{2}} & = &\overset{\cdot }{\Pi }\left( S, E, I, R\right) -\mu \left( 1-\frac{S^{\ast }}{S}\right) \left( \mu N-\mu S-\frac{\beta SI}{N} \right) \\ &&-\beta \left( 1-\frac{S^{\ast }}{S}\right) \left( \frac{\left( \left( \mu N-\mu S-\frac{\beta SI}{N}\right) I+\left( \tau E-\left( \gamma +\mu \right) I\right) S\right) }{N}\right) \\ &&+\left( 1-\frac{E^{\ast }}{E}\right) \beta \left( \frac{\left( \left( \mu N-\mu S-\frac{\beta SI}{N}\right) I+\left( \tau E-\left( \gamma +\mu \right) I\right) S\right) }{N}\right) \\ &&-\left( 1-\frac{E^{\ast }}{E}\right) \left( \tau +\mu \right) \left( \frac{ \beta SI}{N}-\left( \tau +\mu \right) E\right) \\ &&+\left( 1-\frac{I^{\ast }}{I}\right) \tau \left( \frac{\beta SI}{N}-\left( \tau +\mu \right) E\right) -\left( \gamma +\mu \right) \left( 1-\frac{ I^{\ast }}{I}\right) \left( \tau E-\left( \gamma +\mu \right) I\right) \\ &&+\left( 1-\frac{R^{\ast }}{R}\right) \gamma \left( \tau E-\left( \gamma +\mu \right) I\right) -\left( 1-\frac{R^{\ast }}{R}\right) \mu \left( \gamma I-\mu R\right) . \end{eqnarray}$

(2.36)

For simplicity, we can write

$\begin{equation} \frac{d^{2}L}{dt^{2}} = \Omega _{1}-\Omega _{2} \end{equation}$

(2.37)

where

$\begin{eqnarray} \Omega _{1} & = &\overset{\cdot }{\Pi }\left( S, E, I, R\right) +\mu ^{2}S+\mu \frac{\beta SI}{N}+\mu ^{2}\frac{S^{\ast }}{S}N+\beta \frac{S^{\ast }}{S}\mu I \\ &&+\beta \frac{S^{\ast }}{S}\frac{\tau ES}{N}+\frac{\beta ^{2}SI^{2}}{N^{2}} +\beta \frac{\mu SI+\left( \gamma +\mu \right) IS}{N}+\beta \mu I+\beta \frac{\tau ES}{N} \\ &&+\frac{E^{\ast }}{E}\frac{\beta ^{2}SI^{2}}{N^{2}}+\beta \frac{E^{\ast }}{E }\frac{\left( \mu SI+\left( \gamma +\mu \right) IS\right) }{N}+\frac{E^{\ast }}{E}\left( \tau +\mu \right) \frac{\beta SI}{N}+\left( \tau +\mu \right) ^{2}E \\ &&+\tau \frac{\beta SI}{N}+\frac{I^{\ast }}{I}\tau \left( \tau +\mu \right) E+\left( \gamma +\mu \right) ^{2}I+\left( \gamma +\mu \right) \frac{I^{\ast } }{I}\tau E \\ &&+\gamma \tau E+\frac{R^{\ast }}{R}\gamma \left( \gamma +\mu \right) I+ \frac{R^{\ast }}{R}\gamma I\mu +\mu ^{2}R \end{eqnarray}$

(2.38)

and

$\begin{eqnarray} \Omega _{2} & = &\mu ^{2}N+\mu ^{2}S^{\ast }+\mu S^{\ast }\frac{\beta I}{N} +\beta \mu I+\beta \frac{\tau ES}{N}+\beta \frac{S^{\ast }}{S}\left( \frac{ \beta SI^{2}}{N^{2}}+\frac{\mu SI+\left( \gamma +\mu \right) IS}{N}\right) \\ &&+\left( \frac{\beta ^{2}SI^{2}}{N^{2}}+\beta \frac{\mu SI+\left( \gamma +\mu \right) IS}{N}\right) +\frac{E^{\ast }}{E}\beta \left( \mu I+\frac{\tau ES}{N}\right) \\ &&+\left( \tau +\mu \right) \frac{\beta SI}{N}+\frac{E^{\ast }}{E}\left( \tau +\mu \right) \frac{\beta SI}{N}+\tau \left( \tau +\mu \right) E+\frac{ I^{\ast }}{I}\tau \frac{\beta SI}{N} \\ &&+\left( \gamma +\mu \right) \tau E+\frac{I^{\ast }}{I}\left( \gamma +\mu \right) ^{2}I+\left( \gamma \left( \gamma +\mu \right) I\right) \\ &&+\frac{R^{\ast }}{R}\gamma \tau E+\mu \gamma I+\frac{R^{\ast }}{R}\mu ^{2}R. \end{eqnarray}$

(2.39)

It can be seen that

$\begin{eqnarray} {\rm{If }}\ \Omega _{1} & > &\Omega _{2}\ {\rm{ then }}\ \frac{d^{2}L}{dt^{2}} > 0, \\ {\rm{If }}\ \Omega _{1} & < &\Omega _{2}\ {\rm{ then }}\ \frac{d^{2}L}{dt^{2}} < 0, \\ {\rm{If }}\ \Omega _{1} & = &\Omega _{2}\ {\rm{ then }}\ \frac{d^{2}L}{dt^{2}} = 0. \end{eqnarray}$

(2.40)

2.7. Equilibrium points for second order

Since the equilibrium points also help one to get information on curvatures, in this subsection, we aim to find the equilibrium points of the second-order derivative of our solutions. Thus, we write

$\begin{eqnarray} \overset{\cdot \cdot }{S} & = &-\mu ^{2}N+\mu ^{2}S-\beta \mu I+2\beta \mu \frac{SI}{N}+\beta ^{2}S\left( \frac{I}{N}\right) ^{2}-\beta \tau \frac{SE}{N }+\left( \gamma +\mu \right) \beta \frac{SI}{N} \\ \overset{\cdot \cdot }{E} & = &\beta \mu I-\beta \mu \frac{SI}{N}-\beta ^{2}S\left( \frac{I}{N}\right) ^{2}+\beta \frac{S}{N}\tau E-\beta \frac{S}{N} \left( \gamma +\mu \right) I-\left( \tau +\mu \right) \frac{\beta SI}{N} +\left( \tau +\mu \right) ^{2}E \\ \overset{\cdot \cdot }{I} & = &\tau \frac{\beta SI}{N}-\tau \left( \tau +\mu \right) E-\left( \gamma +\mu \right) \tau E+\left( \gamma +\mu \right) ^{2}I \\ \overset{\cdot \cdot }{R} & = &\gamma \tau E-\gamma \left( \gamma +\mu \right) I-\mu \gamma I+\mu ^{2}R. \end{eqnarray}$

(2.41)

The disease-free equilibrium point of the second order derivative of the model is

$\begin{equation} E^{\ast \circ } = \left( N, 0, 0, 0\right) . \end{equation}$

(2.42)

The endemic equilibrium points of the second order derivative of the model are

$\begin{eqnarray} S^{\ast \circ } & = &\frac{N\left( \mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau \right) }{\tau \beta } \\ E^{\ast \circ } & = &-\frac{\left( \mu ^{2}+\mu \tau +\mu \gamma -\tau \beta +\gamma \tau \right) \mu N}{\tau \beta \left( \tau +\mu \right) } \\ I^{\ast \circ } & = &-\frac{N\mu \left( \mu ^{2}+\mu \tau +\mu \gamma -\tau \beta +\gamma \tau \right) }{\beta \left( \mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau \right) } \\ R^{\ast \circ } & = &-\frac{N\gamma \left( \mu ^{2}+\mu \tau +\mu \gamma -\tau \beta +\gamma \tau \right) }{\beta \left( \mu ^{2}+\mu \tau +\mu \gamma +\gamma \tau \right) }. \end{eqnarray}$

(2.43)

For the classes $I\left(t\right), E\left(t\right)$ to increase, we need

$\begin{eqnarray} \overset{.}{E}\left( t\right) & = &\frac{\beta SI}{N}-\left( \tau +\mu \right) E \\ \overset{.}{I}\left( t\right) & = &\tau E-\left( \gamma +\mu \right) I. \end{eqnarray}$

(2.44)

Indeed, for the endemic process, we can consider the following:

$\begin{equation} \left\{ \begin{array}{c} \overset{.}{E}\left( t\right) > 0\Rightarrow \frac{\beta SI}{N} > \left( \tau +\mu \right) E \\ \overset{.}{I}\left( t\right) > 0\Rightarrow \tau E > \left( \gamma +\mu \right) I, \end{array} \right. \end{equation}$

(2.45)

and since $\frac{S}{N} < 1,$

$\begin{equation} \left\{ \begin{array}{c} \beta I > \left( \tau +\mu \right) E \\ \frac{\tau }{\gamma +\mu } > \frac{I}{E} \end{array} \right. . \end{equation}$

(2.46)

Thus,

$\begin{equation} \left\{ \begin{array}{c} \frac{I}{E} > \frac{\left( \tau +\mu \right) }{\beta } \\ \frac{\tau }{\gamma +\mu } > \frac{I}{E} \end{array} \right. \Rightarrow \frac{\left( \tau +\mu \right) }{\beta } < \frac{I}{E} < \frac{\tau }{\gamma +\mu }. \end{equation}$

(2.47)

The $E\left(t\right)$ and $I\left(t\right)$ classes will increase if the following conditions hold:

$\begin{equation} \frac{\left( \tau +\mu \right) }{\beta } < \frac{I}{E} < \frac{\tau }{\gamma +\mu }. \end{equation}$

(2.48)

We are now checking local maximum, local minimum, or inflection points. We have two equilibrium points, including the disease-free and endemic. To achieve this, we consider the second derivative of the infectious classes

$\begin{eqnarray} \overset{..}{E}\left( t\right) & = &\beta \left[ \frac{\overset{.}{S}IN+ \overset{.}{I}SN-\overset{.}{N}SI}{N^{2}}\right] -\left( \tau +\mu \right) \overset{.}{E} \\ \overset{..}{I}\left( t\right) & = &\tau \overset{.}{E}-\left( \gamma +\mu \right) \overset{.}{I}. \end{eqnarray}$

(2.49)

At $E^{0}$ , we have that

$\begin{eqnarray} \overset{..}{E}\left( t\right) & = &0 \\ \overset{..}{I}\left( t\right) & = &0 \end{eqnarray}$

(2.50)

which can be considered as an inflection point of the spread. This means that, at this point, we could see any increase in the number of infectious classes. Also, at $E^{0\ast }$ , endemic equilibrium $\overset{.}{S}\left(t\right) = \overset {.}{E}\left(t\right) = \overset{.}{I}\left(t\right) = 0$ . Thus, $\overset{..}{ E}\left(t\right) = \overset{..}{I}\left(t\right) = 0$ also. The model predicts two inflection points; therefore, it is expected to have a maximum and a minimum point.

We shall evaluate the sign of the $\overset{..}{E}$ and $\overset{..}{I}$ functions.

$\begin{eqnarray} \beta \left[ \frac{\overset{.}{S}IN+\overset{.}{I}SN-\overset{.}{N}SI}{N^{2}} \right] -\left( \tau +\mu \right) \overset{.}{E}{\rm{ }} & > &0 \\ \tau \overset{.}{E}-\left( \gamma +\mu \right) \overset{.}{I}{\rm{ }} & > &0. \end{eqnarray}$

(2.51)

$\begin{eqnarray} \beta \left[ \frac{\left( \mu N-\mu S-\frac{\beta SI}{N}\right) I+\left( \tau E-\left( \gamma +\mu \right) I\right) S}{N}\right] & > &\left( \tau +\mu \right) \frac{\beta SI}{N}-\left( \tau +\mu \right) E{\rm{ }} \\ \tau \frac{\beta SI}{N}-\left( \gamma +\mu \right) ^{2}I & > &\left( \gamma +2\mu +\tau \right) \tau E{\rm{ }}. \end{eqnarray}$

(2.52)

$\begin{eqnarray} \beta \mu I-\beta \mu SI-\frac{\beta ^{2}SI^{2}}{N^{2}}+\frac{\beta S\tau E}{ N}-\frac{\beta \left( \gamma +\mu \right) IS}{N} & > &\left( \tau +\mu \right) \frac{\beta SI}{N}-\left( \tau +\mu \right) ^{2}E{\rm{ }} \\ \left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) I & > &\left( \gamma +2\mu +\tau \right) \tau E{\rm{ }}, \end{eqnarray}$

(2.53)

and

$\begin{eqnarray} \beta \mu I+\frac{\beta S\tau E}{N}+\left( \tau +\mu \right) ^{2}E{\rm{ }} & > &\left( \tau +\mu \right) \frac{\beta SI}{N}+\beta \mu SI+\frac{\beta ^{2}SI^{2}}{N^{2}}+\frac{\beta \left( \gamma +\mu \right) IS}{N} \\ \left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) I & > &\left( \gamma +2\mu +\tau \right) \tau E{\rm{ }}. \end{eqnarray}$

(2.54)

Since $\frac{S}{N} < 1$ , we can write

$\begin{equation} \left\{ \begin{array}{c} \beta \mu I+\beta \tau E+\left( \tau +\mu \right) ^{2}E{\rm{ }} > \frac{\beta S }{N}\left[ \mu I+\frac{\beta I^{2}}{N}+\left( \gamma +\mu \right) I+\left( \tau +\mu \right) I\right] \\ \left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) I > \left( \gamma +2\mu +\tau \right) \tau E \end{array} \right. , \end{equation}$

(2.55)

$\begin{equation} \left\{ \begin{array}{c} \beta \mu \frac{I}{E}+\beta \tau +\left( \tau +\mu \right) ^{2} > \frac{\beta S }{N}\left[ \mu \frac{I}{E}+\frac{\beta I^{2}}{NE}+\left( \gamma +\mu \right) \frac{I}{E}+\left( \tau +\mu \right) \frac{I}{E}\right] \\ \frac{I}{E} > \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) } \end{array} \right. , \end{equation}$

(2.56)

$\begin{equation} \left\{ \begin{array}{c} \beta \tau +\left( \tau +\mu \right) ^{2} > \frac{\beta SI}{NE}\left[ \gamma +3\mu +\tau \right] -\beta \mu \frac{I}{E} \\ \frac{I}{E} > \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) } \end{array} \right. , \end{equation}$

(2.57)

$\begin{equation} \left\{ \begin{array}{c} \beta \tau +\left( \tau +\mu \right) ^{2} > \frac{I}{E}\left[ \frac{\beta S}{N} \left( \gamma +3\mu +\tau \right) -\beta \mu \right] \\ \frac{I}{E} > \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) } \end{array} \right. . \end{equation}$

(2.58)

$\begin{equation} \left\{ \begin{array}{c} \frac{\beta \tau +\left( \tau +\mu \right) ^{2}}{\frac{\beta S}{N}\left( \gamma +3\mu +\tau \right) -\beta \mu } > \frac{I}{E} \\ \frac{I}{E} > \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) } \end{array} \right. \end{equation}$

(2.59)

with $\frac{S}{N} > \frac{\mu }{\left(\gamma +3\mu +\tau \right) }.$ Thus we have the following condition

$\begin{equation} \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \beta \tau +\left( \gamma +\mu \right) ^{2}\right) } < \frac{I}{E} < \frac{\beta \tau +\left( \tau +\mu \right) ^{2}}{\frac{\beta S}{N}\left( \gamma +3\mu +\tau \right) -\beta \mu }. \end{equation}$

(2.60)

On the other hand, we write the following:

$\begin{eqnarray} \overset{..}{E}\left( t\right) & < &0 \\ \overset{..}{I}\left( t\right) & < &0. \end{eqnarray}$

(2.61)

$\begin{eqnarray} \beta \mu I+\left( \tau +\mu \right) ^{2}E{\rm{ }} & < &\left( \tau +\mu \right) \beta I+\beta \mu I+\beta ^{2}I^{2}+\beta \left( \gamma +\mu \right) I \\ \left( \gamma +\mu \right) ^{2}I & < &\left( \gamma +2\mu +\tau \right) \tau E {\rm{ }}, \end{eqnarray}$

(2.62)

and

$\begin{equation} \left\{ \begin{array}{c} \frac{\left( \tau +\mu \right) ^{2}}{\left( \beta ^{2}N+\beta \left( \gamma +\mu \right) \right) } < \frac{I}{E} \\ \frac{I}{E} < \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \gamma +\mu \right) ^{2}} \end{array} \right. . \end{equation}$

(2.63)

Therefore,

$\begin{equation} \left\{ \begin{array}{c} \overset{..}{E}\left( t\right) < 0 \\ \overset{..}{I}\left( t\right) < 0 \end{array} \right. \Rightarrow \frac{\left( \tau +\mu \right) ^{2}}{\left( \beta ^{2}N+\beta \left( \gamma +\mu \right) \right) } < \frac{I}{E} < \frac{\left( \gamma +2\mu +\tau \right) \tau }{\left( \gamma +\mu \right) ^{2}}. \end{equation}$

(2.64)

3. Numerical solution of the model with exponential decay process

This section presents a numerical scheme based on the Newton polynomial for the numerical solution of the considered model with piecewise differential and integral operators ^[16]. While the Caputo Fabrizio differential operator helps to capture processes exhibiting fading memory behaviors, the stochastic, on the other hand, helps capture processes with randomness. There are processes in nature that exhibit both fading memory and randomness. It is worth noting that such methods cannot be accurately modelled either with Caputo-Fabrizio differential operator or with the random, and thus a combination of both is needed. Indeed this comparison is not at all a new calculus or theory but a mathematical tool to replicate this scenario. This behavior has been observed in data collected for different infectious diseases; therefore, in this section, we provide a mathematical tool to capture such behaviors. We start with the piecewise SEIR model with the classical and stochastic Caputo-Fabrizio case, which is given by

$\begin{equation} \left\{ \begin{array}{c} X\prime = F\left( t, X\right) {\rm{, \ if }}\ 0\leq t\leq T_{1} \\ X_{i}\left( 0\right) = x_{i, 0}, \\ _{T_{1}}^{CF}D_{t}^{\alpha }X = F\left( t, X\right) dt+\sigma _{i}X_{i}B_{i}\prime \left( t\right) , \\ {\rm{ if }}\ T_{1}\leq t\leq T_{2} \\ X_{i}\left( T_{1}\right) = x_{i, 0}, 0 < \alpha \leq 1 \end{array} \right. \end{equation}$

(3.1)

where $\sigma _{i}$ and $B_{i}\left(t\right)$ indicate the stochastic constant and Brownian function, respectively. Also, some notation is presented as below.

$\begin{equation} X = \left[ \begin{array}{c} S \\ E \\ I \\ R \end{array} \right] , F\left( t, X\right) = \left[ \begin{array}{c} f_{1}\left( t, X\right) \\ f_{2}\left( t, X\right) \\ f_{3}\left( t, X\right) \\ f_{4}\left( t, X\right) \end{array} \right] = \left[ \begin{array}{c} \mu N-\mu S-\frac{\beta SI}{N} \\ \frac{\beta SI}{N}-\left( \tau +\mu \right) E \\ \tau E-\left( \gamma +\mu \right) I \\ \gamma I-\mu R \end{array} \right] . \end{equation}$

(3.2)

Applying the associated integral, we can have

$\begin{equation} \left\{ \begin{array}{c} X = X_{i}\left( 0\right) +\int_{0}^{t}F\left( \tau , X\right) d\tau {\rm{, \ if }}\ 0\leq t\leq T_{1} \\ X_{i}\left( 0\right) = x_{i, 0}, \\ X = \left\{ \begin{array}{c} X_{i}\left( T_{1}\right) +\frac{1-\alpha }{M\left( a\right) }F\left( t, X\right) +\frac{\alpha }{M\left( a\right) }\int_{0}^{t}F\left( \tau , X\right) d\tau \\ +\frac{1-\alpha }{M\left( a\right) }\sigma _{i}X_{i}B_{i}\prime \left( t\right) +\frac{\alpha }{M\left( a\right) }\int_{0}^{t}\sigma _{i}X_{i}B_{i}\prime \left( \tau \right) d\tau \end{array} \right. {\rm{ }} \\ {\rm{if }}\ T_{1}\leq t\leq T_{2} \\ X_{i}\left( T_{1}\right) = x_{i, 0}, i = 1, .., 4, 0 < \alpha \leq 1 \end{array} \right. . \end{equation}$

(3.3)

We divide $\left[0, T\right]$ in two:

$\begin{equation} 0\leq t_{0}\leq t_{1}\leq ...\leq t_{n_{1}} = T_{1}\leq t_{n_{1}+1}\leq t_{n_{1}+2}\leq ...\leq t_{n_{2}} = T_{2}. \end{equation}$

(3.4)

Interpolating $f_{i}\left(t, X_{i}\right)$ using the Newton polynomial within $\left[t_{n}, t_{n+1}\right]$ yields

$\begin{equation} \left\{ \begin{array}{c} X_{i}^{n_{1}} = X_{i}\left( 0\right) +\sum\nolimits_{j_{1} = 2}^{n_{1}}\left[ \begin{array}{c} \frac{23}{12}F\left( t_{j_{1}}, X_{j_{1}}\right) \\ -\frac{4}{3}F\left( t_{j_{1}-1}, X_{j_{1}-1}\right) \\ +\frac{5}{12}F\left( t_{j_{1}-2}, X_{j_{1}-2}\right) \end{array} \right] \Delta t \\ , 0\leq t\leq T_{1} \\ X_{i}^{n_{2}} = X_{i}\left( T_{1}\right) +\frac{\left( 1-\alpha \right) }{ M\left( \alpha \right) }F\left( t_{n_{2}}, X_{n_{2}}\right) \\ +\frac{\alpha \Delta t}{M\left( \alpha \right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}} \left[ \begin{array}{c} \frac{23}{12}f\left( t_{j_{1}}, X_{j_{1}}\right) \\ +\sigma _{i}X_{j_{1}}\left( B_{i}\left( t_{j_{1}+1}\right) -B_{i}\left( t_{j_{1}}\right) \right) \\ -\frac{4}{3}f\left( t_{j_{1}-1}, X_{j_{1}-1}\right) \\ +\sigma _{i}X_{j_{1}}\left( B_{i}\left( t_{j_{1}}\right) -B_{i}\left( t_{j_{1}-1}\right) \right) \\ +\frac{5}{12}f\left( t_{j_{1}-2}, X_{j_{1}-2}\right) \\ +\sigma _{i}X_{j_{1}}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \end{array} \right] \\ , T_{1}\leq t\leq T_{2} \end{array} \right. . \end{equation}$

(3.5)

The parameters and initial conditions are given as

$\begin{eqnarray} X_{1}\left( 0\right) & = &1000, X_{2}\left( 0\right) = 1, X_{3}\left( 0\right) = 0, X_{4}\left( 0\right) = 0, X_{1}\left( T_{1}\right) = 999.2, \\ X_{2}\left( T_{2}\right) & = &0.4924, X_{1}\left( T_{3}\right) = 0.1847, X_{1}\left( T_{4}\right) = 0.1109, \end{eqnarray}$

(3.6)

and

$\begin{eqnarray} N & = &1000, T_{1} = 50, \beta = 64/30, \mu = 0.5, \tau = 0.3, \gamma = 0.3, \\ \sigma _{1} & = &0.2, \sigma _{2} = 0.12, \sigma _{3} = 0.2, \sigma _{4} = 0.24. \end{eqnarray}$

(3.7)

The simulation of the numerical solution of the model is depicted in Figure 1.

Figure 1. Numerical simulation of SEIR model for

$\alpha = 0.7$ .

DownLoad: Full-Size Img PowerPoint

In this example, some parameters of the model are chosen such that $\frac{ \tau \beta }{\left(\gamma +\mu \right) \left(\tau +\mu \right) } = 1,$ which can be observed from Eq (2.48).

4. Numerical solution of the model with Mittag-Leffler process

In this section, we adopt the numerical scheme for the model where the first part is classical, and the second part is the fractional derivative with the generalized Mittag-Leffler kernel ^[17]. It is important to note that a model with the generalized Mittag-Leffler kernel helps capture processes, with a passage from stretched exponential to power-law with no steady-state.

$\begin{equation} \left\{ \begin{array}{c} X\prime = F\left( t, X\right) {\rm{, \ if }}\ 0\leq t\leq T_{1} \\ X_{i}\left( 0\right) = x_{i, 0}, \\ _{T_{1}}^{ABC}D_{t}^{\alpha }X = F\left( t, X\right) dt+\sigma _{i}X_{i}B_{i}\prime \left( t\right) , \\ {\rm{ if }}\ T_{1}\leq t\leq T_{2} \\ X_{i}\left( T_{1}\right) = x_{i, 0}, 0 < \alpha \leq 1 \end{array} \right. . \end{equation}$

(4.1)

Applying the Atangana-Baleanu integral, we obtain the following equality:

$\begin{equation} \left\{ \begin{array}{c} X = X_{i}\left( 0\right) +\int_{0}^{t}F\left( \tau , X\right) d\tau {\rm{, \ if }}\ 0\leq t\leq T_{1} \\ X_{i}\left( 0\right) = x_{i, 0}, \\ X = \left\{ \begin{array}{c} X_{i}\left( T_{1}\right) +\frac{1-\alpha }{AB\left( a\right) }F\left( t, X\right) +\frac{\alpha }{AB\left( a\right) \Gamma \left( a\right) } \int_{0}^{t}F\left( \tau , X\right) \left( t-\tau \right) ^{\alpha -1}d\tau \\ +\frac{1-\alpha }{AB\left( a\right) }\sigma _{i}X_{i}B_{i}\prime \left( t\right) +\frac{\alpha }{AB\left( a\right) \Gamma \left( a\right) } \int_{0}^{t}\sigma _{i}X_{i}B_{i}\prime \left( \tau \right) \left( t-\tau \right) ^{\alpha -1}d\tau \end{array} \right. , \\ {\rm{ if }}\ T_{1}\leq t\leq T_{2} \\ X_{i}\left( T_{1}\right) = x_{i, 0}, 0 < \alpha \leq 1. \end{array} \right. \end{equation}$

(4.2)

We divide $\left[0, T\right]$ in three:

$\begin{equation} 0\leq t_{0}\leq t_{1}\leq ...\leq t_{n_{1}} = T_{1}\leq t_{n_{1}+1}\leq t_{n_{1}+2}\leq ...\leq t_{n_{2}} = T_{2}. \end{equation}$

(4.3)

Replacing the functions $f_{i}\left(t, X_{i}\right)$ by their Newton polynomials within $\left[t_{n}, t_{n+1}\right],$ the following numerical scheme is obtained:

$\begin{equation} \left\{ \begin{array}{c} X_{i}^{n_{1}} = X_{i}\left( 0\right) +\sum\nolimits_{j_{1} = 2}^{n_{1}}\left[ \begin{array}{c} \frac{23}{12}F\left( t_{j_{1}}, X_{j_{1}}\right) \\ -\frac{4}{3}F\left( t_{j_{1}-1}, X_{j_{1}-1}\right) \\ +\frac{5}{12}F\left( t_{j_{1}-2}, X_{j_{1}-2}\right) \end{array} \right] \Delta t, \\ {\rm{if }}\ 0\leq t\leq T_{1} \\ x^{n_{2}} = x\left( T_{1}\right) +\frac{\left( 1-\alpha \right) \Delta t}{ AB\left( \alpha \right) }f\left( t_{n_{2}}, X_{n_{2}}\right) \\ +\frac{\alpha \Delta t}{AB\left( \alpha \right) \Gamma \left( \alpha +1\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}F\left( t_{j_{2}-2}, X_{j_{2}-2}\right) \Lambda _{1} \\ +\frac{\alpha \Delta t}{AB\left( \alpha \right) \Gamma \left( \alpha +2\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\left[ \begin{array}{c} F\left( t_{j_{2}-1}, X_{j_{2}-1}\right) \\ -F\left( t_{j_{2}-2}, X_{j_{2}-2}\right) \end{array} \right] \Lambda _{2} \\ +\frac{\alpha \Delta t}{2AB\left( \alpha \right) \Gamma \left( \alpha +3\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\left[ \begin{array}{c} F\left( t_{j_{2}}, X_{j_{2}}\right) \\ -2F\left( t_{j_{2}-1}, X_{j_{2}-1}\right) \\ +F\left( t_{j_{2}-2}, X_{j_{2}-2}\right) \end{array} \right] \Lambda _{3} \end{array} \right. \end{equation}$

(4.4)

$\left\{ \begin{array}{c} +\frac{\alpha \Delta t}{AB\left( \alpha \right) \Gamma \left( \alpha +1\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\sigma _{i}X_{j_{1}}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \Lambda _{1} \\ +\frac{\alpha \Delta t}{AB\left( \alpha \right) \Gamma \left( \alpha +2\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\left[ \begin{array}{c} \sigma _{i}X_{j_{2}-1}\left( B_{i}\left( t_{j_{1}}\right) -B_{i}\left( t_{j_{1}-1}\right) \right) \\ -\sigma _{i}X_{j_{2}-2}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \end{array} \right] \Lambda _{2} \\ +\frac{\alpha \Delta t}{2AB\left( \alpha \right) \Gamma \left( \alpha +3\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\left[ \begin{array}{c} \sigma _{i}X_{j_{2}}\left( B_{i}\left( t_{j_{1}+1}\right) -B_{i}\left( t_{j_{1}}\right) \right) \\ -2\sigma _{i}X_{j_{2}-1}\left( B_{i}\left( t_{j_{1}}\right) -B_{i}\left( t_{j_{1}-1}\right) \right) \\ +\sigma _{i}X_{j_{2}-2}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \end{array} \right] \Lambda _{3}\end{array} \right. , \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ {\rm{if }}\ T_{1}\leq t\leq T_{2} \notag$

(4.5)

where

$\begin{eqnarray} \Lambda _{1} & = &\left( n_{2}-j_{2}+1\right) ^{\alpha }-\left( n_{2}-j_{2}\right) ^{\alpha }\\ \Lambda _{2} & = &\left\{ \begin{array}{c} \left( n_{2}-j_{2}+1\right) ^{\alpha }\left( n_{2}-j_{2}+3+2\alpha \right) \\ -\left( n_{2}-j_{2}\right) ^{\alpha }\left( n_{2}-j_{2}+3+3\alpha \right) \end{array} \right\} , \\ \Lambda _{3} & = &\left\{ \begin{array}{c} \left( n_{2}-j_{2}+1\right) ^{\alpha }\left[ \begin{array}{c} 2\left( n_{2}-j_{2}\right) ^{2}+\left( 3\alpha +10\right) \left( n_{2}-j_{2}\right) \\ +2\alpha ^{2}+9\alpha +12 \end{array} \right] \\ -\left( n_{2}-j_{2}\right) ^{\alpha }\left[ \begin{array}{c} 2\left( n_{2}-j_{2}\right) ^{2}+\left( 5\alpha +10\right) \left( n_{2}-j_{2}\right) \\ +6\alpha ^{2}+18\alpha +12 \end{array} \right] \end{array} \right\} . \notag \end{eqnarray}$

(4.6)

For this model, the parameters and initial conditions are considered as follows:

$\begin{eqnarray} X_{1}\left( 0\right) & = &1000, X_{2}\left( 0\right) = 1, X_{3}\left( 0\right) = 0, X_{4}\left( 0\right) = 0, X_{1}\left( T_{1}\right) = 1000, \\ X_{2}\left( T_{1}\right) & = &4.642\times 10^{-5}, X_{1}\left( T_{1}\right) = 2.315\times 10^{-5}, X_{1}\left( T_{1}\right) = 0.0069 \end{eqnarray}$

(4.7)

and

$\begin{eqnarray} N & = &1000, T_{1} = 50, \beta = 0.001, \mu = 0.1, \tau = 0.1, \gamma = 0.3, \\ \sigma _{1} & = &0.0000025, \sigma _{2} = 0.011, \sigma _{3} = 0.0091, \sigma _{4} = 0.017. \end{eqnarray}$

(4.8)

The simulation of the numerical solution of the model is shown in Figure 2.

Figure 2. Numerical simulation of SEIR model for

$\alpha = 0.6$ .

DownLoad: Full-Size Img PowerPoint

In this example, the parameters of the model are chosen such that $\frac{ \tau \beta }{\left(\gamma +\mu \right) \left(\tau +\mu \right) } < 1,$ which is seen from Eq (2.48).

5. Numerical solution of the model with power-law process

In this section, we consider the following model:

$\begin{equation} \left\{ \begin{array}{c} X\prime = F\left( t, X\right) {\rm{, \ if }}\ 0\leq t\leq T_{1} \\ X_{i}\left( 0\right) = x_{i, 0}, \\ _{T_{1}}^{C}D_{t}^{\alpha }X = F\left( t, X\right) dt+\sigma _{i}X_{i}B_{i}\prime \left( t\right) , \\ {\rm{ if }}\ T_{1}\leq t\leq T_{2} \\ X_{i}\left( T_{1}\right) = x_{i, 0}, 0 < \alpha \leq 1 \end{array} \right. \end{equation}$

(5.1)

where the first part is classical, and second part is the fractional derivative with the generalized Mittag-Leffler kernel. Applying the associated integral, we obtain the following equality:

$\begin{equation} \left\{ \begin{array}{c} X = X_{i}\left( 0\right) +\int_{0}^{t}F\left( \tau , X\right) d\tau {\rm{, \ if }}\ 0\leq t\leq T_{1} \\ X_{i}\left( 0\right) = x_{i, 0}, \\ X = \left\{ \begin{array}{c} X_{i}\left( T_{1}\right) +\frac{1}{\Gamma \left( a\right) } \int_{0}^{t}F\left( \tau , X\right) \left( t-\tau \right) ^{\alpha -1}d\tau \\ +\frac{1}{\Gamma \left( a\right) }\int_{0}^{t}\sigma _{i}X_{i}B_{i}\prime \left( \tau \right) \left( t-\tau \right) ^{\alpha -1}d\tau \end{array} \right. , \\ {\rm{ if }}\ T_{1}\leq t\leq T_{2} \\ X_{i}\left( T_{1}\right) = x_{i, 0}, 0 < \alpha \leq 1 \end{array} \right. . \end{equation}$

(5.2)

We divide $\left[0, T\right]$ in three:

$\begin{equation} 0\leq t_{0}\leq t_{1}\leq ...\leq t_{n_{1}} = T_{1}\leq t_{n_{1}+1}\leq t_{n_{1}+2}\leq ...\leq t_{n_{2}} = T_{2}. \end{equation}$

(5.3)

Replacing the functions $f_{i}\left(t, X_{i}\right)$ by their Newton polynomials within $\left[t_{n}, t_{n+1}\right],$ the associated model can be solved by the following algorithm:

$\begin{equation} \left\{ \begin{array}{c} X_{i}^{n_{1}} = X_{i}\left( 0\right) +\sum\nolimits_{j_{1} = 2}^{n_{1}}\left[ \begin{array}{c} \frac{23}{12}F\left( t_{j_{1}}, X_{j_{1}}\right) \\ -\frac{4}{3}F\left( t_{j_{1}-1}, X_{j_{1}-1}\right) \\ +\frac{5}{12}F\left( t_{j_{1}-2}, X_{j_{1}-2}\right) \end{array} \right] \Delta t, \\ {\rm{if }}\ 0\leq t\leq T_{1} \\ X_{i}^{n_{2}} = X_{i}\left( T_{1}\right) +\frac{\Delta t}{\Gamma \left( \alpha +1\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}F\left( t_{j_{2}-2}, X_{j_{2}-2}\right) \Lambda _{1} \\ +\frac{\Delta t}{\Gamma \left( \alpha +2\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}} \left[ \begin{array}{c} F\left( t_{j_{2}-1}, X_{j_{2}-1}\right) \\ -F\left( t_{j_{2}-2}, X_{j_{2}-2}\right) \end{array} \right] \Lambda _{2} \\ +\frac{\Delta t}{2\Gamma \left( \alpha +3\right) } \sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\left[ \begin{array}{c} F\left( t_{j_{2}}, X_{j_{2}}\right) \\ -2F\left( t_{j_{2}-1}, X_{j_{2}-1}\right) \\ +F\left( t_{j_{2}-2}, X_{j_{2}-2}\right) \end{array} \right] \Lambda _{3} \end{array} \right. \end{equation}$

(5.4)

$\begin{gather} \left\{ \begin{array}{c} +\frac{\Delta t}{\Gamma \left( \alpha +1\right) } \sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\sigma _{i}X_{j_{2}}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \Lambda _{1} \\ +\frac{\Delta t}{\Gamma \left( \alpha +2\right) }\sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}} \left[ \begin{array}{c} \sigma _{i}X_{j_{2}-1}\left( B_{i}\left( t_{j_{1}}\right) -B_{i}\left( t_{j_{1}-1}\right) \right) \\ -\sigma _{i}X_{j_{2}-2}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \end{array} \right] \Lambda _{2} \\ +\frac{\Delta t}{2\Gamma \left( \alpha +3\right) } \sum\nolimits_{j_{2} = n_{1}+3}^{n_{2}}\left[ \begin{array}{c} \sigma _{i}X_{j_{2}}\left( B_{i}\left( t_{j_{1}+1}\right) -B_{i}\left( t_{j_{1}}\right) \right) \\ -2\sigma _{i}X_{j_{2}-1}\left( B_{i}\left( t_{j_{1}}\right) -B_{i}\left( t_{j_{1}-1}\right) \right) \\ +\sigma _{i}X_{j_{2}-2}\left( B_{i}\left( t_{j_{1}-1}\right) -B_{i}\left( t_{j_{1}-2}\right) \right) \end{array} \right] \Lambda _{3} \end{array} \right. , \end{gather} \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ {\rm{if }}\ T_{1}\leq t\leq T_{2}.$

(5.5)

The parameters and initial conditions are given as

$\begin{eqnarray} X_{1}\left( 0\right) & = &1000, X_{2}\left( 0\right) = 1, X_{3}\left( 0\right) = X_{4}\left( 0\right) = 0, X_{1}\left( T_{1}\right) = 28.03, \\ X_{2}\left( T_{1}\right) & = &78.01, X_{3}\left( T_{1}\right) = 593.1, X_{4}\left( T_{1}\right) = 301.4, \end{eqnarray}$

(5.6)

and

$\begin{eqnarray} N & = &1000, T_{1} = 50, \beta = 0.6, \mu = 0.01, \tau = 0.2, \gamma = 0.03, \\ \sigma _{1} & = &0.000005, \sigma _{2} = 0.008, \sigma _{3} = 0.0071, \sigma _{4} = 0.017. \end{eqnarray}$

(5.7)

The simulation of the numerical solution of the model is presented in Figure 3.

Figure 3. Numerical simulation of SEIR model for

$\alpha = 0.75$ .

DownLoad: Full-Size Img PowerPoint

In this example, the parameters of the model are chosen such that $\frac{ \tau \beta }{\left(\gamma +\mu \right) \left(\tau +\mu \right) } > 1,$ which is seen from Eq (2.48).

6. Conclusions

In classical calculus, analysis of the first derivative of a given function helps one to understand the change in such a function. In addition, this formula calculates the velocity of a moving object. However, analysis using the first derivative does not always a priori provide a clear indication of the variation of a function. Therefore, an analysis of the second derivative is required. The study of the second derivative includes information on inflection points, local maximums, and minimums. These elementary analyses can better show spread patterns in epidemiological modeling. A new concept was recently introduced and named strength number, which is obtained by considering the second derivative of the non-linear part of a given infectious disease model. Then, the next-generation matrix methodology is used to obtain the strength number. It was argued that such numbers could help detect waves or instability in a model. This work applies such a concept, together with a second derivative analysis, in a simple SEIR problem. The obtained strength number was negative, with one equilibrium point for the model with second derivatives, indicating that such a model could have only one wave and then die out. We have used piecewise differential operators to add stochastic behavior to the model.

Conflict of interest

The authors declare no conflict of interest.

References

[1]	A. Atangana, D. Baleanu, New fractional derivatives with non-local and non-singular kernel: Theory and application to heat transfer model, Therm. Sci., 20 (2016), 763–769. https://doi.org/10.2298/TSCI160111018A doi: 10.2298/TSCI160111018A
[2]	M. Caputo, M. Fabrizio, A new definition of fractional derivative without singular kernel, Progress in Fractional Differentiation and Applications, 1 (2015), 73–85. https://doi.org/10.12785/pfda/010201 doi: 10.12785/pfda/010201
[3]	M. Caputo, Linear model of dissipation whose Q is almost frequency independent-Ⅱ, Geophys. J. Int., 13 (1967), 529–539. https://doi.org/10.1111/j.1365-246X.1967.tb02303.x doi: 10.1111/j.1365-246X.1967.tb02303.x
[4]	A. Atangana, Mathematical model of survival of fractional calculus, critics and their impact: How singular is our world?, Adv. Differ. Equ., 2021 (2021), 403. https://doi.org/10.1186/s13662-021-03494-7 doi: 10.1186/s13662-021-03494-7
[5]	C. Ji, D. Jiang, N. Shi, The behavior of an SIR epidemic model with stochastic perturbation, Stoch. Anal. Appl., 30 (2012), 755–773. https://doi.org/10.1080/07362994.2012.684319 doi: 10.1080/07362994.2012.684319
[6]	Y. Zhao, D. Jiang, The threshold of a stochastic SIRS epidemic model with saturated incidence, Appl. Math. Lett., 34 (2014), 90–93. https://doi.org/10.1016/j.aml.2013.11.002 doi: 10.1016/j.aml.2013.11.002
[7]	M. Bachar, M. A. Khamsi, M. Bounkhel, A mathematical model for the spread of Covid-19 and control mechanisms in Saudi Arabia, Adv. Differ. Equ., 2021 (2021), 253. https://doi.org/10.1186/s13662-021-03410-z doi: 10.1186/s13662-021-03410-z
[8]	Z. B. Dieudonné, Mathematical model for the mitigation of the economic effects of the Covid-19 in the Democratic Republic of the Congo, PLoS ONE, 16 (2021), e0250775. https://doi.org/10.1371/journal.pone.0250775 doi: 10.1371/journal.pone.0250775
[9]	M. Tomochi, M. Kono, A mathematical model for COVID-19 pandemic-SIIR model: Effects of asymptomatic individuals, J. Gen. Fam. Med., 22 (2020), 5–14. https://doi.org/10.1002/jgf2.382 doi: 10.1002/jgf2.382
[10]	P. Sahoo, H. S. Mondal, Z. Hammouch, T. Abdeljawad, D. Mishra, M. Reza, On the necessity of proper quarantine without lock down for 2019-nCoV in the absence of vaccine, Results Phys., 25 (2021), 104063. https://doi.org/10.1016/j.rinp.2021.104063 doi: 10.1016/j.rinp.2021.104063
[11]	A. Babaei, H. Jafari, S. Banihashemi, M. Ahmadi, Mathematical analysis of a stochastic model for spread of Coronavirus, Chaos Soliton. Fract., 145 (2021), 110788. https://doi.org/10.1016/j.chaos.2021.110788 doi: 10.1016/j.chaos.2021.110788
[12]	I. Ahmed, E. F. D. Goufo, A. Yusuf, P. Kumam, P. Chaipanya, K. Nonlaopon, An epidemic prediction from analysis of a combined HIV-COVID-19 co-infection model via ABC-fractional operator, Alex. Eng. J., 60 (2021), 2979–2995. https://doi.org/10.1016/j.aej.2021.01.041 doi: 10.1016/j.aej.2021.01.041
[13]	J. M. Carcione, J. E. Santos, C. Bagaini, J. Ba, A simulation of a COVID-19 epidemic based on a deterministic SEIR model, Front. Public Health, 8 (2020), 230. https://doi.org/10.3389/fpubh.2020.00230 doi: 10.3389/fpubh.2020.00230
[14]	P. van den Driessche, J. Watmough, Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission, Math. Biosci., 180 (2002), 29–48. https://doi.org/10.1016/S0025-5564(02)00108-6 doi: 10.1016/S0025-5564(02)00108-6
[15]	J. P. La Salle, The stability of dynamical systems, SIAM Press, 1976. https://doi.org/10.1137/1.9781611970432
[16]	A. Atangana, S. İğret Araz, New numerical scheme with Newton polynomial: theory, methods and applications, Academic Press, 2021. https://doi.org/10.1016/B978-0-12-775850-3.50017-0
[17]	A. Atangana, S. İğret Araz, New concept in calculus: Piecewise differential and integral operators, Chaos Soliton. Fract., 145 (2021), 110638. https://doi.org/10.1016/j.chaos.2020.110638 doi: 10.1016/j.chaos.2020.110638

This article has been cited by:

1.	Victoria May P. Mendoza, Renier Mendoza, Youngsuk Ko, Jongmin Lee, Eunok Jung, Managing bed capacity and timing of interventions: a COVID-19 model considering behavior and underreporting, 2022, 8, 2473-6988, 2201, 10.3934/math.2023114
2.	Mehmet Akif Cetin, Seda Igret Araz, Prediction of COVID-19 spread with models in different patterns: A case study of Russia, 2024, 22, 2391-5471, 10.1515/phys-2024-0009
3.	Muhammad Umer Saleem, Muhammad Farman, Rabia Sarwar, Parvaiz Ahmad Naik, Perwasha Abbass, Evren Hincal, Zhengxin Huang, Modeling and analysis of a carbon capturing system in forest plantations engineering with Mittag–Leffler positive invariant and global Mittag–Leffler properties, 2025, 11, 2363-6203, 10.1007/s40808-024-02181-2
4.	Muhammad Farman, Rabia Sarwar, Evren Hincal, Dumitru Baleanu, Kottakkaran Sooppy Nisar, Aceng Sambas, Muhammad Manan Akram, Dynamical Evaluation of the Mittag-Leffler Properties-Based Marine Ecosystem Model, 2025, 2509-9426, 10.1007/s41748-025-00617-y

Reader Comments

Your name:*

Email:*
© 2022 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)

通讯作者: 陈斌, bchen63@163.com

1.
沈阳化工大学材料科学与工程学院沈阳 110142

AIMS Mathematics

1.8 3.1

Metrics

Article views(2094) PDF downloads(119) Cited by(4)

Preview PDF

Download XML

Export Citation

Article outline

Show full outline

AIMS Mathematics

Advanced analysis in epidemiological modeling: detection of waves

Related Papers:

Abstract

1. Introduction

2. A SEIR model

2.1. Positiveness and boundedness of solutions

2.2. Analysis of equilibrium points

2.3. Reproduction number

2.4. Strength number

2.5. First derivative of Lyapunov

2.6. Second derivative of Lyapunov

2.7. Equilibrium points for second order

3. Numerical solution of the model with exponential decay process

4. Numerical solution of the model with Mittag-Leffler process

5. Numerical solution of the model with power-law process

6. Conclusions

Conflict of interest

References

This article has been cited by:

Reader Comments

通讯作者: 陈斌, bchen63@163.com

Metrics

Figures and Tables

Other Articles By Authors

Related pages

Tools

Export File

Citation

Format

Content

Catalog

Abstract

1. Introduction

2. A SEIR model

2.1. Positiveness and boundedness of solutions

2.2. Analysis of equilibrium points

2.3. Reproduction number

2.4. Strength number

2.5. First derivative of Lyapunov

2.6. Second derivative of Lyapunov

2.7. Equilibrium points for second order

3. Numerical solution of the model with exponential decay process

4. Numerical solution of the model with Mittag-Leffler process

5. Numerical solution of the model with power-law process

6. Conclusions

Conflict of interest

References