Analysis and analytical simulation for a biophysical fractional diffusive cancer model with virotherapy using the Caputo operator

1.   Introduction

Cancer stands as one of the most lethal diseases worldwide, ranking as the second leading cause of death globally. According to the World Health Organization, approximately 9.6 million lives were claimed by cancer in 2018 alone, accounting for about 1 in 6 deaths. Developing countries bear the highest impact of this burden, contributing nearly 70% of the annual cancer-related fatalities. In response, policymakers have recognized the significance of prioritizing cancer prevention and control to safeguard the nation's productivity [1],[2]. Researchers have also devoted considerable time and effort to discovering effective treatments, enhancing the efficiency of current low-cost therapies and exploring methods to boost patients' immune systems in their battle against cancer [3],[4].

In recent times, significant progress has been made in developing cancer therapies that can effectively target tumors without harming healthy neighboring tissues. Genetic engineering has played a pivotal role in this pursuit, leading to the discovery of a promising cancer treatment involving genetically altered viruses [5]. These engineered viruses, known as oncolytic viruses, have the unique ability to infect cancer cells specifically. They grow within the abnormal tumor cells and destroy them without affecting the surrounding healthy cells or normal tissue. The process involves the oncolytic viruses interacting with the tumor cell, leading to a burst of new oncolytic viruses. The burst size, which measures the number of new viruses resulting from the lysis of an infected tumor cell, serves as an important indicator of the oncolytic virus's replicability [6]. The treatment of cancer varies according to individual situations, and different therapies may be employed alone or in combination, including surgery, radiotherapy, chemotherapy, hormone therapy, immunotherapy and virotherapy. Virotherapy has emerged as a promising approach, wherein a reprogrammed virus, termed an oncolytic virus, is utilized. Oncolytic viruses effectively infect and destroy cancer cells, harnessing the cell's genetic machinery to replicate themselves and spread to neighboring uninfected cells [7]. This innovative treatment holds immense potential in the ongoing fight against cancer.

Biophysics is a branch of physics that is concerned with studying biological phenomena depending on physical theories, laws, methods and physical quantities. The biophysics concept was introduced in 1892 by Karl Pearson and then expanded and developed to solve extremely difficult problems cannot be solved. Biophysics includes different biological and medical branches, which will be listed as the following: Firstly, it is concerned with studying the structure of many biological phenomena, like DNA, and analyzing its data depending on physical laws and theories. Another biophysical field is the construction of mathematical models (for viruses, pathogens, proteins etc.), developing the needed computer codes, and performing the simulation of the studied case to find the required target. Moreover, building models computationally of the human body, like the brain and nervous systems, is an important biophysical field that gives a better understanding of their working and also tells us how molecules, like blood hormones and other liquids and gases, move around the body. Biophysics used physical and geometrical optics to develop medical imaging like PET scans, CT scans and MRIs. It developed medical devices that are used in radiation therapy, cardiology, defibrillators, pacemakers, artificial heart valves and other fields. It's also used in ecosystems that give us detailed information about the human body, like a pregnant woman's follow-up.

Biophysical models that are expressed mathematically have proven valuable in understanding and addressing various biological and medical challenges, encompassing cancer models involving diverse therapies [8],[9], prey-predator models [9], epidemic models [10]–[13], HIV infection models [14]–[16], COVID-19 models [17], food chain models [18], smoking models [19], glycolysis models [20], drinking population models [21] and vector–host disease models [22]. Among these, some researchers have investigated mathematical models concerning cancer treated with virotherapy, investigating the intricate interaction between the virus and the tumor. Wodarz's foundational biophysical model explored tumor growth in the presence of virotherapy treatment [23], and subsequent modifications of the model aimed to explore the correlation between burst size and virus replicability. The findings revealed that larger burst sizes correspond to a decrease in cancer cells [24]. Conversely, other models have focused on understanding the immune system's response to virotherapy. As the immune system perceives viruses as foreign entities and acts to eliminate them, it may manifest a negative response to virotherapy, potentially compromising the effectiveness of the viral treatment [25]. Numerous studies have also explored the dynamics of interactions between uninfected and infected cells, as well as various immune responses [26],[27]. In this context, a review of Wodarz's fundamental model of oncolytic virus replication [28],[29] is utilized. This model investigates tumor growth and the infection term through a system of ordinary differential equations (ODEs). In this work, the branch of biophysics is concerned with developing a mathematical model in medicine, where the analysis and analytical simulation for a diffusive cancer model in the fractional form is considered. The present study aims to conduct a comprehensive analysis and develop analytical simulations for a fractional diffusive cancer model incorporating virotherapy, utilizing the Caputo operator. The objective is to gain a deeper understanding of the dynamics of the cancer-virus interaction and evaluate the potential efficacy of virotherapy as a treatment strategy for cancer.

