Abstract

In this paper, a fractional model in the Caputo sense is used to characterize the dynamics of HPV with cervical cancer. Generalized mean value theorem has been used to examine whether the infection model has a unique positive solution. The model has two equilibrium points: the disease-free point and the endemic point. The examination of the system’s local and global stability is provided in terms of the basic reproductive number . The global stability analysis has been carried out using an appropriate Lyapunov function and the LaSalle invariant principle. The results demonstrate that in the infection model, if , then the solution converges to the disease-free equilibrium, which is both locally and globally asymptotically stable. Whilst , the endemic equilibrium is considered to exist. Simulations are implemented via a finite difference method with Grünwald-Letnikov discretization approach for Caputo derivative operator to define how changes in parameters impact the dynamic behavior of the system using Matlab.

1. Introduction

The most prevalent sexually transmitted infectious agent in both sexes worldwide is the human papillomavirus (HPV), which gets its name from warts (papillomas) [1]. HPVs are small DNA viruses that attack the cutaneous or mucosal epithelium [2]. There are more than 170 different varieties of HPV that have been found and categorized, and more than 40 of these viruses are the most widespread sexually transmitted ailment in the world [3]. One of the main causes of anal cancer, cervical cancer (CC), and other cancers is HPV. Each year, more than 400000 new cases are reported worldwide, and more than 300 women’s lives were lost to HPV-related causes in just 2018 [4]. Based on the severity of clinical manifestations, genital HPV types are divided into high-risk as subtypes (16, 18, 31, 33, 35, 39, 45, 51, 52, 56, 58, 59, 68, 73, 82) associated with premalignant and malignant cervical, penile, vulvar, vaginal, anal, head, and neck cancers; and low-risk such as subtypes (6, 11, 42, 44, 51, 53, 83) that cause warts or benign, highly proliferative lesions on the genitals [5].

It has been determined that a high-risk HPV DNA sequence, notably HPV 16 and 18, which are present in around 70% of invasive cancers, is present in approximately 99.8% of CC [6]. The World Health Organization (WHO) developed a global strategy to hasten the elimination of CC by 2020 and urged for its eradication as a public health problem in 2018 [7]. It classes as the fourth most frequent malignancy in women, with an anticipated 604,000 new cases and 342,000 fatalities worldwide in 2020 [8].

