Fractional dynamics and computational analysis of food chain model with disease in intermediate predator

1.   Introduction

The linkage dynamics between predator and prey explore significant aspects in the field of ecology due to its several applications. Lotka-Volterra model ^[1,2] is the first mathematical framework describing the prey-predator interplay in mathematical ecology. Recently, it has been detected that stage structure models of the population present the interaction dynamics more accurately than other existing models. The recent literature regarding the development of structured models can be seen in ^[3,4,5,6,7]. In the past decade, some models discussed the impact of infectious disease on environmental ecology. These models actually show the spread of disease in populations and their transmission from susceptible to infected species. Kermack and McKendric ^[8] contributed toward the mathematical theory of epidemics through their model. Several researchers investigated mathematical predator-prey models with disease ^[9,10,11]. Freedman and Waltman ^[12] investigated three interacting predator-prey populations. Chattopadhyay et al. ^[13] suggested a predator-prey model. Kar et al. ^[14] modeled a harvested prey-predator system in 2006. In 2007, Dubey ^[15] studied the prey-predator model with a reserved area. In addition, the persistence and extinction of one prey and a two predator system were explored in ^[16].

In the past decade, several mathematical models prepared on the basis of the Beddington-DeAngelis functional response have been derived in ^{[17,18,19,20,21]}. Dubey et al. ^[22] provided the numerical treatment of the fractional food chain problem. In addition, Abdo et al. ^[23] investigated the three-species prey-predator (PP) model pertaining to the Mittag-Leffler kernel. Moreover, Ghanbari et al. ^[24] presented the numerical results of the PP system having a functional response of the Beddington-DeAngelis type. More recently, Liu et al. ^[25] investigated the Turing patterns of the Leslie-Gower Holling type Ⅲ predator-prey model on several different networks with the help of linear stability analysis. Song et al. ^[26] proposed a PP model organized in multiplex networks to investigate the effect of multiplex structure on the diffusion of predator and prey, and furthermore, the influences on the formation of Turing patterns. Alsakaji et al. ^[27] investigated a delay differential model of one-predator two-prey system with Monod-Haldane and Holling type Ⅱ functional responses. Rihan et al. ^[28] studied the dynamics of a fractional-order delayed model of COVID-19 with vaccination efficacy. More recently, Arif et al. ^[29] propounded and discussed a mathematical food chain model (FCM) involving disease in an intermediate predator. This FCM consists of four ordinary differential equations (ODEs) relating four types of species: prey, the intermediate predator (IP), susceptible infected intermediate predator (SIIP), and the apex predator. One of the significant aspects of this model is the inclusion of the Michaelis-Menten type or Holling type Ⅱ functional response ^[30,31] to represent the PP link. A functional response depicts the rate at which the normal predator consumes the prey.

In this work, a fractional order mathematical FCM proposed by Arif et al. ^[29] is investigated and analyzed in view of the Atangana-Baleanu Caputo (ABC) fractional derivative operator (FDO) for the first time. This derivative operator was given by Atangana and Baleanu ^[32] to study the heat transfer model. This derivative operator provides a nonsingular as well as nonlocal kernel carrying the Mittag-Leffler function (MLF) ^[33], which is significantly better than previously established derivative operators. Atangana and Baleanu ^[32] have proposed two versions, the Atangana-Baleanu FDO (ABFDO) in Caputo sense (i.e., ABC derivative), which is a convolution of a local derivative of a given function with the ML function, and the ABFDO in Riemann-Liouville (RL) sense (i.e., Atangana-Baleanu Riemann (ABR) derivative), which is the derivative of a convolution of a given function that is not differentiable with the ML function. The use of a ML kernel in ABFDO is due to its natural appearance in various physical models because the MLF is a joint venture of power-law and exponential-law which induces completely the effect of memory ^[34]. The inspiration behind the selection of ABFDO is the nonlocal characteristic of the kernel which generates the scope of global analysis in those areas where the trends do not follow the power-law. Recently, a number of models were investigated with ABFDO which can be seen in ^{[35,36,37,38,39,40]}. This study also examines the existence and uniqueness related analysis of the solution of the model. In addition, the stability analysis for the food chain problem is also presented utilizing the Ulam-Hyers approach. In the later part of this work, the numerical solution of the model is explored along with simulations.

