key: cord-0148997-22yw926h
authors: Balzotti, Caterina; D'Ovidio, Mirko; Lai, Anna Chiara; Loreti, Paola
title: Effects of fractional derivatives in epidemic models
date: 2021-07-04
journal: nan
DOI: nan
sha: a0bf27e7af3d12122f4c55c06dd090bf5cbce7eb
doc_id: 148997
cord_uid: 22yw926h

We study epidemic Susceptible-Infected-Susceptible models in the fractional setting. The novelty is to consider models in which the susceptible and infected populations evolve according to different fractional orders. We study a model based on Caputo derivative, for which we establish existence results of the solutions. Also, we investigate a model based on Caputo-Fabrizio operator, for which we provide existence of solutions and a study of the equilibria. Numerical simulations for both models and a direct numerical comparison are also provided.

The interest of the scientific community in mathematical modeling for epidemiology has grown considerably in recent years. The study of epidemic models began in the early 1900s with the pioneering work of Kermack and McKendrick [KM27] . Their idea was to divide the population into groups which distinguish the individuals based on their status with respect to the infection, giving rise to the compartmental modeling for epidemics. The evolution in time of the disease is then described by a system of ordinary differential equations for each considered class.

In this work we focus on the SIS (Susceptible-Infected-Susceptible) model [Het89] , describing infections which do not confer immunity to recovery from illness, such as influenza and common cold. Such theory describes compartmental models, where the population is divided into groups depending on the state of individuals, that is with respect to disease, distinguishing infected, susceptible. The use of mathematical models for epidemiology is useful to predict the behavior of an infection and take strategic decisions in emergency situations to limit the spread of the disease which is microscopically modeled by the fractional order of the derivative.

In recent years, the use of fractional derivatives for epidemic models has grown widely. The main advantage of fractional calculus is that it can incorporate memory effects into the model. Moreover, fractional models have an extra degree of freedom compared to classical models, which is particularly useful for fitting real data when available. We refer to [CLYL21] for a recent review of fractional epidemic models.

In this paper we consider two fractional SIS models. One model is based on Caputo derivative, for which we establish existence results of the solutions and provide numerical simulations. The novelty is to allow the susceptible and infected population evolve according to different fractional orders. The other model is based on Caputo-Fabrizio operator. Here we let the susceptible population evolve according to the Caputo-Fabrizio fractional operator, whereas the infected population dynamics is based on ordinary differential equation. In this case, we rewrite the system as a system of ordinary differential equations, we study the equilibria and present some simulations.

More precisely, let α, α 1 , α 2 ∈ (0, 1). We consider the initial-value problem for Caputo derivative with different orders (1) where β, γ > 0, S 0 , I 0 ≥ 0.

Here, in the formula (1) above we denote the Caputo derivative [Cap08] by The Caputo derivative is well-defined for a function u ∈ C(0, T ), u (s)(t − s) −αi ∈ L 1 ((0, t)), ∀t ∈ (0, T ), i = 1, 2.

The Caputo-Fabrizio operator in (2) is defined by

where M (α) is a non-negative scaling factor satisfying M (0) = M (1) = 1 The case α 1 = α 2 = 1 in (1) has been extensively studied in the literature. Beginning with the models introduced in [Ver38, Ver45, Ver47] in which the logistic equation is used to model population dynamics, a numerous researchers work in the field.

The problem to find a solution to the fractional logistic equation attracted many authors and many works on this topic have been written only providing approximations for that solution. A contribution to this discussion is given in [DL18] based on a series representation of the solution which involves Euler's numbers. However, this approach is not applicable in this new context since it is based on the following property of the Caputo derivative C D α t u(t) = 0 if and only if u(t) = const ∀t ≥ 0, α ∈ (0, 1). Indeed, our effort is to study the fractional indices of derivation in which α 1 may differ from α 2 . We mainly study existence of solution and we propose numerical test based on the method presented in [GLM20] that we illustrate by several pictures. The numerical tests are in agreement with the theory developed in the present paper, and they offer new perspectives (e.g. symmetries emerging in the evolution of the total population N (t) = I(t) + S(t)) for future works. This is not only a challenging problem from the mathematical point of view since it is of interest in many applications in which we may act in different ways in order to slowdown the process.