Infectious disease modeling has drawn a lot of interest recently in an effort to find cures for the diseases and calamities that have beset humanity [9]. A system of equations with state variables and parameters is used in the mathematical model to show the interconnectedness of theories and observations and to investigate approximations and the effects of the parameters as well as anticipate the behavior of the problem over a particular period of time [10] and provide policymakers in public health with information on how to carry out efficient infection control interventions [11]. The fundamental quantities like as birth rates, transmission rates, recovery rates, and mortality rates are expressed by parameters, which are constants incorporated into the equations [12]. Fractional-order differential equation models are utilized as an alternative approach since integer-order differential equations can’t adequately describe experimental and field measurement data. More degrees of freedom and the inclusion of memory effects are advantages of fractional-order differential equation systems over conventional differential equation systems. To put it another way, they offer a useful tool for expressing memory and hereditary aspects that were not inclosed in the classical integer-order system [13]. A variety of fractional operator types are introduced to obtain a deeper understanding of the behavior of the models. These operators have various advantages and disadvantages over each other, such as Riemann-Liouville, Caputo, Caputo-Fabrizio, Hadamard, Katugampola, Atangana-Baleanu, and many more [14]. Furthermore, unlike conventional integer-order derivatives, fractional-order derivatives such as the Caputo–Fabrizio derivative; have a non-singular kernel property. This property makes it very evident that when modeling realistically, the model’s future state depends on both its current and past states [15]. Moreover, the fractional-order explains some characteristics of the dynamical system for the entire time and provides a complete description of the system that covers the entire process space, whereas the classical integer-order derivative addresses certain dynamic characteristics at a specific time [16]. The realization that the majority of challenging dynamical systems are found to be non-singular is recently credited with the enormous increase in the treatment of dynamical systems with fractional-order derivatives. Additionally, they have a lengthy memory, which allows them to provide overall systems effectively [17]. The Mittag-Leffler function and the exponential decay functions, respectively, are used as the kernels of non-singular derivatives like the Atangana–Beleanu and Caputo–Fabrizio fractional derivatives, whilst the power function is used as the kernel of singular derivatives like Caputo [18]. It’s hardly surprising that a lot of models have been researched using fractional-order derivatives given their enormous and wealthy properties. As an illustration, Paul et al. [19] created and studied the SEIR model using fractional-order Caputo derivatives, considering two time-lags: the time required to heal sick individuals and the temporary immune period. Jahanshahi et al. [20] investigated a fractional-order SIRD model in Caputo’s sense with time-dependent memory indexes to include the multi-fractional properties of COVID-19. Caputo fractional derivatives were used by Naik et al. [21] to analyze a fractional-order model for the HIV epidemic’s propagation with optimal control. Chu et al. [22] examined the Holling type II form, a vector-host infectious illness compartmental model with nonlinear saturated incidence and therapy functions. Firstly, the model is formulated mathematically as a nonlinear classical integer-order deferential system. To more accurately represent the dynamics of the disease, the subsequently extended the model to the fractional order by utilizing the well-known Caputo–Fabrizio operator with an exponential decay kernel. Atangana-Baleanu Caputo operator has been used in [9] to examine the diabetes mellitus fractional-order model. Qureshi [14] suggested using Caputo fractional-order operator to create a novel epidemiological system for the measles pandemic. In [23], a model of the dynamics of HIV and malaria transmission with optimal control was created using a Caputo fractional derivative. Owolabi and Edson [24] used the Caputo operator and fixed point theory to model and investigate tuberculosis (TB). Utilizing the Caputo derivative, the dynamical analysis of a fractional-order time-delay glucose-insulin model was carried out in [25]. The researcher in [26] studied a variable-order fractional mathematical model driven by Lévy noise describing the model of the Omicron virus using the concept of Caputo derivative. For additional information on several different forms of fractional derivatives, see reference [27] as well as the references cited therein.

Many researchers created mathematical models to represent the ailment’s dynamics, which helped them to propose ailment control strategies and characterize the dynamics of co-infection with other infectious ailments. The dynamics of HPV and CC (HPV-CC) with preventative strategies including screening, vaccination, and re-vaccination were described using an ordinary differential equation model developed in [28]. Chakraborty et al. [29] developed a mathematical model to analyze how vaccination affects the dynamics of HPV infection prevention and looked into the viability of a vaccination strategy, illustrating analytically and numerically how vaccination can ensure a predictable preventive policy against disease transmission among sexually active individuals who are more likely to develop CC from high-risk and low-high-risk papillomavirus strains, whilst the same mathematical model in the fractional-order with Atangana-Baleanu derivatives was studied, and numerical solutions were obtained using the Adams-Bashforth-Moulton method by [30]. The most crucial epidemiological aspects of HPV infection and related cancers were incorporated into a two-sex deterministic mathematical model that was developed by the researchers in [31]. The model included catch-up vaccination for adults and school-based vaccine administration for teenagers to evaluate the population-level effects of HPV immunization programs. The center manifold theorem, normal forms theory, and the next-generation operator were all used to thoroughly study the model’s dynamics. For several conceivable scenarios, they created an optimal control problem to identify the best HPV vaccination deployment method. They proved that there are optimum control issue solutions and used Pontryagin’s Maximum Principle to describe the prerequisites for optimal control solutions. Akimenko and Fajar studied the stability of the age-structured model of CC cells and HPV dynamics [32]. In [33], the author has researched how antiviral treatments affect the spread of cervical cancer. The model includes two pharmacological therapies. The function of the first one is to prevent new infections, while the second’s function is to prevent viral replication. The numerical outcomes have demonstrated the role of the provided medication in regulating HPV infection and cancer cell growth rate. In [34], the mathematical model for the dynamics of HPV transmission in-host in the presence of an immune response represented by Cytotoxic T-lymphocytes cells has been studied; the model presented taking into account the effects of latent HPV infections, and the dynamics of the model were successfully analyzed. A new mathematical model of CC based on the age-structured of cells at the tissue level has been put out by researchers in [35]. Goshu and Alebachew [36] created a mathematical model for the spread of CC transmission disease in the presence of vaccination and therapy. Omame et al. [37] created and presented a coinfection classical integer-order model for syphilis and HPV with cost-effectiveness analysis and optimal control. While the fractional-order in the Caputo–Fabrizio sense of the coinfection model for HPV and Syphilis is investigated using the nonsingular kernel derivative [38]. HPV and Chlamydia trachomatis coinfection mathematical model with cost-effectiveness optimal control analysis has been formulated and examined. It has been investigated how HPV screening, Chlamydia trachomatis therapy, and both diseases’ preventive measures affect the management of their coinfections, the consequent prevention of malignancies, and pelvic inflammatory disease [39], and Nwajeri et al. [40] researched the fractional-order with Caputo derivatives of the codynamics model for HPV and Chlamydia trachomatis, and created the numerical simulations utilizing the fractional predictor-corrector method. In [41] they developed an optimal control strategy for the HPV-Herpes Simplex Virus type 2 codynamics model that minimizes the cost of implementing controls while also minimizing the number of infectious individuals over the intervention interval. The Runge-Kutta forward-backward sweep numerical approximation method was used to implement the optimal control systems.