The rest part of the paper is subdivided as follows: Section 2 provides fundamental definitions and formulae regarding fractional integral and derivative operators. Section 3 gives a basic description of the food chain model. Section 4 discusses the qualitative property of the model. Sections 5 and 6, respectively, present the existence and uniqueness of the obtained solution. In Section 7, Ulam-Hyers stability approach is applied for FCM. In Sections 8 and 9, the numerical solution and simulation are discussed, respectively. Finally, Section 10 concludes the whole work.

2.   Preliminaries

This section presents a quick review of fractional integral and derivative operators.

Definition 2.1. Let $0 \leqslant \theta \leqslant {\theta _0} < \infty$. Then, Banach space $\Omega = E \times E \times E \times E$, where $E = C\left[{0, \, {\theta _0}} \right]$, with the norm

Definition 2.2. ^[41] Let $\Omega$ be a Banach space. The operator $\wp :\Omega \to \Omega$ is Lipschitzian if $\exists$ a constant $m > 0$ such that

where $m$ is the Lipschitz constant for $\wp$. If $m < 1$, $\wp$ is a contraction.

Definition 2.3. ^[32] The Caputo type fractional integral & derivative of a function $M\left(\theta \right)$ of order $\omega$, respectively, are expressed as

where $k \in {Z^ + }$ and $\Gamma$ is a gamma function.

Definition 2.4. ^[32] The AB fractional integral & derivative of $N\left(\theta \right)$ of order $\omega$ are stated as

Here, ${E_\omega }\left(\theta \right)$ signifies the MLF formulated as ^[33]:

and $B\left(\omega \right) = 1 - \omega + \frac{\omega }{{\Gamma \left(\omega \right)}}$ denotes normalized function with $B\left(1 \right) = B\left(0 \right) = 1$.

Definition 2.5. ^[32] The fractional order ABC derivative fulfills the Lipschitz criterion for two functions $M\left(\theta \right)$ and $N\left(\theta \right)$, and the inequality holds as follows:

Proposition 2.6. ^[42] The solution of

is suggested by

3.   Mathematical description of the food chain problem

In this segment, we provide the basic information about the equations and parameters of the FCM.

The mathematical structure of the FCM with three species suggested by Arif et al. ^[29] is represented by means of four nonlinear ODEs in this way:

where $M'\left(\theta \right) > 0$, $N'\left(\theta \right) > 0$, ${N_1}^\prime \left(\theta \right) > 0$, $\Theta '\left(\theta \right) > 0$. Here, $M\left(\theta \right)$, $N\left(\theta \right)$, ${N_1}\left(\theta \right)$ and $\Theta \left(\theta \right)$, respectively, denote functions of time representing population densities of susceptible prey (SP), susceptible IP (SIP), infected IP (IIP), and the top predator (TP), and all parameters are positive constants.

The food chain model given above describes the relation between SP, SIP, IIP, and the TP. ${\alpha _0}\frac{{M\left(\theta \right)\, N\left(\theta \right)}}{{1 + \mu \, M\left(\theta \right)}}$ is the Michaelis-Menten type (or Holling type Ⅱ) functional response and IP becomes infected with relative function $c\frac{{N\left(\theta \right)\, {N_1}\left(\theta \right)}}{{N\left(\theta \right)\, + {N_1}\left(\theta \right)}}$. Parameters are denoted as follows: μ indicates half saturation constant, rate ${\alpha _0}$ is the per capita rate of predation of the IP, rate ${\alpha _1}$ measures the efficiency of biomass conversion from prey to IP, rate ${\alpha _2}$ is the per capita rate of predation of the TP, rate ${\alpha _3}$ measures the efficiency of biomass conversion from IP to TP, rate ${\alpha _4}$ is the per capita rate of predation of the prey, and rate ${\alpha _5}$ measures the efficiency of biomass conversion from infected IP to TP. Furthermore, $r$ is the intrinsic rate of growth of SP. Moreover, $\mathrm{c}$ measures the rate of contact between SIP and IIP while rate $\mathrm{k}$ represents the transformation from IIP to SIP, as this model is known as the SIS model. In this model, ${\mathrm{d}}_{1}$ and ${\mathrm{d}}_{2}$, respectively, stand for natural deaths of intermediate & TPs. Finally, $\mathrm{q}$ denotes harvesting of an IIP.