Recently the Caputo-Fabrizio operator [CF15] has attracted many researchers. The peculiarity of this fractional operator is the presence of non-singular kernel in contrast with the Caputo derivative in which singular kernel appears in the definition. However, the Caputo-Fabrizio operator can be used to model processes with memory.

Our work is to consider (2) in which we are able to reduce the problem to an ordinary differential equation. Here we obtain the solution of the differential equation and their equilibria. The analysis is concluded by numerical simulations. Finally, in order to point out the effects of the Caputo and the Caputo-Fabrizio differential operators, we show a direct comparison between the evolutions of the numerical solutions to the system (1) and the system (2).

We conclude with the plan of the paper. In Section 2 we introduce the SIS model with fractional Caputo derivatives with different orders. In Section 2.1, we investigate the existence of solutions. To this end, we adopt a constructive approach whose ideas are borrowed from the Carathéodory existence theorem for ordinary differential equations, see [CL55, Chapter 2] . In Section 2.2 we complete the study with some numerical simulation. In Section 3 we address a SIS model with fractional Caputo-Fabrizio operator. In Section 3.1 we establish existence of solution by rewriting the system as a system of ordinary differential equations, and we characterize the associated equilibria. Section 3.2 is devoted to numerical simulations, also including numerical tests directly comparing the proposed Caputo and Caputo-Fabrizio models. Finally in Section 4 we draw our conclusions.

In this section we investigate the properties of the Caputo fractional SIS model with different fractional orders (1), which we rewrite here for reader's convenience.

(3)

with α 1 , α 2 ∈ [0, 1]. The particular case of fractional SIS models of Caputo type with α 1 = α 2 has been studied in several works, see e.g. [BDL20, ES13, HOEK18] and references therein. The novelty of our work is the use of different fractional indices. This approach has been already proposed in [LWLT19] , where the authors study an inverse problem to calibrate the parameters for the dengue fever. The results obtained fit well the real data, suggesting that mixed order fractional epidemic models are needed in applications.

Here we establish existence results for the solutions of (3) and we present some related numerical simulations.

We introduce the notation

2.1. Solutions. The existence of a solution to (3) in the time interval [0, T ) (with T ≤ +∞) with (positive) initial data (S 0 , I 0 ) is implied by the existence of a couple of absolutely continuous functions (S(t), I(t)) satisfying the Volterra-type fractional integral equation

f (S(s), I(s))(t − s) α1−1 ds and

for all t ∈ (0, T ). Fix α 1 , α 2 ∈ [0, 1]. In the following we make use of the following function

In the next two results we establish some invariance properties for the field f in the cases γ ≥ β and γ < β, respectively.

Note that below we often make use of the symbols (x 0 , y 0 ) to denote couples of positive real numbers and x(t), y(t) to denote continuous real valued functions. This choice is meant to lighten the notation and we point out that, in the search of solutions of (3), (x 0 , y 0 ) plays the role of initial data (S 0 , I 0 ) whereas x(t) and y(t) represent the (approximated) evolution of the susceptible population S(t) and infected population I(t), respectively. Lemma 1. Let α 1 , α 2 ∈ [0, 1] and let 0 < β ≤ γ and fix ε ∈ (0, 1). Let T = T ε,γ be the positive real number such that Then for any couple of absolutely continuous functions x, y such that (x(t), y(t)) ∈ Q x0,y0 for all t ∈ [0, T ], the associated functions

Proof. Let (x(t), y(t)) ∈ Q x0,y0 for all t ∈ [0, T ], In particular one has x(t) > 0 and y(t)

and this implies

On the other hand,

and this concludes the proof.

The next result deals with the case β > γ and it posits some invariance properties of f , similar to those in Lemma 1, in a time interval [0, T ]. The main difference with Lemma 1 is that in this case the time T depends not only on system parameters (β ≤ γ, α 1 and α 2 ) but also on initial data.

Then for any couple of continuous functions

Proof. Fix x 0 , y 0 > 0. Let x(t), y(t) be two continuous functions satisfying (x(t), y(t)) ∈ Q x0,y0 for all t ∈ [0, T ] with T = T β,γ,x0,y0 satisfying (5). Since

then one has

Arguing as above, one also deduces

We are finally in position to state the main result of the section, namely a global existence for the solutions of (3) in the case γ ≥ β and a local existence result in the case β > γ.

Theorem 1. Let α 1 , α 2 ∈ [0, 1] and let 0 < β ≤ γ. Then for every initial data S 0 , I 0 > 0 the system (3) admists a positive solution in (0, +∞).