In light of these accomplishments, we are inspired to investigate the HPV-CC model in this study with the Caputo fractional-order operator, which is best appropriate for modeling biological and physical facts. The motivation to continue our investigation with the Caputo derivative which is a modification of the Riemann–Liouville definition is that this type of fractional operator has advantages concerning other present derivatives; the Caputo derivative of a constant function yields zero, which is usual in mathematics. The Caputo operator first solves an ordinary differential equation, then it takes a fractional integral to get the desired fractional derivative order. More crucially, Caputo’s fractional differential equation allows the use of local initial conditions to be included in the derivation of the model, and finally, the inclusion of the memory effect on the Caputo derivative.

The purpose of this study is to investigate the dynamics of the fractional-order of the HPV-CC model, using the Caputo derivative. The organization of the paper is as follows: Some essential definitions and fractional calculus properties are given in Section 2. The fractional HPV-CC model is introduced in Section 3, along with proof of the solutions’ existence and uniqueness. In Section 4, it is shown that the suggested model is locally stable at the disease-free and endemic equilibrium points. The global stability of the suggested model at the disease-free and endemic equilibrium points is examined in Section 5. In Section 6, we construct a finite difference scheme for the suggested model and show that it maintains the boundedness and positivity of the solutions to the investigated model; here, we will discuss the asymptotic stability of this scheme. To demonstrate the applicability and effectiveness of the finite difference method, some numerical simulations for the suggested model are shown in Section 7, and the final Section is devoted to the conclusion.

2. Preliminary Concepts

Some basic results pertaining to fractional calculus are presented in this section.

Definition 1. Laplace transform of the function is defined by the following improper integral,

Definition 2. MittagLeffler function is defined as

Definition 3. Mittag-Leffler function of two parameters (generalization of Mittag-Leffler function) is given by

Lemma 4. For and , we acquire

Definition 5. Riemann-Liouville fractional integral of order of a function is defined by

Definition 6. Caputo fractional derivative of order for a function where and is given by

Definition 7. Caputo fractional derivative on the half axis of order of a function and is defined as followsWhen , (7) takes the the following form

3. HPV-CC Model Formulation

This work examines a modified version of the classical integer-order HPV-CC model that was previously examined in [42]. The model’s classes have been defined as follows: The total individual population who are sexually active at time , denoted by , is divided into five compartments: susceptible individuals , individuals exposed to HPV , infectious individuals showing no symptoms of HPV , infectious individuals showing symptoms of HPV , and individuals infected with CC . Thus,