We present a brief presentation of the model which may indicate the biological relevance of it. Behavior of the entire biological community is assumed to arise from the coupling of the interacting species $M$, $N$, ${N_1}$, and $\Theta$, where the top predator $\Theta$ prey on intermediate predators $N$ and ${N_1}$, and intermediate predators prey on $M$. This is the practical assumption from both mathematical and biological points of view. A specific feature of these food chain systems is that if one species dies out, all the species at higher trophic levels die out as well. It is also assumed that in the absence of the predators the prey population density grows according to a logistic curve with carrying capacity and with an intrinsic growth rate constant $r$$\left({r > 0} \right)$. The consideration of functional response provides motivation to study a food chain model under the framework of nonlinear ODEs.

Now, we replace the classical derivatives in the model (3.1) with ABC fractional derivatives ${}_0^{ABC}D_\theta ^\omega$ of order $\omega$ to capture memory effect in the model in this way:

along with following conditions:

4.   Nonnegative solutions: qualitative properties of the model

In this segment, we discuss qualitative properties of the nonnegative solutions of the fractional FCM (3.2).

For biological reasons, each variable in model (3.2) must be a nonnegative real-valued function. In other words, $\left({M, \, N, \, {N_1}, \, \Theta } \right) \in R_ + ^4$, where $R_ + ^4 = \left\{ {u = \left({{u_1}, \, {u_2}, {u_3}, {u_4}} \right):\, {u_i} \geqslant 0, \, \forall \, i = 1, \, 2, 3, \, 4} \right\}$. Now, we demonstrate that all the solutions of the model (3.2) with (3.3) are absolutely nonnegative.

Lemma 4.1. All solutions of the model (3.2) lie in $R_ + ^4$.

Proof. We define

Then from the FCM (3.2), we attain

${\left. {{}_0^{ABC}D_\theta ^\omega N\left(\theta \right)\, } \right|_{\varphi \left({N)} \right)}} = k\, {N_1} \geqslant 0$, since $0 < k \leqslant 1$ and ${N_1} \geqslant 0$,

Thus, $\left({M(\theta), \, \, N(\theta), \, \, {N_1}\left(\theta \right), \, \Theta \left(\theta \right)} \right) \in R_ + ^4.$ Hence, the lemma is proved.

Theorem 4.2. Consider the subsequent initial value problem (IVP)

where ${}_0^{ABC}D_\theta ^\omega$ signifies the ABC fractional derivative operator and $V\left({\theta, \, U\left(\theta \right)} \right):{\Re ^ + } \times {\Re ^m} \to {\Re ^m}$ denotes a vector field. This system (4.3) definitely has a unique solution on $\left[{0, \, \infty } \right)$ if

(a) $V\left({\theta, \, U\left(\theta \right)} \right)$ and all of its partial derivatives are continuous $\forall$$U \in {\Re ^m}$.

(b) $\left\| {V\left({\theta, \, U\left(\theta \right)} \right)\, } \right\| \leqslant {a_1} + {a_2}\left\| {U\left(\theta \right)\, } \right\|$ for each $U \in {R^m}$ for ${a_1}, \, {a_2} > 0$.

Now, it is easy to establish that the above-mentioned criteria are fulfilled by the set of equations of the model (3.2) with (3.3), and, thus, it confirms the existence of unique nonnegative solutions for the model (3.2) with (3.3).