The authors in [30] have reformulated the model of oncolytic virus replication that studies the tumor growth and the infection term common with virotherapy and immune response and its impact in both uninfected and infected cells by considering five populations, as follows:

subject to the initial population data ϕi(0,x)=ϕi0, i = 1, …, 5. In model (1), the parameters a_i, i = 1, …, 5 denote the diffusion terms. The five-cell population are ϕ₁(t, x) uninfected cancer cells, ϕ₂(t, x) infected cancer cells, ϕ₃(t, x) free virus, ϕ₄(t, x) naive immune cells and ϕ₅(t, x) activated immune cells. The description of parameters that appear in (1) can be arranged as, see [31],[32].

In this work, we introduce a moderation for model (1) by substituting the time derivative with the Caputo fractional derivative. This modification makes the sides of the equations in the model (1) have not the same dimension, which requires using an auxiliary parameter σ, see [33]. Based on this discussion, we present the following fractional diffusive cancer model with virotherapy:

where Dc,tα denotes the Caputo fractional derivative concerning the time t of order α > 0. The urgent need to progress cancer research serves as the inspiration for our paper. Understanding the complicated dynamics of cancer growth and the effectiveness of treatment approaches is crucial since cancer creates significant issues around the world. Our research is driven by the objective to advance cancer research through illuminating cancer modeling and therapeutic strategies. We specifically look at how fractional calculus can be used, especially in the Caputo sense, to improve our understanding of how cancer progresses. In addition, our research explores the potential topic of virotherapy, which uses viruses to target and treat cancer cells. We examine the potential of virotherapy and assess its potency as a cancer treatment strategy within a mathematical modeling framework.

The existence and uniqueness of the governing model (2) are examined in this paper. In order to solve the cancer model, it applies the Laplace residual power series method (LRPSM) [34] as an analytical tool, demonstrating its applicability and usefulness in solving real-world issues, which combine the power series method and the Laplace transform method [35]

The paper is organized as an introduction in the first section. In Section 2, we present the definition of the Caputo derivative and some of its properties. Section 3 is devoted to investigating the existence and uniqueness of the proposed model. The LRPSM is employed to obtain analytical solutions in Section 4. In Section 5, we present some obtained numerical results in tables and figures. Discussion about the obtained results is presented in Section 6. Finally, some conclusions are provided in Section 7.

2.   Preliminaries

This section is devoted to introducing some basic concepts about fractional calculus and the essential properties and results that will be useful in our work.

Definition 1. [36] For an integrable function ϕ, the Caputo derivative of fractional order α ∈ (0,1) is given by:

where t ≥ 0, m = [α] + 1. The Riemann-Liouville fractional integral of order α, Re (α) > 0 is given by:

Lemma 1. [36] For 0 < α < 1 and t ≥ 0, we get the following property for the Caputo derivative:

Definition 2. [36] Let ϕ(t, x) be a piecewise continuous function on [0, ∞) × I and of exponential order η. The Laplace transformation (LT) of ϕ(t, x) is defined as:

while the inverse Laplace transformation of Φ(s, x) is given by:

where ν₀ lies in the right half plane of the absolute convergence of the Laplace integral.

Lemma 2. [37] Let ϕ(t, x) be a piecewise continuous function on [0, ∞) × I and of exponential order η. Then the Laplace transformation of Dtαϕ(t,x),0<α<1 is given by:

3.   Existence and uniqueness analysis

This section is devoted to show that the fractional model (2) has a unique solution using the fixed point theorem. Toachive our goal, the model (2) can be written as:

such that

Apply the Riemann-Liouville fractional integral to both sides of the model (9) with the aid of Lemma 1 and (10) to get:

To show that the kernels Fi(t,ϕi),i=1,…,5 satisfy the Lipschitz condition we introduce the following assumption:

C1: For the continuous functions ϕi and ϕ ¯i,i=1,…,5 belong to L[0,1], there exist constants ωi,ω ¯i,ω*∈𝔑,i=1,…,5, such that the following hold true:

for i = 1, …, 5. Now, we can introduce the following result:

Theorem 1. The kernels Fi(t,ϕi),i=1,…,5 satisfy the Lipschitz conditions and contraction, provided that the assumption C1 and the following inequalities are satisfied:

Proof. Using the definition of F1(t,ϕ1) in (10) with aid of assumption C1, we get:

Now, we show that the kernel F2(t,ϕ2) satisfies the Lipschitz condition.

We are going to prove that F3(t,ϕ3) satisfies the Lipschitz condition.

Regarding the kernel F4(t,ϕ4), we can see that it achieves the Lipschitz condition through the following:

Finally, for the kernel F5(t,ϕ5), we obtain:

From the discussion in (20)–(24) with aid (15)–(19), we obtain that the kernels Fi(t,ϕi),i=1,…,5 satisfy the Lipschitz conditions and contraction.

According to (11), we define the following recursions:

We take the norm for the sides of (25) with the help of Theorem 1 and we obtain:

for i = 1, …, 5. Now, we can present the following result:

Theorem 2. The fractional diffusive cancer model (2) has a solution provided that the assumption C1 and the following inequality holds the true:

Proof. Firstly, we define:

Using the results obtained in (26), we obtain:

for i = 1, …, 5. This means the functions Qi,q(t)→0 as q → ∞ for i = 1, …, 5. This completes the proof.

For the uniqueness analysis for the solution of the model (2), we introduce the following Theorem:

Theorem 3. The fractional diffusive cancer model (2) has a unique solution provided that the assumption C1 and the following inequality hold:

for i = 1, …, 5.

Proof. Assume that the model (2) has another pair of solutions ϕi*(t,x),i=1,…,5. Then, ϕi* satisfies the integral system:

Consequently, with the help of (11) and (31), one has:

for i = 1, …, 5. Thus, we conclude:

Therefore, using the assumption (30) and the obtained result (33), we observe that ∥ϕi(t,x)−ϕi*(t,x)=0∥ for i = 1, …, 5, which implies ϕi(t,x)=ϕi*(t,x).

4.   Laplace residual power series solutions

In this section, we try to utilize the LRPSM to construct an analytical solution for the fractional diffusive cancer model (2). This method was investigated in many works in the literature. As a first step, we take the Laplace transformation with respect to temporal variable t to both sides for the equations in the biophysical model (2), and we get:

Assume that L{ϕi(t,x)}=Φi(s,x) for i = 1, …, 5. Using Lemma 2, then we can rewrite Eqs (34)–(38) as:

The proposed technique assumes that the solution of the system (39)–(43) can be given in the following form:

for i = 1, …, 5. Hence, the nth-truncated solution can be written as:

for i = 1, …, 5. It is important to note here that the discussion of series convergence (44) has been presented in detail in the literature. To determine the initial guess b_i,0(x), we present the following Lemma:

Lemma 3. [37] Let ϕ(t, x) be a piecewise continuous function on [0, ∞) × I and suppose that Φ(s,x)=L{ϕ(t,x)}. Then,

According to this lemma, we obtain the initial guess bi,0(x)=ϕi0, for i = 1, …, 5. Therefore, the nth-truncated solution (45) is:

for i = 1, …, 5. Now, define the Laplace residual functions, LRes_i(s, x), related to (39), (43) as:

Consequently, the nth-truncated Laplace residual functions, LResin(s,x) gives:

Substitute the nth-truncated solution (47) into (53)–(57) and multiply both sides by snα+1,n=1,2,…, then solve the obtained algebraic system:

recursively, we can construct the desired coefficients bi,r(x),r=1,2,…,n,i=1,…,5. Substitute the obtained coefficients bi,r(x),r=1,2,…,n,i=1,…,5, into (47) and take the inverse Laplace transformation to get the nth-approximating solution for the model (2). Mathematica software packages have been used to carry out numerical operations and implement these steps, and good results have been obtained, which are shown in the tables and graphics that follow.

5.   Numerical results

Mathematica software packages have been used to carry out numerical operations and implement the LRPSM steps, and good results have been obtained, which are shown in the following tables and graphics. We consider the parameters values as:

where the diffusion terms have been ignored. The initial conditions are considered to be:

Tables 1–5 show the constructed 12th-approximating solutions with the consecutive's errors Cons. error=|ϕi12(t,x)−ϕi10(t,x)| for i = 1, …, 5, while the fractional derivative order α was considered at different values, namely, α = 0.7, 0.8 and 0.9.

In Tables 6–10, we present the obtained approximating solution, absolute error and relative error at derivative order α = 1. The obtained LRPS approximating solutions were compared with the results introduced in [1] to get the absolute error and relative error

For more demonstration, we present Figure 1 that shows the LRPS approximating solutions at different fractional derivative orders. In addition, we show the comparative results obtained at different approximate levels, namely, n = 100 and n = 12.

6.   Discussion

The present study investigates the dynamics of a fractional diffusive cancer model incorporating virotherapy, a promising approach that utilizes viruses to combat cancer cells. The model utilizes the Caputo fractional derivative to describe the non-local dynamics of the cancer cells and the virotherapy treatment. The study employs analytical techniques and simulation methods to gain insights into the behavior of the model under various scenarios. The findings shed light on the potential efficacy of virotherapy in controlling cancer growth and highlight the significance of fractional calculus in modeling complex biological systems.