Susceptible individuals acquire HPV infection with an infection force ofwhere, represents the rate of HPV infection transmission, represents effective contact rate, , and represent the relative infectiousness of and respectively with .

With probability where , individuals in class develop to individuals in class . individuals in class develop to individuals in class with probability . After experiencing an HPV symptom, individuals in class proceed to individuals in class with a rate . Individuals in class compartments may have CC with the rate of development . Based on the formulations and presumptions mentioned above, the HPV-CC model takes the following form (associated biological parameters of the model are shown in Table 1) [42]:

With the following non-negative initial conditions,

The reason for examining the fractional-order case is the noteworthy uniqueness of these fractional-order systems with hereditary qualities and non-local features (memory) that have not been observed with the integer-order differential operators commonly found in biology. Additionally, mistakes resulting from overlooked parameters can be minimized by modeling real-life processes with fractional-order differential equations [21]. As a result, our suggested fractional-order model for HPV-CC using the Caputo derivatives is:

Subject to the initial conditions (11), where , , and . If then Model (12) reduces to Model (10).

Biological parameters in Model (12) have been modified to make sure that the left and right-hand sides of the model have the same dimension as follows (the schematic diagram of Model (13) is shown in Figure 1):here,

3.1. Existence and Uniqueness of Non-Negative Solutions

Consider the following initial value problem of autonomous nonlinear Caputo fractional-order derivative systemcombined with the initial conditionwhere and the function is called vector filed.

Theorem 8 (see [43]). Suppose that attains the following conditions:(i) and are continuous,(ii), ,For almost each and every . Then, the solution of (14) and (15) is not only existent but also unique on .

Lemma 9 (see [44]). (Generalized mean value theorem) Assume that and with are continuous for each and respectively, then, one has

Lemma 10 (see [45]). Assume that and with .

(i)If . then is increasing .(ii)If . then is decreasing .

Theorem 11. For any initial condition fulfilling (15), System (14) has a unique solution on . Additionally, is still in and bounded.

Theorem 12. Model (13) has a unique solution on for positive initial conditions (11) and this solution is still in for every . Furthermore,

Proof. Let us reformulate Model (13) in the form of Caputo fractional derivative system of order , as followswhere , , ,Here , , and .
It is clear that vector function satisfies the first condition of Theorem 11. The second one needs to be proven. From System (17), we get,Hence, the second condition of Theorem 11 is proven. Next, we need to prove that this solution is nonnegative. From (13), we have , for all , for all , for all , for all .
We conclude that the solution of Model (13) will remain in based on Lemmas 9 and 10.
Finally, we demonstrate that the solution is bounded. Adding Model’s (13), results inTaking Laplace transform in (19) into account, we get
Consequently, we have
It follows that as , This completes the proof.

4. Local Stability

In this section, the local stability analysis of equilibrium points will be covered.

Theorem 13. The equilibrium points of Model (13) are locally asymptotically stable if the following Matignon condition is satisfiedwhere stands for the class of all eigenvalues of the Jacobian matrix of Model (13).

4.1. Local Stability Analysis of Disease-Free Equilibrium Point

Model’s (13) equilibrium points are derived by setting Model’s (13) right side to zero. The disease-free equilibrium point is given by .

The threshold value of Model (13) which is known as a basic reproduction number can be computed using the next-generation matrix indicated in Van den Driessche and Watmough [46]. The transmission matrix and transition matrix of Model (13) are obtained respectively as follows

Then, we compute the derivative of matrices and with respect to infected classes at , and we get Jacobian of and , respectively so that

The aforementioned matrices are utilized to determine for Model (13) using the spectral radius . That is, , and is given by

Theorem 14. is locally asymptotically stable if and unstable if .

Proof. The linearization matrix of Model (13) around , is:here, .
Then the characteristic equation of the linearized Model (24) iswhere , , and .
It is clear that are negative. The three remaining eigenvalues are the roots of the characteristic equationThe necessary and sufficient conditions on three-order Routh-Hurwitz determinants such that the zeros of (26) have negative real parts or (20) to be satisfied are and .
It is easy to show thatandAll stability conditions have been met, then, is locally asymptotically stable if and unstable if .