5.   Analysis of existence of the solution of the FCM

Here, we investigate the existence of a solution of arbitrary order FCM with disease in the intermediate predator. Now exerting AB integral operator of fractional order in the system of Eq (3.2) in the following way:

Making use of the definition of AB fractional integral operator in the system of Eq (5.1), we acquire

For simplified presentation of the above system of Eq (5.2), we express

Theorem 5.1. The kernels ${\rho _1}$, ${\rho _2}$, ${\rho _3}$, and ${\rho _4}$ fulfill the Lipschitz criterion, and contraction of the conditions $0 \leqslant {\lambda _1} < 1$, $0 \leqslant {\lambda _2} < 1$, $0 \leqslant {\lambda _3} < 1$, and $0 \leqslant {\lambda _4} < 1$ are satisfied.

Proof. First, we initiate with the kernel ${\rho _1}$. Suppose $M\left(\theta \right)$ and ${M^ * }\left(\theta \right)$ are two functions for the kernel ${\rho _1}$ fulfilling the conditions $\left\| M \right\| \leqslant {l_1}$ and $\left\| {{M^ * }} \right\| \leqslant l_1^ *$. Similarly, $N\left(\theta \right)$ and ${N^ * }\left(\theta \right)$ are assumed to be the functions for the kernel ${\rho _2}$ satisfying the criteria $\left\| N \right\| \leqslant {l_2}$ and $\left\| {{N^ * }} \right\| \leqslant l_2^ *$, ${N_1}\left(\theta \right)$ and ${N_1}^ * \left(\theta \right)$ are the two functions for the kernel ${\rho _3}$ satisfying the conditions $\left\| {{N_1}} \right\| \leqslant {l_3}$ and $\left\| {{N_1}^ * } \right\| \leqslant l_3^ *$, $\Theta \left(\theta \right)$ and ${\Theta ^ * }\left(\theta \right)$ are the two functions for the kernel ${\rho _4}$ satisfying the conditions $\left\| \Theta \right\| \leqslant {l_4}$ and $\left\| {{\Theta ^ * }} \right\| \leqslant l_4^ *$.

Let

Here, $M\left(\theta \right)$ signifies the bounded function, then

Thus, the kernel ${\rho _1}$ satisfies the Lipschitz criterion. In addition, if ${\lambda _1} \in \left[{0, \, 1} \right)$, then it will also be a contraction for ${\rho _1}$. Similarly,

Let

Here, $N\left(t \right)$ signifies the bounded function, then

Therefore, the Lipschitz criterion is fulfilled for the kernel ${\rho _2}$. Further, if $0 \leqslant {\lambda _2} < 1$, then it is also a contraction for ${\rho _2}$.

Similarly, the Lipschitz criterion is satisfied for the kernels ${\rho _3}$ and ${\rho _4}$ as follows:

where

Moreover, if $0 \leqslant {\lambda _3} < 1$ and $0 \leqslant {\lambda _4} < 1$, it is also a contraction for ${\rho _3}$ and ${\rho _4}$.

In view of the system of Eq (5.3), the system of Eq (5.2) takes the following form:

Now, the subsequent recursive formulae are constructed in this way:

with the following conditions:

Next, we consider difference of the successive terms as

It is worth noting that

Now, utilizing the set of Eq (5.17) along with the triangular inequality, we attain

It is already proved that the kernels ${\rho _1}$, ${\rho _2}$, ${\rho _3}$, and ${\rho _4}$ satisfy the Lipschitz condition, so the set of Eq (5.19) reduces to

Consequently, we attain the following result:

Applying the same procedure, we attain other results as follows:

Taking Eqs (5.21) and (5.22) into account, we acquire the existence of the solution of the considered model. This establishes the theorem.

Theorem 5.2. The fractional FCM involving Holling type Ⅱ functional response expressed by the system of Eq (3.2) possesses a solution if $\exists$${\theta _0}$ in this way:

Proof. It is assumed that $M\left(\theta \right)$, $N\left(\theta \right)$, ${N_1}\left(\theta \right)$, and $\Theta \left(\theta \right)$ are functions of bounded nature. Now, using Eqs (5.21) and (5.22) along with the recursive algorithm, we get

Evidently, Eq (5.18) assures about the existence and smoothness of the functions. Thus, the solution of system (3.2) exists as well as is continuous. Furthermore, to examine that the system (5.14) is a solution of the FCM model (3.2), it is presumed that

Thus, we find

Now, we have

Making the process recursively, we achieve

For $\theta = {\theta _0}$, we obtain

Applying the similar methodology, we further get

Now, employing the limit on inequality (5.29) as $n \to \infty$, we find $\left\| {{F_{1n}}\left(\theta \right)\, } \right\| \to 0$. Implementing the same procedure, we have $\left\| {{F_{in}}\left(\theta \right)\, } \right\| \to 0$, $i = 2, \, 3, \, 4$, and this establishes the theorem.

6.   Uniqueness of system of solutions of the model

Here, we show the uniqueness of the solution of the fractional food chain model (3.2). We assume that ${M^ * }\left(\theta \right)$, ${N^ * }\left(\theta \right)$, $N_1^ * \left(\theta \right)$, and ${\Theta ^ * }\left(\theta \right)$ is another set of solutions for the ABC fractional order model (3.2), then

Taking the norm of equations of system (6.1) provides

Now, employing the results presented in (5.6) and (5.9)–(5.11) in the set of inequalities (6.2), we have

After simplification, we achieve

which produces

Theorem 6.1. The fractional food chain model (3.2) will possess the unique solution if

Proof. If the conditions presented in (6.6) hold, then the set of inequalities (6.5) provides

In view of properties of norm, the set of conditions (6.7) implies that

Thus,

Therefore, the food chain model (3.2) has a unique solution. Hence, the theorem proved.

7.   Ulam-Hyers stability analysis for the model

The Ulam-Hyers stability approach ^[43,44] has been used for problems with fractional derivatives ^{[45,46,47,48]}. Now, we apply this approach to discuss stability for model (3.2) by virtue of nonlinear functional analysis.

Definition 7.1. System (3.2) with (3.3) possesses Ulam-Hyers stability if there exist $\sigma = \max \left({{\sigma _1}, \, {\sigma _2}, {\sigma _3}, {\sigma _4}} \right) > 0$ and $\aleph = \max \left({{\aleph _1}, \, {\aleph _2}, {\aleph _3}, {\aleph _4}} \right) > 0,$for each $\tilde M, \, \, \tilde N, \, \, {\tilde N_1}, \, \tilde \Theta \in E \times E \times E \times E,$ having the subsequent inequalities

then there exists $\left({M, \, N, \, {N_1}, \, \tilde \Theta } \right) \in E \times E \times E \times E$ fulfilling system (3.2) with conditions given by

such that

Remark 7.2. Consider the small perturbations ${h_1}, \, {h_2}, \, {h_3}, {h_4}\, \in C\left[{0, \, {\theta _0}} \right]$ that depend only on the solutions such that ${h_1}\left(0 \right) = 0, \, {\text{ }}{h_2}\left(0 \right) = 0, \, {\text{ }}{h_3}\left(0 \right) = 0, {\text{ }}{h_4}\left(0 \right) = 0$, with the following properties:

(1)

$\left| {{h_1}\left(\theta \right)} \right| \leqslant {\aleph _1}, {\text{ }}\, \left| {{h_2}\left(\theta \right)} \right| \leqslant {\aleph _2}, {\text{ }}\, \left| {{h_3}\left(\theta \right)} \right| \leqslant {\aleph _3}, {\text{ }}\, \left| {{h_4}\left(\theta \right)} \right| \leqslant {\aleph _4}$, for $\theta \in \left[{0, \, {\theta _0}} \right]$, ${\aleph _i} > 0$, $i = 1, \, 2, \, 3, \, 4$. (7.4)