If otherwise 0 < γ < β then for every initial data S 0 , I 0 > 0 the system (3) admists a positive solution in (0, T ] where T = T x0,y0,β,γ > 0 satisfies

Proof. In order to have a more light notation, let x 0 := S 0 > 0 and y 0 := I 0 > 0. Also fix ε := 1/2. Define

Note that, if γ < β then T meets the hypothesis of Lemma 2.

Consider Q x0,y0 as in Lemma 1 if γ ≥ β and Q x0,y0 as in Lemma 2 if γ < β. We first prove that (3) admits a positive solution in [0, T ] -note that T is independent from x 0 , y 0 when γ ≥ β. To this end, consider the sequence of functions for n ≥ 1

Set for brevityz n (t) := (x n (t),ỹ n (t)). Claim 1. The sequence of functions {z n (t)} n≥1 is well defined, equicontinuous and equibounded in [0, T ]. In particularz n (t) ∈ Q x0,y0 for all t ∈ [0, T ] and for all n ≥ 1. We prove the claim by showing by induction that for all n and for all k = 1, . . . , n,z n (t) is well defined

where ω is the modulus of continuity of the function G(·; α 1 , α 2 ). If k = 1, then by the definitions in (7) and (8) it follows that z n (t) = (x 0 , y 0 ) for all t ∈ [0, T /n] and the base of the induction readily follows. Fix k such that 1 < k < n and assume thatz n (t) is well defined on [0, kT /n],z n (t) ∈ Q x0,y0 for every t ∈ [0, kT /n], and that (9) holds. Then using the second lines of (7) and (8) we have that the definition ofz n (t) continuously extends to [0, (k + 1)T /n] by setting for t ∈ (kT /n, (k + 1)T /n]

By applying Lemma 1 or Lemma 2 (according to the cases γ ≥ β and γ < β, respectively) to (x(s), y(s)) = z(s) :=z n (s + (k − 1)T /n) (which satisfies (x(s), y(s)) ∈ Q x0,y0 for all s ∈ [0, T /n] by inductive hypothesis) we have thatz n (t) ∈ Q x0,y0 for all t ∈ [kT, (k + 1)T /n]. It is left to show that

The case in which t 1 , t 2 ∈ [0, T /n] is trivial, becausez n is constant in that interval. It is left to discuss the case in which t 1 ∈ [0, T /n] and t 2 ∈ [T /n, (k + 1)T /n]. In this case |t 2 − t 1 | = t 2 − t 1 ≥ t 2 − T /n = |t 2 − T /n|. Then arguing as above one gets

This concludes the proof of the inductive step and, consequently, of the Claim 1.

Now, by Claim 1 and by Ascoli-Arzela's theorem, there exists a subsequence {z n k (t)} converging uniformly in [0, T ] to a continuous limit function z(t) = (S(t), I(t)) satisfying z(t) ∈ Q x0,y0 -recall indeed that Q x0,y0 ⊂ (0, +∞) × (0, +∞) is a compact set. This implies in particular S(t) > 0 and

Then, by Lebesgue's dominated convergence theorem, for every fixed t ∈ (0, T ] and for i = 1, 2

Therefore, for all t ∈ (0, T ], choosing k sufficiently large to have t > T /n k , one hasz n k (t) = (x n k (t),ỹ n k (t)) wherex n k andỹ n k satisfy

On the other hand (x n k (t),ỹ n k (t)) → (S(t), I(t)) as k → ∞ for all t ∈ (0, T ], and recalling x 0 = S 0 and y 0 = I 0 one deduces that

Then (S(t), I(t)) is the required solution of (3) in (0, T ]. If β > γ then we are done, because we proved the local existence of a solution of (3). If β ≤ γ, then we can iteratively extend (S(t), I(t)) to (0, +∞). For instance by applying the result to the initial datum (x 0 , y 0 ) := (S(T /2), I(T /2)) we then obtain a solution defined in [0, T /2 + T ] and so on.

2.2. Numerical simulations. The numerical discretization of system (3) follows [GLM20]. Let us consider a general equation f (u(t) ). and consider a numerical grid which uniformly divides the time interval [0, T ] into N t steps of length ∆t. We denote by u n = u(t n ), with t n = n∆t. Let α ∈ (0, 1), then the Caputo derivative can be approximated as

