Solving a fractional order model of HIV infection of CD4+ T cells

Solving a fractional order model of HIV infection of CD4+ T cells

Mathematical and Computer Modelling 54 (2011) 2132–2138 Contents lists available at ScienceDirect Mathematical and Computer Modelling journal homepa...

560KB Sizes 3 Downloads 38 Views

Mathematical and Computer Modelling 54 (2011) 2132–2138

Contents lists available at ScienceDirect

Mathematical and Computer Modelling journal homepage: www.elsevier.com/locate/mcm

Solving a fractional order model of HIV infection of CD4+ T cells Ahmet Gökdoˇgan a , Ahmet Yildirim b , Mehmet Merdan a,∗ a

Gümüşhane University, Department of Mathematics Engineering, 29100-Gümüşhane, Turkey

b

Ege University, Department of Mathematics, 35000 İzmir, Turkey

article

abstract

info

Article history: Received 27 February 2011 Received in revised form 12 May 2011 Accepted 12 May 2011

In this paper, a multi-step differential transform method (MsDTM) is performed to give approximate and analytical solutions of nonlinear fractional order ordinary differential equation systems such as a model for HIV infection of CD4+ T cells. The numerical solutions obtained from the proposed method indicate that the approach is easy to implement and accurate when applied to systems of fractional differential equations. Some figures are presented to show the reliability and simplicity of the methods. © 2011 Elsevier Ltd. All rights reserved.

Keywords: A model for HIV infection of CD4+ T cells Multi-step differential transform method Fractional order ordinary differential equation systems

1. Introduction In this paper, dynamics of a model for HIV infection of CD4+ T cells is examined [1]. The components of the basic three component models are the concentration of susceptible CD4+ T cells. CD4+ T cells infected by the HIV viruses and free HIV virus particles in the blood are denoted respectively by T (t ), I (t ) and V (t ). CD4+ T cells are also named as leukocytes or T helper cells. These with order cells in human immunity systems fight against diseases. HIV uses cells in order to propagate. 800 The number of CD4+ T cells of a healthy person is 1200 mm3 . Here, a fractional order model of for HIV infection of CD4+ T cells is reviewed. These quantities satisfy α1

D T = q − α T + rT Dα2 I = k∗ VT − β I Dα3 V = N β I − γ V

 1−

T +I Tmax



− k∗ VT (1)

with the initial conditions: T (0) = r1 , I (0) = r2 and V (0) = r3 . Throughout this paper, we set q = 0.1, α = 0.02, β = 0.3, r = 3, γ = 2.4, k∗ = 0.0027, N = 10, Tmax = 1500. The logistic growth of the healthy CD4+ T cells is now described by (1 − TT +I ), and proliferation of infected CD4+ T cells is neglected. The term kVT describes the incidence of HIV infection of max

healthy CD4+ T cells, where k > 0 is the infection rate. Each infected CD4+ T cell is assumed to produce N virus particles during its lifetime, including any of its daughter cells. The body is believed to produce CD4+ T cells from precursors in the bone marrow and thymus at a constant rate q. When stimulated by antigen or mitogen, T cells multiply through mitosis with a rate r. Tmax is the maximum level of CD4+ T cell concentration in the body [2–5]. Differential equations of fractional order have been the subject of many studies owing to their implementation in various applications in fluid mechanics, viscoelasticity, biology, physics and engineering. Lately, a large amount of studies developed concerning the application of fractional differential equations in nonlinear dynamics [6–11]. As most fractional differential



Corresponding author. E-mail addresses: [email protected] (A. Gökdoˇgan), [email protected] (A. Yildirim), [email protected], [email protected] (M. Merdan). 0895-7177/$ – see front matter © 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.mcm.2011.05.022

A. Gökdoˇgan et al. / Mathematical and Computer Modelling 54 (2011) 2132–2138

2133