(2) Furthermore, one has

Now,

Lemma 7.3. The solution of perturbed problems

fulfills the relations

where ${\tilde M_{{h_1}}}\left(\theta \right)$, ${\tilde N_{{h_2}}}\left(\theta \right)$, ${\tilde N_{1, {h_3}}}\left(\theta \right)$, ${\tilde \Theta _{{h_4}}}\left(\theta \right)$ are solutions of Eqs (7.7)–(7.10), respectively. Here, $\tilde M$, $\tilde N,$ ${\tilde N_1}$, $\tilde \Theta$ satisfy the set of conditions (7.1) and $m = \frac{{\Gamma \left(\omega \right) - \Gamma \left({\omega + 1} \right) + \theta _0^\omega }}{{B\left(\omega \right)\, \Gamma \left(\omega \right)}}$.

Proof. As suggested by Remark 7.2 and Lemma 7.3, the solutions of Eqs (7.7)–(7.10) are, respectively, given by

Also, we find

It follows from the Remark 7.2 that

Similarly,

Similarly,

Hence, the lemma is proved.

Theorem 7.4. Under the assumptions of Theorems 5.1 and 5.2 and conditions (5.6), (5.9)–(5.11), the system (3.2) with (3.3) will possess Ulam-Hyers stability in $\Omega$.

Proof. Let $\tilde M, \, \tilde N, \, {\tilde N_1}, \, \tilde \Theta \in E$ be the solutions of inequalities (7.1) and the functions $M, N, \, {N_1}, \, \Theta \in E$ be unique solutions of Eq (3.2) with the conditions (7.2). Now, due to the set of conditions (7.2), we obtain

That is,

Hence, the set of Eq (7.19) transforms to

Now, making use of inequality (5.6) and Lemma 7.3, we get

which yields that

where

For ${\sigma _1} = \frac{{2m}}{{1 - {\varphi _1}}}$, we obtain

Similarly, by condition (5.9) and Lemma 7.3, we obtain

which implies that

where

For ${\sigma _2} = \frac{{2m}}{{1 - {\varphi _2}}}$, we obtain

Similarly, we conclude that

where

Thus, for some $\aleph, \, \, \sigma > 0$,

Therefore, the Ulam-Hyers stability of the model (3.2) with (3.3) is established.

8.   Numerical solution of the ABC type FCM model

Here, a numerical approach for the solution of the model is discussed. For this, we consider an IVP with the ABC fractional derivative as follows:

Employing the AB fractional integral operator on Eq (8.1), we get

Taking $\theta = {\theta _n} = n\hbar$ in Eq (8.2), we obtain

Now, the linear Lagrange interpolation of $V\left({t, \, U\left(t \right)} \right)$ provides

where ${U_i} = U\left({{\theta _i}} \right)$.

By using Eq (8.4) in Eq (8.3), the estimated solution of Eq (8.1) is obtained as

where

Using the numerical method (8.5)–(8.7), the solution of the model (3.2) is generated recursively in this way:

The solution of the model (3.2) is achieved by means of the above obtained iterative numerical schemes (8.8).

9.   Numerical simulation

In this part, the above obtained numerical iterative schemes represented by Eqs (8.6)–(8.8) are utilized to perform the simulations for the fractional food chain model (3.2). The following values for the parameters of the discussed model are considered for the simulation purpose ^[29]: $r = 0.6$, ${\alpha _0} = 0.44$, $\mu = 0.406$, ${\alpha _4} = 0.479$, ${\alpha _1} = 0.309$, ${\alpha _2} = 0.45$, $c = 0.202$, $k = 0.095$, ${d_1} = 0.126$, ${d_2} = 0.07$, ${\alpha _5} = 0.292$, $q = 0.5$, ${\alpha _3} = 0.35,$and a time step size $\hbar = 1.0 \times {10^{ - 3}}$.

