ABSTRACT
We present two slightly different constructions of a SIR model in which both the time taken to remove the individual from the infectious compartment and the infectivity have a memory according to Mittag-Leffler distributions. The second construction clearly points out where the proposed generalizations are acting, starting from the classic SIR model. Using integrodifferential techniques, we state and demonstrate novel analytically results on positivity, monotonicity in limiting case, and equilibrium points. The results are also verified numerically.
Keywords: Fractional SIR model; Mittag-Leffler functions; positivity; integrodifferential equations; equilibrium
1 INTRODUCTION
Fractional differential equations are extensively applied in memory-dependent mathematical models. In fact, the Fractional Calculus includes the possibility of considering the memory of the studied phenomenon, that is, the dependence on previous stages, once its operators are non local. Most often, a fractional model is designed in an ad hoc conception: a classical model is generalized by allowing the integer order of the derivative to be arbitrary. However, this practice can lead to some physical misinterpretation, and previous works study manners to construct fractional models with biological meaning 1), (2), (3), (17.
This work aims to complete previous discussions of fractional SIR models with Mittag-Leffler memory, presenting the detailed proofs of our results in 16. To this aim, in Subsection 2.1, we present an overview of the work 3, a recent and unusual SIR model constructed by Continuous Time Random Walks (CTRW). Both the time of removal of the individual from the infectious compartment and the infectivity are age-dependent and follow a Mittag-Leffler distribution, causing the Riemann-Liouville fractional derivative to arise throughout Laplace Transform techniques. In Subsection 2.2, we show how to obtain from the theory presented in 15 the same fractional model constructed with CTRW. It shows clearly where the generalizations proposed are acting, from the viewpoint of the classical SIR model. Then, we propose and demonstrate new results about positivity, monotonicity (Section 3), and equilibrium points (Section 4).
MATLAB simulations with codes proposed in 17 and other references verify the theory numerically. We use basic definitions and results of Fractional Calculus from classic references 6), (12), (14), (20), (22), (23. We denote I 1−α the Riemann-Liouville integral and D 1−α the Riemann-Liouville derivative operators with initial point zero.
2 THE MODEL
The construction proposed in 3 is studied in details in 17), (19 and correlated references. Here, we present a brief overview, for completeness.
2.1 An overview
The main idea is that, if there are S(t) susceptible people at time t in a population of size N, then the expected number of new infections per infected individual is given by σ(t,t' in a step time ∆T . Here, σ (t, t') is the transmission function of an individual first infected at time t' < t.
The survival function Φ(t, t') gives the probability that a person first infected at time t' < t stays infected at time t. Considering that there is no disease before t = 0, the flow of individuals into infectious compartment I in time t is recursively given by
where i 0 is the initial condition: the quantity of infectious that emerged in time 0.
The transmission function σ (t, t') is dependent of both the present t, due, e.g., to cultural or containment measures, and the age of infection t - t'. The dependence of t is expressed by an extrinsic infectivity ω, independent of the individual. On the other hand, there is an intrinsic infectivity ρ dependent of age. So,
For the survival function, it is assumed that natural death and recuperation/death by disease of an infected individual are independent. Then,
where ϕ(t - t') is the probability that an individual first infected in t' is still with the disease at time t. Furthermore, θ(t, t') is the probability that an individual first infected in t' has not yet died a natural death by time t, being given by
where γ is the vital rate.
Individuals in compartment I at time t must have entered this compartment at some previous time and remained there until t. Therefore, we can express the number of infected individuals as follows:
Writing ψ(t) = −dϕ(t)/dt, we derive Eq. (2.5) through Leibniz’s Rule, obtaining
The dependency on q +(I, t') is removed by defining memory kernels for infectivity and recovery. Henceforth we consider i 0 = 1 for simplicity. Eq. (2.5) can be rewritten as
Using the Laplace transform in the first integral of (2.6) we can write
where it is defined the infectivity memory kernel
In the same way, based in the second integral of Eq. (2.6),
where it is defined the recovery memory kernel by
Fractional derivatives are incorporated into the model by choosing ψ(t) with potential law and ρ(t) related to the choice of ψ(t):
for 0 < α ≤ 1, where τ is a scale parameter. The corresponding survival function is
We can calculate the Laplace transform of the removal memory kernel defined by the Mittag-Leffler function, obtaining:
The Laplace transform of the kernel implies that a convolution with the recovery memory kernel can be written as
Another fractional derivative is incorporated into the infectivity memory kernel by considering
Since ρ(t) ≥ 0 is necessary, we take 0 < α ≤ β ≤ 1. Using Eq. (2.16), we obtain the Laplace transform of the infectiousness kernel:
Substituting Eqs. (2.15) and (2.17) into the master equation Eq.(2.6) and rewriting, we obtain a fractional SIR model:
Once there is no disease before time 0, the initial condition must be in the form S(0) = N − i 0 , I(0) = i 0 , R(0) = 0.
2.2 Another approach for constructing a fractional SIR model
The function ψ(t) can be related with the continuous random variable X that gives the individual’s waiting time in infectious compartment despite the vital dynamics. To see this, note that the cumulative distribution of X is F(t) = P(X ≤ t) = 1 − ϕ(t). Thus, the probability density function of X is ψ(t) = −dϕ(t)/dt. The probability distribution F is a Mittag-Leffler distribution . If α < 1, the expectation is infinite. This can represent, for instance, cases in which long waiting times are not really improbable. However, if α = 1, it represents an exponential distribution. For waiting times exponentially distributed, the first moment of the random variable X exists, and τ is exactly the mean waiting time in infectious compartment despite the vital dynamics. It is worth to note that, if γ(t) ≡ γ > 0, the effective mean waiting time in infectious compartment, described by both the effects of recovery and natural death, is finite and given by γ α-1 /(γ α + τ -α ) for any α ∈ (0, 1] 15. In this Section, we show how to obtain from the presented theory in 15 the same fractional model constructed with CTRW in the previous Section.
The compartmental SIR (Susceptible-Infected-Removed) model was introduced by Kermack and McKendrick in 1927 11. With constants rates of vital dynamics (γ), infectivity (σ) and mean removal time (τ), the SIR model can be written as
From the second equation, we can write
We generalize e −t/τ , which indicates an exponentially distributed waiting time in the infectious compartment, with a Mittag-Leffler probability .
This generalization constructs a fractional SIR model with general remotion time as in 2. However, as seen in previous Subsection, in 3 there is also a generalization in the transmission function: in fact, the infectivity σ is seen as a function with an extrinsic term ω(t), depending on the present time (for instance, cultural and political measures), and an intrinsic term ρ(t - t'), depending on the age of infection:
If there is no memory, i.e., if the force of infectivity is the same along the days of infection, so ρ is a constant.
Note that a people that becomes infectious in a time t' is still infectious in t if:
-
It is alive (probability e −γ(t−t′) ), and
-
It is not removed (probability ).
In this case, its infectivity is ω(t)ρ(t - t'). So, the flux of new infectious is given recursively by
Aiming to use the Laplace Transform of Mittag-Leffler functions in Eq.(2.24) and considering that the intrinsic infectivity is positive and monotonically decreasing with time, it is proposed
with 0 < α ≤ β ≤ 1. Then,
Note that
where
Note that Y indicates the summation of people that entered in I compartment in each 0 < t' ≤ t and are still alive and infectious at t. So,
This implies that
Follows then from (2.26) and (2.35) that
Finally, deriving (2.33) by Leibniz rule and using
we get
This is the master equation of the fractional SIR model (2.18)-(2.20) with constant γ(t) ≡ γ:
In next Section, we demonstrate the non-negativity of the solution. Furthermore, we show that, in the limiting case where the vital dynamics is null, the compartments S and R are monotone.
3 MONOTONICITY IN THE LIMITING CASE AND NON-NEGATIVITY
The solution of the fractional SIR model constructed in previous Section is expected to be non-negative, since the compartments reproduce groups of people. In fact, we demonstrate the following result:
Proposition 1.
Let γ(t) ≡ γ ≥ 0. Admitting that the solution of the system (2.18)-(2.20) exists uniquely with initial conditions S(0) ≥ 0, I(0) = i 0 > 0 and R(0) = N − S(0) − I(0) ≥ 0, and is continuous for t > 0, then:
-
and are non-negative functions;
-
S(t), I(t), R(t) ≤ N for all t.
Proof. First, we prove the positivity of S(t) for both γ > 0 and γ = 0. Initially, consider γ > 0. Suppose absurdly that the set {t > 0; S(t) < 0} is not empty. Let T > 0 be the infimum of this set. Once the solution is continuous, we have S(T) = 0. From Eq. (2.18), we have S'(T) = γN > 0. Thus, there is an interval (T, T + ε) such that S'(T ) > 0 and S is increasing in this interval, absurd, since we define T the infimum of the set {t > 0; S(t) < 0}, and so there is a sequence in {t > 0; S(t) < 0} that converges to T from the right. Therefore, S(t) > 0 for all t.
On the other hand, if γ = 0, we can write
where
Once S(0) ≥ 0, Eq.(3.4) implies that S(t) ≥ 0 for all t.
Given the positivity of S(t), as σ and Φ are positive functions, Eq. (2.1) or Eq. (2.26) gives q +(I, t) ≥ 0 for all t. From Eq. (2.5) or Eq. (2.33) it follows directly that I(t) > 0 for all t. Now, applying the inverse Laplace Transform to Eq. (2.8), we have
It follows from the positivity of ρ, ϕ, θ and q +(I, t) that
for all t. Therefore, from Eq. (2.17), we have , for all t. The same can be concluded for .
Finally, suppose absurdly that the set {t > 0; R(t) < 0} is not empty. Let T ≥ 0 be the infimum of this set. Once the solution is continuous and R(0) ≥ 0, we must have R(T) = 0. Then,
Thus, there is an interval (T, T + ε) such that R'(t) > 0 and R is increasing in this interval, absurd, since there is a sequence in {t > 0; R(t) < 0} that converges to T from the right. Therefore, R(t) ≥ 0 for all t.
Since S(t) + I(t) + R(t) = N, it follows from S(t), I(t), R(t) ≥ 0 that S(t), I(t), R(t) ≤ N, completing the proof. □
For the limit case in which γ = 0, we also have the following result:
Proposition 2.Considering the limit system (2.18)-(2.20) for which γ(t) ≡ 0, the functions S(t) and R(t) are monotone decreasing and increasing, respectively.
Proof. It follows directly from S(t), and be non-negative that, if γ(t) ≡ 0, then S'(t) ≤ 0 and R'(t) ≥ 0. Thus S and R are monotonically decreasing and increasing, respectively. □
We note that this condition has good epidemiological agreement, but we cannot guarantee this for ad hoc SIR fractional models, where the integer order derivative is replaced by the Caputo derivative:
In fact, Figure 1 indicates the solution of the system (3.9). We illustrate the non-monotonicity of the compartments S and R in the model without vital dynamics. This problem is not solved by balancing the units of the parameters initially considered. With the same orders, the corrections would transform them numerically into other constants, but the same in each compartment. On the other hand, in Figure 2 the system (2.18)-(2.20) is numerically solved. We observe that monotonicity remains in the model with arbitrary orders, as observed in Proposition 2, even if α ≠ β. Thus, the oscillations that can appear in this model are due exclusively to vital dynamics. We use a population N 0 = 106, initial conditions S(0) = N 0 − 1, I(0) = 1 and R(0) = 0 and dt = 0.1 for both simulations. Our previous works 13), (18 present a deep discussion about the ad hoc fractionalization (3.9).
Next Section deals with equilibrium points. In particular, we study its relationship with the basic reproduction number and approaches to stability analysis.
4 EQUILIBRIUM POINTS AND STABILITY
Classical epidemiological ODEs’ model as (2.21) are time invariant 7. Their local behavior allows the parameters to define the disease dynamics independently of the time in which one begins the modeling. However, this is not valid in model (2.18)-(2.20). Different pasts modify the future dynamics. Its non-local behavior also implies unexpected behavior of the reproduction number ℜ(t). The reproduction number in t can be defined as the expected number of individuals infected by a person first infected in time t. The basic reproduction number is this number for t = 0, i.e., the average number of secondary infections that occur when an infectious individual is introduced into a completely susceptible population 9. The definition of the reproduction numbers in the proposed model is non intuitive. We mention that, in 3, it is proposed an integral construction for the basic reproduction number. In 19, we extend the reproduction number proposal to any time t and also propose the S-variable reproduction numbers. These discussions have important implications in equilibrium analysis and peak conditions. For instance, the peak does not occur when ℜ(t) = 1, nor ℜS (t) = 1 19.
For this work proposal, we only use the basic reproduction number 3:
Returning to the equilibrium calculation, we consider γ(t) ≡ γ constant, so θ(t, 0) = e −γt . Taking the limit when t → ∞, we calculate the limits of the form according to 19.
Assuming , we obtain a disease-free state,
and, if ω* > 0, an endemic state:
Particularly, when ω(t) ≡ ω, this result is the same as that obtained in 3. Note that the viability of the endemic equilibrium (4.3) requires
which is related to the basic reproduction number (4.1).
In fact, if ω(t) ≡ ω, the criterion for the viability given in Eq. (4.4) can be rewritten as ℜ0 > 1. Furthermore, the value S ∗ of the endemic state given in (4.3) is of the form S ∗ = N/ℜ0. Thus, ℜ0 has an essential relationship with the final size of the infection.
Until now, we have proven that if there are asymptotically stable equilibria for the case γ > 0, then they are given by Eq.(4.2)-(4.3). The limiting case in which γ = 0 is studied in 16:
Theorem 1.If ω(t) ≡ ω, γ(t) ≡ 0, i 0 > 0 and α = β in the system (2.18)-(2.20), then the solution asymptotically approaches equilibrium (S, I, R)∞ , where
and
W0represents the main branch of Lambert’s function W5, and S0 = S(0).
Proof. See 16. □
For γ > 0, the disease-free state is expected to be an asymptotically stable equilibrium when ω* < (τγ)β−α + (τγ)β , while the endemic state is expected to be asymptotically stable if ω* > (τγ)β−α + (τγ)β . These hypothesis are formulated, for constant ω(t), in 3.
To develop our results, based on 10, we need the integral forms of the infectious and recovered compartments:
where F(t) = 1 − ϕ(t).
The following results are demonstrated:
Theorem 2.If ω(t) is bounded with, γ(t) ≡ γ and β = 1 in the system (2.18)-(2.20), the disease-free equilibrium of Eq. (4.2) is globally asymptotically stable if ω* < (τγ)1−α + (τγ).
Proof. We define the constant A such that
Let and suppose absurdly that J > 0. Thus, it is possible to choose ε small enough such that 2ε + A(1 + ε)(J + ε) < J. Indeed, defining the continuously increasing function E(t) = 2t + A(1 + t)(J + t), we have E(0) = AJ < J and , so there is ε > 0 such that E(ε J ) = J. Therefore, we just need to choose ε ∈ (0, ε J ).
We have . Furthermore, we know that and . Thus, choosing ε, we can choose a time t 1 > 0 sufficiently large such that Φ(t, 0)i 0 < ε/2, I(t) < J + ε and ω(t) < (1 + ε)ω* for t > t 1. Choosing t 1, since lim θ(t, 0) = 0, we can choose t 2 large enough such that Nt 1 θ(t, 0)ω m /τ < ε for t > t 2, where ω m is the maximum value reached by ω(t). So, for t > t 1 + t 2, we have t − t 1 > t 2 and
Now, we remember that S(t), I(t) ≤ N. In the first integral we use Φ(t, 0) = ϕ(t) θ(t, 0) ≤ θ(t, 0). In the second integral, we use that I(t) < J + ε and ω(t) < (1 + ε)ω* for t > t 1, obtaining
Finally, as t − t 1 > t 2, we have θ(t, t 1) < θ(t 2). Furthermore, by the definition of A, we have . Therefore, by choosing ε, t 1 and t 2, it follows that, for all t > t 1 + t 2, we have
contradiction. Therefore, J = 0 and .
Remembering that, under the conditions of the statement, we have
then implies that . Therefore, . □
Theorem 3.If ω(t) ≡ ω, γ(t) ≡ γ and β = 1 in the system (2.18)-(2.20), the endemic equilibrium of Eq. (4.3) is locally asymptotically stable whenever it is feasible, that is, when ω > (τγ)1−α + (τγ).
Proof. We have . We remember that if β = 1 and ω is constant, then σ(t) ≡ ω/τ. Thus,
Now, we observe that
i.e,
Returning to Eq. (4.3), we note that R ∗ = (τγ)−α I ∗. So, we write
We do V = I − I ∗ and W = R − R ∗ to translate the equilibrium to the origin. Then, by Eqs. (4.7), (4.8), (4.14) and (4.17) and observing that SI − S ∗ I ∗ = S ∗ V − I ∗(V + W) − V (V + W), we get
on what
The nonlinear Volterra integral system of Eq. (4.18) can be written in matrix form as:
on what
The characteristic equation of the linearization of Eq. (4.20) is
where J is the Jacobian of G evaluated at 0.
The stability analysis of the origin uses a result cited in [10]: If the solutions of Eq. (4.20) exist in [0, ∞) and are bounded, F(t) ∈ C[0, ∞), F(t) → 0 when t → ∞, A(t) ∈ L 1[0, ∞), G(X ) ∈ C 1(ℝ2), G(0) = 0, J is non-singular and Eq. (4.23) does not have roots with non-negative real parts, then the origin is a locally asymptotically stable equilibrium.
Most of the conditions are straightforward. In particular, F(t) → 0 due to the existence of the Laplace Transform of the Mittag-Leffler function. We need to prove that there are no characteristic roots of Eq. (4.23) with a non-negative real part. First, we calculate J, obtaining
Thus, the characteristic equation is given by
We observe that the integrals have the form of a Laplace transform and, in 10, is requested the hypothesis that the probability of permanence ϕ is dominated by an exponential decay, in order to guarantee the existence of the transform. The Mittag-Leffler function is not dominated by exponential decay, but the Laplace transform converges if Re(λ + γ) > 0, and |λ + γ| > 1/τ8.
Eq. (4.25) can be written as
Which resumes to
Finally, since Φ(t, 0) = θ(t, 0)ϕ(t) = e −γt (1 − F(t)), we have
The condition ω > (τγ)1−α + (τγ) is equivalent to
Supposing absurdly that Re(λ) ≥ 0, then
Thus,
and, to fulfill Eq. (4.29), it will be necessary to have
But, see that, if Re(λ) > −γ, we have
Since λ = a + bi with a ≥ 0, we have
following that
contradiction with (4.33). Therefore, the equilibrium (S ∗ , I ∗ , R ∗) is locally asymptotically stable. □
5 DISCUSSION AND CONCLUSION
We discuss the SIR model from 3 with the aim of expanding on previous results and sharing the findings. This model is of interest to us because it is based on a meaningful construction that naturally involves derivatives of arbitrary order. We also demonstrate how the same fractional model can be derived from the theory presented in 15. This approach helps to clarify the role of the proposed generalizations in the context of the classical SIR model.
In addition to introducing the model, we establish non-negativity and investigate the behavior of the compartments S and R, showing monotonicity in the limiting case where γ ≡ 0 is constant. This aligns with biological interpretations.
Equilibrium points are analyzed, and we prove global asymptotic stability for the disease-free equilibrium and local asymptotic stability for the endemic equilibrium, assuming simplifications such as γ(t) ≡ γ (what is reasonable) and β = 1 (subject to further investigation due to the neglect of time-dependent infectivity). We hypothesize that the region with i 0 > 0 is a stable region for (S ∗ , I ∗ , R ∗), even for β < 1. Figures 3 and 4 illustrate Theorems 2 and 3, with different initial conditions and β < 1.
Furthermore, a recent study by 24 provides a global analysis of the equilibria of an SIS model related to the model in equations (2.18)-(2.20) with β = 1, with the main difference being the existence of an entry depending on I in the S compartment. Their geometric approach to global stability problems may extend the results of Theorem 3 globally, although this extension is not straightforward. Also, 24 uses an approximation for the fractional derivative of I that needs to be further explored in the case of the SIR model. Other memory kernels, as those with α ≤ 1 < β , offer opportunities for generalizing the model to various types of diseases.
Finally, we emphasize that the model in (2.40) uses the tempered fractional Riemann-Liouville derivative with an exponential tempering function. This is a notable application of Tempered Fractional Calculus, an extension of classical Fractional Calculus 21. Additionally, Eq. (2.12) relates to an α−exponential function, the usefulness of which has been recently discussed 4.
REFERENCES
- 1 C.N. Angstmann, B.I. Henry & A.V. McGann. A fractional-order infectivity SIR model. Physica A: Statistical Mechanics and its Applications, 452 (2016), 86-93.
- 2 C.N. Angstmann, B.I. Henry & A.V. McGann. A fractional order recovery SIR model from a stochastic process. Bulletin of mathematical biology, 78(3) (2016), 468-499.
- 3 C.N. Angstmann, B.I. Henry & A.V. McGann. A fractional-order infectivity and recovery SIR model. Fractal and Fractional, 1(1) (2017), 11.
- 4 G. Bengochea, M. Ortigueira & L. Verde-Star. The causal α-exponential and the solution of fractional linear time-invariant systems. International Journal of Systems Science, (2024), 1-17.
- 5 R.M. Corless, G.H. Gonnet, D.E. Hare, D.J. Jeffrey & D.E. Knuth. On the Lambert W function. Advances in Computational mathematics, 5(1) (1996), 329-359.
- 6 K. Diethelm. “The analysis of fractional differential equations: An application-oriented exposition using differential operators of Caputo type”. Springer Science & Business Media, Heidelberg (2004).
- 7 C.I. Doering & A.O. Lopes. “Equações diferenciais ordinárias”. Number 517.2 in Coleção Matemática Universitária. IMPA, Rio de Janeiro (2008).
- 8 H.J. Haubold, A.M. Mathai, R.K. Saxena et al. Mittag-Leffler functions and their applications. Journal of applied mathematics, 2011 (2011).
- 9 H.W. Hethcote. The mathematics of infectious diseases. SIAM review, 42(4) (2000), 599-653.
- 10 H.W. Hethcote & D.W. Tudor. Integral equation models for endemic infectious diseases. Journal of mathematical biology, 9(1) (1980), 37-47.
- 11 W.O. Kermack & A.G. McKendrick. Contributions to the mathematical theory of epidemics-I. 1927. Bulletin of mathematical biology , 53(1-2) (1991), 33-55.
- 12 A.A. Kilbas, H.M. Srivastava & J.J. Trujillo. “Theory and applications of fractional differential equations”, volume 204. Elsevier, Amsterdam (2006).
- 13 S.R. Mazorche & N.Z. Monteiro. Modelos epidemiológicos fracionários: o que se perde, o que se ganha, o que se transforma? Proceeding Series of the Brazilian Society of Computational and Applied Mathematics, 8(1) (2021).
- 14 K.S. Miller & B. Ross. “An introduction to the fractional calculus and fractional differential equations”. Wiley, New York (1993).
- 15 N.Z. Monteiro, R.W. dos Santos & S.R. Mazorche. Constructive fractional models through MittagLeffler functions. Computational and Applied Mathematics, 43(4) (2024), 177.
- 16 N.Z. Monteiro & S.R. Mazorche. Analysis and application of a fractional SIR model constructed with Mittag-Leffler distribution. In “Proceedings of the XLII Ibero-Latin-American Congress on Computational Methods in Engineering, CILAMCE, Online” (2021).
- 17 N.Z. Monteiro & S.R. Mazorche. Fractional derivatives applied to epidemiology. Trends in Computational and Applied Mathematics , 22 (2021), 157-177.
- 18 N.Z. Monteiro & S.R. Mazorche. Estudos Numéricos em um Modelo SIR Fracionário. Proceeding Series of the Brazilian Society of Computational and Applied Mathematics, 9(1) (2022).
- 19 N.Z. Monteiro & S.R. Mazorche. Some remarks on an arbitrary-order SIR model constructed with Mittag-Leffler distribution. Matemática Contemporânea, 51 (2022), 25-42.
- 20 K. Oldham & J. Spanier. “The Fractional Calculus theory and applications of differentiation and integration to arbitrary order”. Elsevier, New York (1974).
- 21 M.D. Ortigueira, G. Bengochea & J.T. Machado. Substantial, tempered, and shifted fractional derivatives: Three faces of a tetrahedron. Mathematical Methods in the Applied Sciences, 44(11) (2021), 9191-9209.
- 22 I. Podlubny. “Fractional differential equations: an introduction to fractional derivatives, fractional differential equations, to methods of their solution and some of their applications”. Elsevier, San Diego (1998).
- 23 S.G. Samko, A.A. Kilbas, O.I. Marichev et al. “Fractional integrals and derivatives”, volume 1. Gordon and breach science publishers, Yverdon Yverdon-les-Bains, Switzerland (1993).
- 24 Z. Wu, Y. Cai, Z. Wang & W. Wang. Global stability of a fractional order SIS epidemic model. Journal of Differential Equations, 352 (2023), 221-248.