equations do not have the exact analytical solutions, the present solutions are obtained from approximation and numerical techniques used extensively. Recently, the Adomian decomposition method [12–17], homotopy perturbation method [18–29], homotopy analysis method [30] and variational iteration method [31–33] have been applied for solving a wide range of problems. The differential transform method (DTM) is a numerical as well as analytical method for solving integral equations, ordinary, partial differential equations and differential equation systems. The method provides the solution in terms of convergent series with easily computable components. The concept of the differential transform was first proposed by Zhou [34] and its main application concern with both linear and nonlinear initial value problems in electrical circuit analysis. The DTM gives exact values of the nth derivative of an analytic function at a point in terms of known and unknown boundary conditions in a fast manner. This method constructs, for differential equations, an analytical solution in the form of a polynomial. It is different from the traditional high order Taylor series method, which requires symbolic computations of the necessary derivatives of the data functions. The Taylor series method is computationally taken long time for large orders. The DTM is an iterative procedure for obtaining analytic Taylor series solutions of differential equations. Recently, the application of differential transformation method is successfully extended to obtain analytical approximate solutions to linear and nonlinear differential equations of fractional order [35–37]. However, DTM has some drawbacks. By using the DTM, we obtain a series solution, actually a truncated series solution. This series solution does not exhibit the real behaviors of the problem but gives a good approximation to the true solution in a very small region. To overcome the shortcoming, MsDTM was presented in [38–40]. Lately, Ma and Fan [41] analyzed a specific sub-class of N-soliton solutions, formed by linear combinations of exponential traveling waves, for Hirota bilinear equations. The starting point is to solve a system of nonlinear algebraic equations for the wave related numbers, called the N-wave solution condition. The resulting system tells what Hirota bilinear equations the linear superposition principle of exponential waves will apply to. Higher dimensional Hirota bilinear equations have a better opportunity to satisfy the N-wave solution condition since there are more parameters to choose from. Applications were made for the 3 + 1 dimensional KP, Jimbo–Miwa and BKP equations, thereby presenting their particular exact multiple wave solutions. Moreover, an opposite procedure was proposed for conversely generating Hirota bilinear equations possessing multiple wave solutions formed by linear combinations of exponential waves, and a few such examples were computed. Also, a direct and systematic solution procedure for constructing multiple wave solutions to nonlinear partial differential equations is proposed [42]. The presented method is oriented toward the ease of use and capability of computer algebra systems, allowing us to carry out the involved computations conveniently through powerful computer algebra systems. It is by the use of computer algebra systems that in each of the cases of two-wave and three-wave solutions. The method can also be easily applied to other nonlinear evolution and wave equations in mathematical physics. We cannot use exp-function or other traveling wave solution method for fractional problems. Because fractional derivatives include complex integration computations. This complexity does not allow us to use traveling wave solution methods. The aim of this paper is to extend the application of the multi-step differential transformation method [43] to solve a fractional order model for HIV infection of CD4+ T cells (1).

2. Fractional differential transform method Consider a general system of fractional differential equations: Dα∗ 1 x1 (t ) + h1 (t , x1 , x2 , . . . , xm ) = g1 (t ), Dα∗ 2 x2 (t ) + h2 (t , x1 , x2 , . . . , xm ) = g2 (t ),

.. .

(2)

Dα∗ m xm (t ) + hm (t , x1 , x2 , . . . , xm ) = gm (t ), α

where D∗ i is the derivative of xi of order αi in the sense of Caputo and 0 < αi ≤ 1, subject to the initial conditions x1 (t0 ) = d1 ,

x2 (t0 ) = d2 , . . . , xm (t0 ) = dm .

(3)

In this paper, we introduce the multi-step fractional differential transform method used in this paper to obtain approximate analytical solutions for the system of fractional differential equations (1). This method has been developed in [36] as follows: Dqx0 f

(x) =

1

dm

Γ (m − q) dxm

[∫

x x0

f (t )

(x − t )1+q−m

]

dt ,

(4)

for m − 1 ≤ q < m, m ∈ Z + , x > x0 . Let us expand the analytical and continuous function f (x) in terms of fractional power series as follows:

2134

A. Gökdoˇgan et al. / Mathematical and Computer Modelling 54 (2011) 2132–2138