Initially, the present study focused on and refined the model of oncolytic virus replication. This model investigates the intricate relationship between tumor growth, virotherapy's infection dynamics, immune response and their collective impact on both infected and uninfected cells. To achieve a comprehensive understanding, we have precisely considered five distinct populations within the system as speculated in [1]. We proposed a modification (1) to the model by incorporating the Caputo fractional derivative in place of the traditional time derivative. This alteration leads to a difference in dimensions on the sides of the equations in the model (1). Consequently, to address this disparity, we introduce an auxiliary parameter σ, as referenced in [3] for diffusive cancer models with virotherapy.

The results from equation [9] show the uniqueness of our model for virotherapy's infection dynamics. This section interprets the findings from both the analytical stability analysis and the numerical simulations. It elaborates on the implications of the results for virotherapy as a potential cancer treatment strategy and discusses the potential challenges and limitations of the model in capturing the complexity of cancer dynamics. Furthermore, we also analyzed the LRPSM to derive an analytical solution for the fractional diffusive cancer model (2). Numerical operations using Mathematica software packages were used to implement the LRPSM steps, and we obtained very significant results. The results are demonstrated in Tables 1–10 and Figures 1a–1e.

The results in Tables 1–10 and Figures 1a–1e showed that during virotherapy, the number of infected cells reduces but the number of healthy cells also decrease. But according to our results, if there is a high concentration of Za, there will be an increase of healthy cells. For centuries, researchers have extensively explored the immune response against cancer cells. It has been established that the immune system actively surveils the body to identify and target abnormal cells. In our investigation, we focused on the diffusion coefficient of the immune system, denoted as d5, in the context of cancer treatment. To study the effect of the treatment on cancer growth, we employed the tanh-expansion method to obtain solutions for model (6). The solution, represented by equation (23), reveals that the diffusion coefficient d5 is influenced by d1, d2, d3 and d4. Specifically, we observed that the diffusion of activated immune system cells (d5) increases when the diffusion coefficients of the virus (d2) or infected cells (d3) decrease. Conversely, the diffusion of activated immune system cells (d5) decreases as the diffusion of the virus (d2) or infected cells (d3) increases. These findings shed light on the intricate interactions between different components within the system and provide valuable insights into the role of diffusion coefficients in the treatment's impact on cancer growth and immune response.

Based on the results obtained from the analytical simulation of the diffusive cancer model with virotherapy using the Caputo operator, we predict that this approach can potentially help reduce the side effects of virotherapy. The use of the Caputo fractional derivative allows for a more accurate representation of the underlying dynamics, which could lead to improved treatment strategies and better management of the therapy's impact on both infected and uninfected cells. By gaining deeper insights into the system's behavior through analytical methods, it may be possible to optimize the treatment and enhance its efficacy while minimizing undesirable side effects, thereby advancing the field of virotherapy for cancer treatment. The combination of analytical techniques and numerical simulations enhances our understanding of the system's behavior, paving the way for more effective and targeted cancer treatments in the future. The study's findings contribute to the growing field of fractional calculus applications in modeling complex biological systems and hold promise for advancements in cancer therapy.

7.   Conclusions

We have studied the existence and uniqueness of the solution to this proposed and important model by presenting the main theories in this work. In addition, the LRPSM was applied to derive approximate solutions for the presented fractional model. The steps of this method were implemented using Mathematica software packages and good results were extracted. The results obtained were presented through tables showing the approximate solutions of the model. Moreover, a comparison of the derived solutions with other results presented in the literature was introduced to evaluate the efficiency of this method and its smoothness of use to solve such systems. As a future work, this fractional model can be studied by considering another fractional operator, such as the Atangana-Baleanu fractional operator, and comparing the results to clarify the extent to which the results are affected by changing the fractional operator. Also, other numerical techniques can be used to derive approximate solutions for this model and compare them with the results presented in this work, considering different values of the parameters. To provide a thorough understanding of cancer dynamics, future research may investigate more intricate and realistic cancer models, taking into account elements like tumor heterogeneity, geographical impacts and varying immune responses. The mathematical model can be improved through experimental validation, which includes observations of cancer cell and immune system behavior in response to virotherapy. Performing a sensitivity analysis to look at parameter variations might give insights into how robust the results are and help pinpoint important elements affecting results. Through clinical trials and case studies, it is possible to examine the clinical applicability of virotherapy for particular cancer types and determine its viability and effectiveness in actual treatment situations.

Use of AI tools declaration

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

Conflicts of interest

The authors declare that they have no conflicts of interest.

Authors' contributions

All authors conceived of the study, participated in its design and coordination, drafted the manuscript, participated in the sequence alignment and read and approved the final manuscript.