with C n,0 = g(n), C n,j = g(n − j) − g(n − (j − 1)) for j = 1, . . . , n − 1 and g(r) = r 1−α − (r − 1) 1−α for r ≥ 1. The numerical scheme to solve (11) is then given by

Let us denote by S n = S(t n ) and I n = I(t n ). By applying this discretization to system (3) we obtain We show the evolution in time of S(t), I(t) and their sum S(t) + I(t) as α 1 , α 2 changes: in Figures 1 we considered a case in which β > γ and in 2 γ < β. Note in plots (a) and (b) that the steepness of the solutions S(t) and I(t) in the long run are related to the size of α 1 and α 2 , respectively: this can be interpreted as a time delay effect of the Caputo fractional operator. An interesting phenomenon is also the lack of monotonicity of the solutions, in contrast with ordinary SIS models.

The sum of the two classes, see plots (c), shows that N (t) = S(t) + I(t) is in general not monotone. For instance, in Figure 1 , N (t) first decreases and then increases when α 1 > α 2 , while it first increases and then decreases when α 1 < α 2 (a symmetrical behavior emerges in Figure 2 ). In case α 1 = α 2 then the sum is constant and we recover the theory developed in [BDL20, HOEK18] . A rigorous study of the symmetries emerging in the simulations and, in particular the intersections of all the functions N (t) at the same time, is still under investigation.

In this section we are concerned with the SIS model using Caputo-Fabrizio type fractional derivatives. We refer to [HA20, KSK + 21, MSK19, UKF + 20] for some examples of epidemic models based on the Caputo-Fabrizio fractional operator.

More precisely, here we study the the fractionary SIS system (2), which we recall to be (12) where β, γ > 0, S 0 , I 0 ≥ 0, I 0 + S 0 > 0. Also recall that the Caputo-Fabrizio operator of order α ∈ [0, 1] for a function u ∈ H 1 ((a, b)), a < b is

where M (α) is a non-negative scaling factor satisfying M (0) = M (1) = 1. The indentity

is crucial to reconduct (12) to a system of ordinary differential equations, which is the first step of our investigation. We make the following key assumption, relating the order of derivation with the system parameters:

Note that, for any γ > 0 the condition (14) is satisfied by choosing α sufficiently close to 1 whereas (15) trivially holds choosing M (α) ≡ 1, which is a setting earlier explored and motivated in [LN15] .

3.1. Solutions and equilibria. As anticipated above, we begin by rewriting (12) as a system of ordinary equations for α ∈ [0, 1).

Theorem 2. Let β, γ > 0, S 0 , I 0 ≥ 0, I 0 + S 0 > 0 and assume the conditions (14) and (15). Let

Also define the function g α :

Then any couple non-negative absolutely functions continuous (S(t), I(t)) is a solution of fractional system (12) if and only if

∀t ≥ 0 (16) and I(t) solves

Proof. Preliminarly remark that, by (13) and by the first equation of (12)

By differentiating both sides of the first equation of (12), we deduce that (12) is equivalent to the ordinary system

and remark that (18) implies P (t) = 0 and consequently P (t) = P (0) = P α (S 0 , I 0 ) for all t ≥ 0. Then, if I(t) and S(t) are solutions of (18) (or equivalently, of (12)) then

I(t) I(t) = P α (S 0 , I 0 ).

Solving above equation with respect to S(t) and selecting the only possibly non-negative solution, we obtain for all t ≥ 0

Incidentally notice that if α = 1 then P 1 (I 0 ) = I 0 + S 0 =: N , B 1 = C 1 = 1 we recover from above relation the classical identity S(t) = N − I(t). We check that S(t), as a function g α of I(t), is well defined. To this end set note that g α (x) is defined in R, because

More precisely, above inequality holds because the discriminant of above polynomial reads

Plugging the identity S(t) = g α (I(t)) in the second equation of (12) we obtain (18). Finally we check that S(t) and I(t) are non-negative. By (16), S(t) ≥ 0 if and only if I(t) ∈ [0, P α (S 0 , I 0 )/(α − (1 − α)γ)] for all t ≥ 0. Remark that this condition is satisfied for t = 0, indeed we have