4.2. Local Stability Analysis of Endemic Equilibrium Point

The endemic equilibrium point of Model (13) denoted by , where , , , , .

Theorem 15. is locally asymptotically stable if and unstable if .

Proof. The linearization matrix of Model (13) around , ishere , , , , and .
Matrix (30) eigenvalues yield the following fifth-degree polynomial equationwhere , , , and .
The following eigenvalue is easily obtained . The remaining eigenvalues are the roots of the characteristic equationThe necessary and sufficient conditions on fourth-order Routh-Hurwitz determinants for the zeros of (32) to have negative real parts or (20) to be fulfilled are and .
It is simple to establish that , and
Finally, we show that .where . has just been demonstrated to be locally asymptotically stable if and unstable if .

5. Global Stability

In this section, we examine the global stability for model (13) by building appropriate Lyapunov functions. In order to show global stability, we take the function into consideration.

Lemma 16 (see [47]). Assume be a continuous function. Then, for any time

Theorem 17. is globally asymptotically stable whenever .

Proof. Suppose the following Lyapunov function:Clearly, and .
Calculate along the solution of model (13) asTherefore, guarantees for all that
Moreover, it is easy to confirm that . Hence, is the largest compact invariant subset of . We conclude that is globally asymptotically stable provided as a result of LaSalle’s invariance principle [48].

Theorem 18. is globally asymptotically stable when .

Proof. Let us construct the following Lyapunov function:Obviously, and . Evaluate along the trajectories of model (13), yieldsAt , we possess where Thenwhich according to the arithmetic mean-geometric mean inequality, is less than or equal to zero. Thus, ensures for all thatFurther, it is clear that . Therefore, is the largest compact invariant subset of . According to LaSalle’s invariance principle [48], is globally asymptotically stable so long as .

6. Construction of Finite Difference Scheme

To acquire numerical solutions for HPV-CC fractional-order model (13), we create the finite difference scheme in this section. To estimate the fractional derivative, one can utilize Grünwald-Letnikov approach. The suggested finite difference scheme maintains the positivity and stability of the equilibrium points of the solution corresponding to the continuous epidemiological models of fractional order.

Let , a finite interval. The usual notation and are used, where and denotes the numerical solution of the analytical solution of at the grid points .

Grünwald-Letnikov discretization approach based on Caputo derivative is defined aswhere and

Lemma 19 (see [49]). Suppose that , therefore the coefficients and fulfil for the properties

The finite difference scheme of Model (13) using Grünwald-Letnikov discretization approach (41) is as followshere

Because each of these equations is linear in , and , hence through some calculations the following explicit expressions can be obtained

6.1. Non-Negativity and Boundedness of Finite Difference Scheme

This subsection examines a few characteristics of proposed Scheme (44). Keeping in mind that system (13) have unique non-negative solutions and also all the parameters are positive.

Theorem 20 (Non-negativity). Assume that in (44) , and , then , and is satisfied for all .

Proof. The induction will be used to prove this theorem. For in Scheme (44), we getWe presume that in system (44), , and . So, for

Theorem 21 (boundedness). Assume that in Scheme (44) . Then there exists a constantsuch that

Proof. When each equation in Scheme (44) is multiplied by its denominator, the result ishere and
Applying the induction principle and utilizing Lemma (19), it follows that for For , we haveFor , we havehere, and .
Now, we assume that for , isThus for , we gethere, and Thus,

6.2. Stability of Finite Difference Scheme

This subsection investigates the stability of Scheme (44).

Definition 22. Scheme (44) is said to be asymptotically stable, if there exist constants as , such that Leads for any arbitrary initial values , and .

Theorem 23. Assume that the supposition of Theorems (20) and (21) are provided. Hence, Scheme (44) is asymptotically stable.

Proof. According to Theorem (21), we deduce that Scheme (44) is asymptotically stable.

7. Numerical Simulations