Fig. 1. Plots of approximations obtained 5 term-MsDTM with dt = 0.001 (a1)–(b1)–(c1) for α1 = α2 = α3 = 1, (a2)–(b2)–(c2) for α1 = α2 = α3 = 0.9, (a3)–(b3)–(c3) for α1 = α2 = α3 = 0.8, (a4)–(b4)–(c4) for α1 = α2 = α3 = 0.7 and (a5)–(b5)–(c5) for α1 = α2 = α3 = 0.6.

f (x) =

∞ −

k

F (k)(x − x0 ) α ,

k=0

where α is the order of fraction and F (k) is the fractional differential transform of f (x).

(5)

A. Gökdoˇgan et al. / Mathematical and Computer Modelling 54 (2011) 2132–2138

2135

In order to avoid fractional initial and boundary conditions, we define the fractional derivative in the Caputo sense. The relation between the Riemann–Liouville operator and Caputo operator is given by

 D∗x0 f (x) = q

Setting f (x) − as follows:

 − 1 k (k) f (x) − (x − x0 ) f (x0 ) . k! k=0

Dqx0

∑m−1

1 k=0 k!

m−1

(6)

(x − x0 )k f (k) (x0 ) in Eq. (4) and using Eq. (6), we obtain fractional derivative in the Caputo sense [44]

   m −1 ∑   1 k (k)   f (t ) −  (t − x0 ) f (x0 )    k! m ∫ x  1 d   k=0 Dq∗x0 f (x) = dt     Γ (m − q) dxm  x0  (x − t )1−q−m       

(7)

since the initial conditions are implemented are implemented to the integer order derivatives, the transformation of the initial conditions are defined as follows:

 k  +   if α ∈ Z , F (k) =    if k ̸∈ Z + α

1



k

d α f ( x)

k ! α



k

dx α

x =x 0

for k = 0, 1, 2, . . . , (qα − 1)

(8)

0,

where, q is the order of fractional differential equation considered. The following theorems that can be deduced from Eqs. (4) and (5) are given below, for proofs and detailed see [36]. Theorem 1. If z (t ) = x(t ) ± y(t ), then Z (k) = X (k) ± Y (k). Theorem 2. If z (t ) = cy(t ), then Z (k) = cY (k). Theorem 3. If z (t ) = x(t )y(t ), then Z (k) =

∑k

k1 =0

X (k1 )Y (k − k1 ).

Theorem 4. If z (t ) = (t − t0 )n , then Z (k) = δ(k − α p) where,

δ(k) =



1

k=0

0

k ̸= 0. q

Theorem 5. If z (t ) = Dt0 [g (t )] then Z (k) =

Γ (q+1+ αk ) Γ (1+ αk )

G(k + α q).

According to fractional DTM, by taking differential transformed both sides of the systems of equations given Eqs. (2) and (3) is transformed as follows:

Γ (α1 (k + 1) + 1) X1 (k + 1) + H1 (k) = G1 (k), Γ (α1 k + 1) Γ (α2 (k + 1) + 1) X2 (k + 1) + H2 (k) = G2 (k), Γ (α2 k + 1) .. .

(9)

Γ (αm (k + 1) + 1) Xm (k + 1) + Hm (k) = Gm (k), Γ (αm k + 1) X1 (0) = d1 ,

X2 (0) = d2 , . . . ,

Xm (0) = dm .

(10)

Therefore, according to DTM the N-term approximations for (2) can be expressed as

ϕ1,n (t ) = x1 (t ) =

N −

X1 (k)t k ,

k=1 N

ϕ2,n (t ) = x2 (t ) =



.. . ϕm,n (t ) = xm (t ) =

X2 (k)t k ,

k=1

N − k=1

Xm (k)t k .

(11)

2136

A. Gökdoˇgan et al. / Mathematical and Computer Modelling 54 (2011) 2132–2138

Fig. 2. Phase portraits obtained 5 term-MsDTM with dt = 0.001 (d1)–(e1)–(f1) for α1 = α2 = α3 = 1, (d2)–(e2)–(f2) for α1 = α2 = α3 = 0.9, (d3)–(e3)–(f3) for α1 = α2 = α3 = 0.8, (d4)–(e4)–(f4) for α1 = α2 = α3 = 0.7 and (d5)–(e5)–(f5) for α1 = α2 = α3 = 0.6.