Moreover if I(t) = P α (S 0 , I 0 )/(α − (1 − α)γ) then I (t) = −γI(t) < 0 and, consequently I(t) ≤ P α (S 0 , I 0 )/(α − (1 − α)γ). Finally, remark that 0 is an equilibrium for (18) and that the velocity field f α (I) = β − γ − βI g α (I) + I I is locally Lipschitz continuous. Therefore, by the local uniqueness of the solutions of (18) if I 0 > 0 then I(t) > 0 for all t.

3.1.1. Equilibria. We now characterize the equilibria and we study the asymptotic behavior of (12).

Proposition 1. Assume conditions (14) and (15) and let S 0 , I 0 ≥ 0. The equilibria of the system (17) are 0 and

Setting R := β/γ, E α (I 0 ) > 0 if and only if R > 1. In this case

Proof. The first part of the claim follows by a direct computation, using in particular (15) and the fact that β, γ > 0 implies R > 0. Let x R = E α (S 0 , I 0 ) if R > 1 and x R = 0 if R ≤ 0. We proved that in Theorem 2 that I(t) is non-negative, and consequently [0, +∞) is an invariant set for the dynamics (18). Moreover, by a direct computation one can check that the function V (x) := (x − x R ) 2 is a Lyapunov function for (18) in [0, +∞) and, consequently, x R is a globally asymptocally stable equilibrium for (18) and this concludes the proof.

Note that, as in the classical case, the qualitative properties of the system strongly depend on the reproduction number R = β/γ. The value E α (S 0 , I 0 ) can be viewed as a fractional generalization of the endemic equilibrium of classical SIS models, which is a stable equilibrium when R > 1.

The next result deals with the monotonicity and the asymptotic behaviour of the function N (t) := S(t) + I(t), which is not constant unless we are in the ordinary case α = 1. As we show below, the asymptotic behaviour of N (t) can be used for inverse problems, i.e., the goal to reconstruct the fractional order α of (12) from the observed data on the long range.

Proposition 2. Assume S 0 , I 0 > 0 and α ∈ [0, 1). If R > 1 and if I 0 ∈ (0, E(S 0 , I 0 )) then N (t) := I(t) + S(t) is a strictly increasing function converging to

as t → +∞. If otherwise either R > 1 and I 0 > E α (S 0 , I 0 ) or R ≤ 1, then N (t) is strictly decreasing and tends to P α (S 0 , I 0 )/M (α) as t → +∞.

Proof. We preliminary remark that g α (x) is a convex function, indeed

Moreover, by a direct computation,

and this, together with the convexity of g α , implies

Also remark that, using the assumption (15), i.e., M (α) ≥ α, we obtain

We conclude this preliminary study on g α by noticing that