There are no general techniques for solving systems of fractional differential equations analytically, similar to the classical theory of differential equations. The fractional case is much more challenging to treat even approximately [50]. There are various techniques based on a continuous expansion formula for the fractional derivative that, in some circumstances, can be used to approximately solve the original fractional-order model [51].

In this section, we’ll simulate the solution to HPV-CC fractional-order Model (13) using suggested Scheme (44). The presented model’s parameter values are either extrapolated from earlier model studies or based on assumptions. Table 2 lists the values of these parameters for the two scenarios and .

Proposed Scheme (44) is implemented with a time step size , and figures of the numerical solutions are presented using varied initial conditions that satisfy Theorem 22 and various values of and derivative order .

The convergence of all solutions toward the equilibrium point is depicted in Figures 24 across a range of different values for and . This result was attained when . The infected classes converge to zero over time while the susceptible compartment rises initially and then converges to . is thereby shown to be globally asymptotically stable. When , all solutions converge toward the equilibrium point which is globally asymptotically stable. This result is shown in Figures 57 for a variety of and values indicating that upon the introduction of HPV, the susceptible compartment over time gradually reduces and stays constant and can’t be treated. While the exposed compartment , the asymptomatic compartment , the infected compartment , and the compartment affected by cervical cancer gradually increase with time .

Consider HPV-CC fractional-order Model (13) with initial condition and the values of parameters when given in Table 2 with . In this scenario, and hence, Model (13) has one disease-free equilibrium and it is globally asymptotically stable, which means that the disease disappears. This result is illustrated in Figure 3. Similarly to the second scenario when with . We get . Consequently, Model (13) has a unique endemic equilibrium . Additionally, it is globally asymptotically stable. This result is sketched in Figure 6 and illustrates the continued prevalence of the disease among the population.

Figures 4 and 7 provide the numerical results for the distinction between fractional and integer orders. Figures 4 and 7 demonstrate that, in comparison to classical integer-order models, differential equations with fractional-order derivatives have wealthy dynamics and more accurately depict biological systems.

Drawing from the preceding discussion and the numerical results presented in Figures 27, we deduce that may be influenced by an individual’s prior experience with or understanding of the illness. Consequently, the numerical results validate that differential equations with fractional-order derivatives explain biological systems more accurately than classical integer-order models.

Depending on the values of the parameters mentioned in Table 2, Tables 3 and 4 show how has an impact on decreasing and increasing , respectively. The prevalence of HPV-CC is influenced by the level of interpersonal contact within a community .

By analyzing these variations, we may predict more about how fractional orders affect the model’s dynamics and how crucial the parameter is in influencing population dynamics. Moreover, Figures 27 present results that provide a basis for additional investigation and refinement of the model, which could result in the development of more potent disease control and prevention techniques.

8. Conclusion and Future Directions

In this article, we proposed a fractional-order HPV with cervical cancer (HPV-CC) model with Caputo derivative of order . This dynamical model was more compatible for describing biological phenomena with memory than the integer-order model. Some mathematical results related to the model are presented. The local stability of the model for and was given. The model had two equilibrium points: the disease-free point and the endemic point . The examination of the system’s local and global stability was provided in terms of the basic reproductive number . Using LaSalle invariant principle and the appropriate Lyapunov function, the analysis of global stability had been conducted. The results demonstrated that in the infection model, if , then the solution converged to , which was both locally and globally asymptotically stable. Whilst , was considered to exist. To determine how changes in parameters affect the dynamic behavior of the proposed system, simulations were constructed using a finite difference scheme with Grünwald-Letnikov discretization approach for Caputo derivative operator. The outcomes acquired demonstrated that the finite difference method is a precise and effectual technique to obtain the numerical solution of the suggested nonlinear fractional-order HPV-CC model.

In future work, the presented model after being modified, will be integrated with HIV-AIDS; employing nonlocal and nonsingular kernel operators such such Caputo-Fabrizio. The generalized Adams-Bashforth Moulton method will be applied for the numerical simulations.

Data Availability

No data were used to support this study.

Conflicts of Interest

There authors declare that there are no conflicts of interest.