3. Multi-step fractional differential transform method The approximate solutions (2) are generally, as will be shown in the numerical experiments of this paper, not valid for large t. A simple way of ensuring validity of the approximations for large t is to treat (9)–(10) as an algorithm for approximating the solutions of (2)–(3) in a sequence of intervals choosing the initial approximations as

A. Gökdoˇgan et al. / Mathematical and Computer Modelling 54 (2011) 2132–2138

2137

x1,0 (t ) = x1 (t ∗ ) = d∗1 , x2,0 (t ) = x2 (t ∗ ) = d∗2 ,

.. .

(12)

xm,0 (t ) = xm (t ∗ ) = d∗m . In order to carry out the iterations in every sub-interval [0, t1 ), [t1 , t2 ), [t2 , t3 ), . . . , [tj−1 , t ) of equal length h, we would need to know the values of the following [43], x∗1,0 (t ) = x1 (t ∗ ),

x∗2,0 (t ) = x2 (t ∗ ),

...,

x∗m,0 (t ) = xm (t ∗ ).

(13) ∗

But, in general, we do not have these information at our clearance except at the initial point t = t0 . A simple way for obtaining the necessary values could be by means of the previous n-term approximations ϕ1,n , ϕ2,n , . . . , ϕm,n of the preceding sub-interval, i.e., x∗1,0 ∼ = ϕ1,n (t ∗ ),

x∗2,0 ∼ = ϕ2,n (t ∗ ), . . . ,

x∗m,0 ∼ = ϕm,n (t ∗ ).

(14)

4. Numerical results We will apply classic DTM and the adaptive MsDTM to nonlinear fractional order ordinary differential equation systems (1). Applying classic DTM for Eq. (1)

Γ (α1 k + 1) T (k + 1) = − Γ (α1 (k + 1) + 1)



r Tmax

k −

 qδ(k) + (r − α)T (k) − k −

T (k1 )I (k − k1 ) − k



k 1 =0

k −

r Tmax

T (k1 )T (k − k1 )

k1 =0

 V (k1 )T (k − k1 )

k1 =0

Γ (α2 k + 1) I (k + 1) = Γ (α2 (k + 1) + 1)



k −



k

(15)

 V (k1 )T (k − k1 ) − β I (k)

k1 =0

Γ (α3 k + 1) v(k + 1) = {N β I (k) − γ V (k)} Γ (α3 (k + 1) + 1) T (0) = 0.05, I (0) = 0.1, V (0) = 0.5. By applying the multi-step DTM to Eq. (15), Ti (n), Ii (n), Vi (n) satisfy the following recurrence relations for n = 1, 2, . . . , N − 1.

Γ (α1 k + 1) Ti (k + 1) = Γ (α1 (k + 1) + 1)



r Tmax

k −

 qδ(k) + (r − α)Ti (k) −

Ti (k1 )Ii (k − k1 ) − k



k1 =0

Γ (α2 k + 1) Ii (k + 1) = Γ (α2 (k + 1) + 1)

k −

r Tmax

k −

Ti (k1 )Ti (k − k1 )

k1 =0

 Vi (k1 )Ti (k − k1 )

k1 =0

 k



k −

 Vi (k1 )Ti (k − k1 ) − β Ii (k)

(16)

k1 =0

Γ (α3 k + 1) {N β Ii (k) − γ Vi (k)} Γ (α3 (k + 1) + 1) T0 (0) = 0.05, I0 (0) = 0.1, V0 (0) = 0.5. Ti (0) = Ti−1 (ti ), Ii (0) = Ii−1 (ti ), Vi (0) = Ti−1 (ti ),

vi (k + 1) =

i = 1, 2, . . . , M .