Now, by Theorem 2 (24) N (t) = I (t) + S (t) = I (t)(1 + g α (I(t)).

Assume R > 1 and I 0 ∈ (0, E(S 0 , I 0 )). Then 

In particular I(t) ≥ I 0 for all t ≥ 0 and, in view of (21) and of (24), we deduce N (t) > 0 for all t > 0, hence N (t) is strictly increasing. Assume now that either R ≤ 1 or I 0 ≥ E(S 0 , I 0 ).

First we discuss the case in which R > 1 and I 0 > E(S 0 , I 0 ). Then I (t) < 0 and, since E α (S 0 , I 0 ) is a (globally asymptotically stable) equilibrium by Proposition 1, then I(t) > E α (S 0 , I 0 ) for all t ≥ 0. Then, in view of (22) and of (24), we deduce N (t) < 0 for all t ≥ 0.

Finally, if R ≤ 1 then I (t) < 0. Since I(t) > 0 for all t > 0 then, in view of (23) and of (24), we deduce N (t) < 0 for all t ≥ 0 and this concludes the proof.

From above result, we can estimate the fractional order α from the asymptotic behaviour of the total population N (t). Indeed, if R > 1 and N (t) → N ∞ then the corresponding fractional order α is the solution of the equation

Assuming as in [LN15] M (α) ≡ 1, and setting N 0 := I 0 + S 0 , above equation reduces to the explicit formula

Note in particular that α = 1 if and only if N ∞ = N 0 , confirming the fact that, in the proposed mixed fractional model, N (t) is constant if and only if α = 1.

3.2. Numerical simulations. The numerical discretization of (12) is based on the results of Theorem 2. Let us consider again the numerical grid introduced in Section 2.2. To compute the discrete evolution in time of I and S we first solve system (17) by a proper ODE solver and then we discretize (16). Specifically, we use the MATLAB tool ode23t to compute I n+1 and then, following (16), we obtain S n+1 = g α (I n+1 ). The parameter β, instead, is set to β = 0.7 in Figure 3 , which implies R > 1, while β = 0.1 in Figure  4 , and thus R < 1. The dashed lines represent the equilibria estimated in Proposition 1. Both S and I monotonically converge to their equilibria. Note that the time delay induced by the fractional order emerges in the fact that the smallest is α, the slowest is the convergence of the related solutions to equilibria. Finally remark that for α = 1 the sum N (t) = S(t) + I(t) is constant, since (12) coincides with the classical SIS model, where if α < 1 then the monotonicity of N (t) is in agreement with Theorem 2. In Figures 5 and 6 we fix α = 0.5 and we set S 0 = 10 − I 0 and let vary the initial data I 0 in {2, 4, 6, 8}. We note that the qualitative behavior of the solutions is in agreement with the theory developed here and it is not much affected by initial data.

3.2.1. Comparison between Caputo SIS model and Caputo-Fabrizio SIS model. We conclude this section with some tests pointing out the effect of the choice of a particular fractional operator on SIS models. We directly compare the numerical solution to the system with Caputo fractional SIS model (3) and Caputo-Fabrizio fractional SIS model (12). To this end we fix the intial data S 0 = 6 and I 0 = 4. We consider the set the fractional orders α 1 and α 2 in (3), the fractional order α in (12) and we set α 1 = α, α 2 = 1, letting α vary in the set {0.2, 0.5, 0.8}. Figure 7 depicts the results in a case in which γ < β, in particular γ = 0.2 and β = 0.7, whereas in Figure 8 we test a case in which β < γ, in particular γ = 0.2 and β = 0.1. Note that the evolution of susceptible population S(t) is mostly affected by the change of differential operator and this effect is augmented by choosing small fractional orders. Finally, note that the Caputo-Fabrizio operator preserves the monotonicity properties of the ordinary SIS case, while the Caputo operator used in (3) yields a more complex structure. 

We explored the effects of fractional differential operators on SIS epidemic models. The main novelty of the present paper consists in letting the susceptible and the infected population evolve with different orders of fractional differential operators. We presented two fractional SIS models with mixed fractional orders. One of them involves the Caputo fractional derivative, characterized by a singular, power law kernel. The other proposed model relies on the recently introduced Caputo-Fabrizio differential operator, which has a non-singular, exponential kernel. For both models we established existence results for the solutions and we conducted a qualitative analysis by means of numerical simulations. In the case of the Caputo-Fabrizio operator, under the assumption that the fractional behavior is restricted to the susceptible population whereas the infected population evolves according to an ordinary differential equation, we were able to move some further step forward in the analysis of the system. Indeed we characterized the equilibria, noticed their strong dependence on the reproduction number according to classical theory, and proposed a method for inverse problems. Finally, we numerically, directly compared the proposed Caputo and Caputo-Fabrizio SIS models, in order to let emerge the effects of each particular differential operator on the system. We believe that the possibility of tuning the memory effects in a single compartment of the population (susceptible/infected) by means of ad hoc fractional orders, may provide finer and effective tools in data fitting and mathematical modeling of epidemic dynamics.

Further possible extensions of the present work include the extension of the proposed methods to other epidemic models, for instance the SIR model, and to controlled, mixed fractional epidemic dynamics. 

Fractional sis epidemic models

Linear models of dissipation whose Q is almost frequency independent

A new definition of fractional derivative without singular kernel

Theory of ordinary differential equations. Tata McGraw-Hill Education

Review of fractional epidemic models

Solutions of fractional logistic equations by Euler's numbers. Phys. A

On the solution of fractional order SIS epidemic model

A contribution to the mathematical theory of epidemics

Study of mathematical model of Hepatitis B under Caputo-Fabrizo derivative

Properties of a new fractional derivative without singular kernel

Novel parameter estimation techniques for a multi-term fractional dynamical epidemic model of dengue fever

A Caputo-Fabrizio fractional differential equation model for HIV/AIDS with treatment compartment

A fractional model for the dynamics of tuberculosis infection using Caputo-Fabrizio derivative

Notice sur la loi que la population suit dans son accroissement

Recherches mathématiques sur la loi d'accroissement de la population

Deuxième mémoire sur la loi d'accroissement de la population