The graphs presented here show the behavior of solutions for different values of some parameters. The numerical solution of the model is plotted through Figures 1–10 for various initial values and fractional order $\omega$. Figure 1 shows the impact of distinct values of fractional order $\omega$on the dynamics of solutions $M\left(\theta \right)$, $N\left(\theta \right)$, ${N_1}\left(\theta \right)$, and $\Theta \left(\theta \right)$. Figure 2 presents the effect of various values of $\mu$ on the nature of solutions $M$, $N$, ${N_1}$, and $\Theta$. Figure 3 describes the behavior of solutions $M$, $N$, ${N_1}$, and $\Theta$ for distinct values of $r$. Figures 4 and 5, respectively, depict the dynamics of solutions $M$, $N$, ${N_1}$ and $\Theta$ for distinct values of ${\alpha _0}$ and ${\alpha _1}$. Figures 6 and 7, respectively, show the impact of distinct values of ${\alpha _2}$ and ${\alpha _3}$ on the dynamics of solutions $M$, $N$, ${N_1}$ and $\Theta$. Figures 8 and 9 elucidate the nature of solutions $M$, $N$, ${N_1}$, and $\Theta$ for various values of ${d_1}$ and ${d_2}$. Finally, Figure 10 explains the impact of various values of $k$ on the dynamics of solutions $M$, $N$, ${N_1}$, and $\Theta$. It is easy to examine that the achieved numerical results are compatible with the conceptual conjectures in the foregoing sections.

10.   Concluding remarks and observations

This paper examines the fractional order FCM having Holling type Ⅱ functional response. The ABC fractional derivative operator providing Mittag-Leffler kernel is used for this purpose. The best feature of this derivative operator is that it provides scope for inclusion of memory-related properties very efficiently along with utilization of the information. This FDO possesses all the features of the Caputo-Fabrizio derivative with extra features of nonsingular and nonlocal character of a kernel. The main reason behind the use of FDO is the better description of memory characteristics which can be possible through the inclusion of nonlocal features. The use of an ML kernel in ABFDO is due to its natural appearance in various physical and biological models because the ML function is a joint venture of power-law and exponential-law which induces the memory effect completely. The inspiration behind the selection of the AB derivative is the nonlocal characteristic of the kernel which generates the scope of global analysis in those areas where the trends do not follow the power-law. Fractional order modeling also enhances the accuracy of the analysis and provides an extended degree of freedom in the model. However, the most important attribute of the fractional order modeling is to provide an excellent tool for the description of memory and hereditary characteristics which are generally ignored by systems of integer-order derivatives. In addition, fractional order derivatives are also useful to model the behavior of real systems and in the investigation of dynamical systems. The numerical approach is implemented to solve the model to get an approximate solution. The existence and uniqueness related analysis for the model are also presented. In addition, the model is also discussed regarding the Ulam-Hyers stability approach. Moreover, the graphical presentations for numerical solutions of the model confirm the authenticity of the numerical scheme utilized in this paper. It is clear that the achieved numerical results for the model correspond well with the conceptual findings. As a future research scope of the work, the presented numerical approach and stability analysis can also be applied to other physical and biological models.

Author contributions

Jagdev Singh and Devendra Kumar: Conceptualization; Behzad Ghanbari and Kottakkaran Sooppy Nisar: Methodology; Jagdev Singh, Behzad Ghanbari, Ved Prakash Dubey and Kottakkaran Sooppy Nisar: Software; Jagdev Singh, Ved Prakash Dubey and Devendra Kumar: Investigation; Jagdev Singh, Behzad Ghanbari and Kottakkaran Sooppy Nisar: Validation; Jagdev Singh, Behzad Ghanbari, Ved Prakash Dubey, Devendra Kumar and Kottakkaran Sooppy Nisar: Writing-original draft. All authors have read and approved the final version of the manuscript for publication.

Use of AI tools declaration

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

Acknowledgments

The authors extend their appreciation to Prince Sattam bin Abdulaziz University for funding this research work through the project number (PSAU/2024/01/78916).

Conflict of interest

The authors declare no conflicts of interest.