Fig. 1 shows the approximate solutions for system (1) obtained for the different values of αi . Fig. 1 indicates periodic behavior of the approximate solutions for system (1) obtained for the values of α1 = α2 = α3 = 1, 0.9, 0.8, 0.7 and 0.6. A number of oscillation may be obtained by increasing and then decreasing αi values of fractional order as in the system (1) above. In Fig. 2, T − I , T − V , I − V and T − I − V phase portraits obtained using 5-term MsDTM solutions are given. Consider the trajectory passing through the point (200, 1200) in Fig. 2(d(1)–d(5)). At the point the ratio of CD4+ T cells infected by the HIV viruses to the concentration of susceptible CD4+ T cells is high; as a result of the population of CD4+ T cells infected by the HIV viruses drops. The ratio of I (t ) to T (t ) drops, and so the population of T (t ) increases. The resulting cyclic behavior is repeated over and over and is shown as the largest closed trajectory in Fig. 2.

2138

A. Gökdoˇgan et al. / Mathematical and Computer Modelling 54 (2011) 2132–2138

5. Conclusions In this paper, MsDTM method was used for finding approximate analytical solutions of nonlinear ordinary differential equation systems of fractional order such as a model for HIV infection of CD4+ T cells. We demonstrated the accuracy and efficiency of these methods by solving some ordinary differential equation systems of fractional order. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44]

W. Liancheng, Y.L. Michael, Mathematical analysis of the global dynamics of a model for HIV infection of CD4+ T cells, Math. Biosci. 200 (2006) 44–57. A.S. Perelson, P.W. Nelson, Mathematical analysis of HIV — I dynamics in vivo, SIAM Rev. 41 (1) (1999) 3–44. L. Wang, M.Y. Li, Mathematical analysis of the global dynamics of a model for HIV infection of CD4+ T cells, Math. Biosci. 200 (2006) 44–57. B. Asquith, C.R.M. Bangham, The dynamics of T-cell fratricide: application of a robust approach to mathematical modelling in immunology, J. Theoret. Biol. 222 (2003) 53–69. M. Nowak, R. May, Mathematical biology of HIV infections: antigenic variation and diversity threshold, Math. Biosci. 106 (1991) 1–21. V. Daftardar-Gejji, H. Jafari, An iterative method for solving nonlinear functional equations, J. Math. Anal. Appl. 316 (2) (2006) 753–763. J.H. He, Approximate analytical solution for seepage flow with fractional derivatives in porous media, Comput. Methods Appl. Mech. Engng. 167 (1998) 57–68. S. Momani, An explicit and numerical solutions of the fractional KdV equation, Math. Comput. Simulation 70 (2) (2005) 110–118. S. Momani, Non-perturbative analytical solutions of the space- and time-fractional Burgers equations, Chaos Solitons Fractals 28 (4) (2006) 930–937. S. Momani, Z. Odibat, Analytical solution of a time-fractional Navier–Stokes equation by Adomian decomposition method, Appl. Math. Comput. 177 (2) (2006) 488–494. S. Momani, Z. Odibat, Analytical approach to linear fractional partial differential equations arising in fluid mechanics, Phys. Lett. A 355 (2006) 271–279. G. Adomian, Nonlinear stochastic differential equations, J. Math. Anal. Appl. 55 (1976) 441–452. G. Adomian, A review of the decomposition method in applied mathematics, J. Math. Anal. Appl. 135 (1988) 501–544. G. Adomian, Solving Frontier Problems of Physics: The Decomposition Method, Kluewr Academic Publishers, Boston, 1994. N.T. Shawagfeh, The decomposition method for fractional differential equations, J. Fract. Calc. 16 (1999) 27–33. A.M. Wazwaz, The modified decomposition method for analytic treatment of differential equations, Appl. Math. Comput. 173 (2006) 165–176. S. Momani, N.T. Shawagfeh, Decomposition method for solving fractional Riccati differential equations, Appl. Math. Comput. 182 (2006) 1083–1092. S. Abbasbandy, Homotopy perturbation method for quadratic Riccati differential equation and comparison with Adomian’s decomposition method, Appl. Math. Comput. 173 (2006) 493–500. Z. Odibat, S. Momani, A reliable treatment of homotopy perturbation method for Klein–Gordan equations, Phys. Lett. A 365 (5–6) (2007) 351–357. Z. Odibat, S. Momani, Modified homotopy perturbation method: application to quadratic Riccati differential equation of fractional order, Chaos Solitons Fractals 36 (1) (2008) 167–174. S. Momani, Z. Odibat, Numerical comparison of methods for solving linear differential equations of fractional order, Chaos Solitons Fractals 31 (2007) 1248–1255. Z. Odibat, S. Momani, Numerical methods for nonlinear partial differential equations of fractional order, Appl. Math. Modelling 32 (1) (2008) 28–39. S. Momani, Z. Odibat, Numerical approach to differential equations of fractional order, J. Comput. Appl. Math. 207 (2007) 96–110. Z. Odibat, A new modification of the homotopy perturbation method for linear and nonlinear operators, Appl. Math. Comput. 189 (1) (2007) 746–753. J.H. He, Homotopy perturbation technique, Comput. Methods Appl. Eng. 178 (1999) 257–262. J.H. He, Application of homotopy perturbation method to nonlinear wave equations, Chaos Solitons Fractals 26 (3) (2005) 695–700. J.H. He, The homotopy perturbation method for nonlinear oscillators with discontinuities, Appl. Math. Comput. 151 (2004) 287–292. J.H. He, Homotopy perturbation method for bifurcation of nonlinear problems, Int. J. Nonlinear Sci. Numer. Simul. 6 (2) (2005) 207–208. M. Merdan, Homotopy perturbation method for solving a model for hiv infection of CD4+ T cells, Istanb. Commerce Uni. J. Sci. 12 (2007) 39–52. A.E. Ajou, Z. Odibat, S. Momani, A. Alawneh, Construction of analytical solutions to fractional differential equations using homotopy analysis method, IAENG Int. J. Appl. Math. 40 (2) (2010) IJAM_40_2_01. Z. Odibat, S. Momani, Application of variation iteration method to nonlinear differential equations of fractional order, Int. J. Nonlinear Sci. Numer. Simul. 1 (7) (2006) 15–27. J.H. He, Variational iteration method — a kind of non-linear analytic technique: some examples, Int. J. Nonlinear Mech. 34 (1999) 699–708. J.H. He, Variational iteration method for autonomous ordinary differential systems, Appl. Math. Comput. 114 (2000) 115–123. J.K. Zhou, Differential Transformation and its Applications for Electrical Circuits, Huazhong University Press, Wuhan, China, 1986 (in Chinese). V.S. Ertürk, S. Momani, Solving systems of fractional differential equations using differential transform method, J. Comput. Appl. Math. 215 (2008) 142–151. A. Arikoglu, I. Ozkol, Solution of fractional differential equations by using differential transform method, Chaos Solitons Fractals 34 (5) (2007) 1473–1481. A. Al-rabtah, V.S. Ertürk, S. Momani, Solutions of a fractional oscillator by using differential transform method, Comput. Math. Appl. 59 (2010) 1356–1362. M.J. Jang, C.L. Chen, Y.C. Liu, On solving the initial value problems using the differential transformation method, Appl. Math. Comput. 115 (2000) 145–160. M.J. Jang, C.L. Chen, Y.C. Liu, Analysis of the response of a strongly nonlinear damped system using a differential transformation technique, Appl. Math. Comput. 88 (1997) 137–151. I.H. Abdel-Halim Hassan, On solving some eigenvalue problems by using a differential Transformation, Appl. Math. Comput. 127 (2002) 1–22. W.X. Ma, E. Fan, Linear superposition principle applying to Hirota bilinear equations, Comput. Math. Appl. 61 (2011) 950–959. W.X. Ma, T. Huang, Y. Zhang, A multiple exp-function method for nonlinear differential equations and its application, Phys. Scr. 82 (2010) 065003 (8 pp.). Z.M. Odibat, C. Bertelle, M.A. Aziz-Alaoui, G.H.E. Duchamp, A multi-step differential transform method and application to non-chaotic or chaotic systems, Comput. Math. Appl. 59 (2010) 1462–1472. M. Caputo, Linear models of dissipation whose Q is almost frequency independent II, Geophys. J. Roy. Astronom. Soc. 13 (1967) 529–539.