# Analysis of Fractional Order Mathematical Model of Hematopoietic Stem Cell Gene-Based Therapy.

1. IntroductionAn immune system is body's primary defense system that has a function to fight against microbes. Humoral immunity consists of innate immunity and adaptive immunity. The main principle of innate immunity is an initial defense which responds quickly against microbes. Meanwhile, adaptive immunity would be activated when innate immunity failed to eradicate microbes. Adaptive immunity consists of B lymphocytes cells (B cells) and T lymphocytes cells (T cells). T cells consist of CD4+ T cells and cytotoxic T lymphocytes (CTL and CD8+ T cells). Adaptive immunity has some capabilities such as specificity, diversity, and memory, so it can eradicate microbes effectively [1].

Viruses are one of the microorganisms that infect the immune system. A virus will be detected and then eradicated by innate immunity. However, there are several viruses that can escape from innate immunity and infect CD4+ T cells. Infected CD4+ T cells will produce cytokines to stimulate proliferation and differentiation of precursor CTL cells to become effector CTL cells that would eradicate infected CD4+ T cells [2].

Hematopoiesis is the process to derive blood cells that are located in bone marrow. The first stage of hematopoietic is a stem cell. Hematopoietic stem cells (HSCs) have an ability to multipotency and self-renewing. While in embryonic phase, HSCs that migrate to thymus would be differentiated to mature T cells [1]. Gene-based therapy has been discussed to cure immune impairing infections, such as HIV. The therapy's concept is specifying technic and immune responses to face viruses by using genetically modified HSCs. In 2012, Kitchen et al. performed an in vivo experiment with HSCs gene-based therapy. Based on the experiment, engineered HSCs have the ability to establish functional antivirus responses that suppressed HIV replication. Hence, it will allow suppression of infected CD4+ T cells and prevent suppression of uninfected CD4+ T cells [3].

Korpusik [4] developed a mathematical model to study the influence of a constant influx of CTLs on the dynamic of the virus and immune system interactions, using a modification of the basic mathematical model of virus-induced impairment of help [2]. In addition, Korpusik and Kolev [5] developed a mathematical model to study the dynamical virus-immune system interactions after a single injection of CD8+ T cells derived from HSCs.

Ordinary differential equation (ODE) system that consists of first-order differential equations can be generalized to fractional differential equation (FDE) system that consists of fractional order differential equations a, with a fractional order parameter 0 < [alpha] [less than or equal to] 1 [6]. In most biological systems, FDE are naturally connected to systems with memory [7]. In addition, memory effect has an important role in the disease spread. The presence of memory effects on past events will affect the disease spread in the future. The distance of memory effect indicates the history of disease spread. Thus, memory effects on the spread of an infectious disease can be investigated using fractional derivatives [8-13].

In this paper, we extend the ODE model from Korpusik and Kolev [14] into a fractional differential equation system model. We also determined equilibria and stability of the equilibria from the proposed model. Finally, we perform numerical simulation to support the mathematical model interpretation.

2. Mathematical Model Formulation

In this section, we proposed a mathematical model of HSCs gene-based therapy based on [5]. The model was constructed under the following assumptions:

(1) The model consists of five compartments, namely, concentration of uninfected CD4+ T cells (X), concentration of infected CD4+ T cells (Y), concentration of precursor CTL (W), concentration of effector CTL (Z), and concentration of CD8+ T lymphocytes derived from HSCs (H).

(2) Uninfected CD4+ T cells are produced with constant rate.

(3) Viral infections occur only in a human body.

(4) Free virus particle only infected uninfected CD4+ T cells.

(5) Effector CTL only killed infected CD4+ T cells.

(6) Injected HSCs will be differentiated into a precursor CTL.

(7) The proliferation rate of precursor CTL cells depends on susceptible CD4+ T cells, infected CD4+ T cells, and CTL precursor cells at the current time.

The transmission diagram of the model is shown in Figure 1.

The basic model from Korpusik and Kolev [5] was given as follows:

dX/dt = [LAMBDA] - [theta]X - [beta]XY (1a)

dY/dt = [beta]XY - aY - pYZ (1b)

dW/dt = m[gamma]H + cWXY - qWY - [b.sub.1]W (1c)

dZ/dt = qWY - [b.sub.2]Z (1d)

dH/dt = [gamma]H (1e)

Here, all parameters [LAMBDA], [theta], [beta], a, p, c, [b.sub.1], [b.sub.2], [gamma] > 0, 0 [less than or equal to] m, q [less than or equal to] 1 and X(0),Y(0),W(0),Z(0),H(0) [greater than or equal to] 0. The description of the parameter for the model is presented in Table 1.

The differential equation (1a) describes the dynamic of uninfected CD4+ T cells concentration. Concentrations of CD4+ T cells increase by production of bone marrow and decreased caused by natural death rate and infected by free virus particle. The differential equation (1b) shows the dynamic of infected CD4+ T cells concentrations. Concentrations of infected CD4+ T cells increase, caused by infection of free virus particle. Infected CD4+ T cells will be decreased because of natural death rate and it was killed by effector CTL cells.

The differential equation (1c) describes the dynamic of precursor CTL cells concentrations. Concentrations of precursor CTL cells increased by CTL proliferation process and HSCs were differentiated into CD8+ T cells that are parts of precursor CTL cells. Concentrations of precursor CTL cells are decreased, caused by differentiating into effector CTL, and decreased by natural death rate. Equation (1d) shows the dynamic of concentrations of effector CTL cells at time. Concentrations of effector CTL cells are increased, caused by differentiating of precursor CTL into effector CTL phase, and decreased, caused by natural death rate. The last equation (1e) shows the dynamic of CD8+ T cells concentrations that derived from HSCs. Concentrations of CD8+ T cells are decreased, caused by successfully passing the thymic selection and differentiating into functional precursor CTL cells.

Next, we consider a fractional order model of (1a)-(1e) above. The fractional model corresponding to system (1a)-(1e) is as follows:

[mathematical expression not reproducible] (2)

where the fractional order derivative 0 < [alpha] [less than or equal to] 1. Fractional derivative of the model (2) is adopted from Caputo's definition. The main advantages of Caputo approach are the initial values for fractional differential equations with the Caputo derivatives taking on the same form as for integer order differential equations [15]. The Caputo fractional derivative is defined as follows.

Definition 1 (see [15]). Let [alpha] > 0, t > 0, and n [member of] N. Caputo fractional derivative [D.sup.[alpha]] [??] d[t.sup.[alpha]]/d[t.sup.[alpha]], with fractional order [alpha], of function f(t) is defined by

[mathematical expression not reproducible] (3)

where [GAMMA](*) is the gamma function.

3. Stability Analysis

In this section, we study stability of the equilibriums of the fractional order model (2) above. We begin by computed the basic reproduction number [R.sub.0] of model (2). Basic reproduction number [R.sub.0] is defined as the average number of new cases of an infection caused by one typical infected individual, in a population consisting of susceptible only [16]. If [R.sub.0] < 1, then the infections will die out, while if [R.sub.0] > 1,then there is an epidemic case [17].

The stability theorem on fractional order system was given in the following theorem.

Theorem 2 (see [18]). Consider a nonlinear fractional order system

[D.sup.[alpha]] x (t) = f (x) (4)

where 0 < [alpha] [less than or equal to] 1, x [member of] [R.sup.n], and f [member of] [R.sup.n]. The equilibrium points [x.sup.*] of system (4) are calculated by solving equation f(x) = 0. The equilibrium points [x.sup.*] are locally asymptotically stable if all the eigenvalues [[lambda].sub.j] (j = 1,2, ..., n) of the Jacobian matrix A = [partial derivative]f/[partial derivative]x evaluated at the equilibrium points [x.sup.*] satisfy the following condition:

[absolute value of (arg [[lambda].sub.j])] > [alpha][pi]/2. (5)

The equilibriums of the fractional mathematical model in (2) satisfy the following equations:

[LAMBDA] - [theta]X - [beta]XY = 0 (6)

[beta]XY - aY - pYZ = 0 (7)

m[gamma]H + cWXY - qWY - [b.sub.1]W = 0 (8)

qWY - [b.sub.2]Z = 0 (9)

-[gamma]H = 0 (10)

The fractional model in (2) has three equilibriums, namely, virus-free equilibrium ([E.sub.0]), CTL-Exhaustion Equilibrium ([E.sub.1]), and control immune equilibrium [E.sub.2]. The virus-free equilibrium is condition when there is CD4+ T cells in human body (X [not equal to] 0), no infected CD4+ T cells (Y = 0), and CTL precursor is not activated (W = 0). The virus-free equilibrium of model (2) above is given by [E.sub.0] = ([X.sub.0],[Y.sub.0],[W.sub.0],[Z.sub.0],[H.sub.0]) = ([LAMBDA]/[theta], 0, 0, 0, 0).

Basic reproduction number ([R.sub.0]) is an important parameter in epidemiological cases. Basic reproduction number [R.sub.0] is defined by average secondary infections caused by one primary infection in susceptible population. In this paper, we use Next Generation Matrix (NGM) that developed by [19] to determine [R.sub.0]. By using the Next Generation method, we obtained the basic reproduction number [R.sub.01], = [beta][LAMBDA]/[theta]a. Stability of virusfree equilibrium is presented in Theorem 3.

Theorem 3. The virus-free equilibrium [E.sub.0] of model (2) is locally asymptotically stable if and only if [R.sub.01], < 1. Proof. The Jacobi matrix that evaluated at the virus-free equilibrium ([E.sub.0]) of the model (2) above is as follows:

[mathematical expression not reproducible]. (11)

The eigenvalues of the Jacobi matrix J(E0) are A1 = -0, A2 = [beta][LAMBDA]/[theta]-a, [[lambda].sub.3] = -[b.sub.1], [[lambda].sub.4] = -[b.sub.2], and [[lambda].sub.5] = -[gamma]. Hence, we have [absolute value of (arg([[lambda].sub.1]))] = [absolute value of (arg([[lambda].sub.3]))] = [absolute value of (arg([[lambda].sub.4]))] = [absolute value of (arg([[lambda].sub.5]))] = [pi] > [pi]/2. Meanwhile, [absolute value of (arg([[lambda].sub.2]))] = [pi] > [pi]/2 if and only if [R.sub.01] = [beta][LAMBDA]/a[theta] < 1. If the condition [R.sub.01] = [beta][LAMBDA]/a[theta] < 1 is satisfied, then it is clear that all the eigenvalues satisfy the condition [absolute value of (arg [[lambda].sub.j])] > [alpha][pi]/2, for j = 1,2, 3,4, 5. Hence, the virus-free equilibrium ([E.sub.0]) of model (2) is locally asymptotically stable if [R.sub.01] < 1.

We continue with the stability analysis of the CTL-Exhaustion Equilibrium ([E.sub.1]) of the model (2). The CTL-Exhaustion Equilibrium is condition when there is CD4+ T cells in human body (X [not equal to] 0), there are infected CD4+ T cells (Y [not equal to] 0), but CTL precursor is not activated yet (W = 0). The CTL-Exhaustion Equilibrium of the model in (2) is given by [E.sub.1] = (a/[beta], ([beta][LAMBDA]-A[theta])/[beta]a, 0, 0, 0). This equilibrium will exist whenever [R.sub.01] > 1.

Theorem 4. Let q[beta]-ca < 0 and [R.sub.02] = [b.sub.1]a[[beta].sup.2]/(a[theta]-[beta]]LAMBDA])(q[beta]-ca). The CTL-Exhaustion Equilibrium ([E.sub.1]) of the model in (2) is locally asymptotically stable if and only if [R.sub.02] > 1.

Proof. The Jacobi matrix evaluated at the free virus equilibrium [E.sub.1] = (a/[beta], ([beta][LAMBDA]-a[theta])/[beta]a, 0, 0, 0) of the model in (2) above is as follows:

[mathematical expression not reproducible]. (12)

The eigenvalues of the Jacobi matrix J([E.sub.1]) above are [[lambda].sub.1] = -[gamma], [[lambda].sub.2] = -[b.sub.2], [[lambda].sub.3] = (a[theta]-[beta][LAMBDA]) ((-ca + q[beta])/a[[beta].sup.2])-[b.sub.1], and [[lambda].sub.4,5] is root of polynomial [[lambda].sup.2] + b[lambda] + c = 0, where b = [beta][LAMBDA]/[theta], and c = [beta][LAMBDA]-a[theta]. We have [absolute value of (arg([[lambda].sub.1])] = [absolute value of (arg([[lambda].sub.2])] = [pi] > [pi]/2. Then we determined the argument of [[lambda].sub.4] and [[lambda].sub.5]. Based on the condition of existence of [E.sub.1] that [R.sub.01] = [beta][LAMBDA]/a[theta] > 1, we obtained c = [beta][LAMBDA]-a[theta] > 0. By using Routh-Hurwitz criterion, we obtained [absolute value of (arg([[lambda].sub.4])] = [absolute value of (arg([[lambda].sub.5])] > [pi]/2.

Since all parameters are assumed to have positive value, the third eigenvalue [[lambda].sub.3] is a real number. We obtained that [absolute value of (arg([[lambda].sub.3])] = [pi] > [pi]/2 if and only if the condition ((a[theta]-[beta][LAMBDA]) (Q[beta]-ca) < [b.sub.1]a[[beta].sup.2]) is fulfilled. It is clear that all the eigenvalues satisfy the condition [absolute value of (arg([[lambda].sub.j])] > [alpha][pi]/2, for j = (1, 2, 3, 4, 5). Hence, the CTL-Exhaustion Equilibrium [E.sub.1] of model (2) is locally asymptotically stable if and only if (a[theta]-[beta][LAMBDA])(q[beta]-ca) < [b.sub.1]a[[beta].sup.2] or [R.sub.02] > 1. This completes the proof

Last, we analyze the stability of the control immune equilibrium ([E.sub.2]). The control immune equilibrium is the condition when there are CD4+ T cells in human body (X [not equal to] 0), there is infected CD4+ T cells (Y [not equal to] 0), and CTL precursor is activated (W [not equal to] 0) so that can give the immune response. The control immune equilibrium is given by [E.sub.2]([X.sup.*], [Y.sup.*], [W.sup.*], [Z.sup.*], [H.sup.*]), where

[mathematical expression not reproducible]. (13)

The equilibriums will exist if the conditions A = [([LAMBDA]c-[theta]q).sup.2]-2[beta][b.sub.1] ([LAMBDA]c + [theta]q) + [([beta][b.sub.1]).sup.2] [greater than or equal to] 0, [theta] > [beta][Y.sup.*] and [R.sub.01] + [beta][Y.sup.*]/[theta] > 1 are satisfied.

The Jacobi matrix of the model in (2) that evaluated in the equilibrium point [E.sub.2] is

[mathematical expression not reproducible]. (14)

By using the control immune equilibrium condition, the Jacobi matrix J([E.sub.2]) could be simplified into

[mathematical expression not reproducible]. (15)

From the Jacobi matrix J([E.sub.2]) above, we get the polynomial characteristic equations as follows:

([lambda] + [gamma]) ([lambda] + [b.sub.2]) ([[lambda].sup.3] + [a.sub.1][[lambda].sup.2] + [a.sub.2][lambda] + [a.sub.3]) = 0 (16)

where [a.sub.1] = [LAMBDA]/[X.sup.*], [a.sub.2] = [[beta].sup.2] [X.sup.*][Y.sup.*] + p[Y.sup.*][W.sup.*](c[X.sup.*]-q), and [a.sub.3] = p[Y.sup.*][W.sup.*]([LAMBDA]c-[LAMBDA]q/[X.sup.*]-c[beta][X.sup.*][Y.sup.*]).

Based on the characteristic equations (16) above, we obtain eigen values [[lambda].sub.1] = -[gamma], [[lambda].sub.2] = -[b.sub.2], and [[lambda].sub.3,4,5] is the roots of the following equation:

[[lambda].sup.3] + [a.sub.1][[lambda].sup.2] + [a.sub.2][lambda] + [a.sub.3] = 0. (17)

The argument of the first and second eigen values is n, [absolute value of (arg([[lambda].sub.1]))] = [absolute value of (arg([[lambda].sub.2]))] = [pi]. Next, we will determine the argument of the third, fourth, and fifth eigen values. The third, fourth, and fifth eigen values are roots of polynomial with degree of three (17), satisfying the following conditions:

(i) [[lambda].sub.3] + [[lambda].sub.4] + [[lambda].sub.5] = -[a.sub.1]

(ii) [[lambda].sub.3] [[lambda].sub.4] + [[lambda].sub.4] [[lambda].sub.5] + [[lambda].sub.3] [[lambda].sub.5] = [a.sub.2]

(iii) [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] = -[a.sub.3]

Then, we will analyze two conditions when [a.sub.3] < 0 and [a.sub.3] > 0. If [a.sub.3] < 0, then we have Theorem 5. On the other hand, if [a.sub.3] > 0 then we obtain Theorem 6 below.

Theorem 5. If [LAMBDA](c[X.sup.*] -q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] < 1, then the immune control equilibrium [E.sub.2] of model (2) is unstable.

Proof. Let [a.sub.3] = p[Y.sup.*][W.sup.*]([LAMBDA]c-[LAMBDA]q/[X.sup.*]-c[beta][X.sup.*][Y.sup.*]) < 0; we obtained [LAMBDA](c[X.sup.*]-q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] < 1. Since [a.sub.3] < 0, then based on the third condition (iii) above we have [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] > 0. There are two possibilities, either all the eigen values of the polynomial (17) are real number or one of the eigen values is real number and two other eigen values are complex numbers.

(1) Let all the eigen values [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] be real number. Because of [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] > 0, then there are two possibili ties. First, all the eigen values are real positive number. If all the eigen values are real positive numbers, then it is clear that [absolute value of (arg([[lambda].sub.3]))] = [absolute value of (arg([[lambda].sub.4]))] = [absolute value of (arg([[lambda].sub.5]))] = 0 < [alpha][pi]/2. Second, one of the eigen values is real positive number and two others are real negative number. Let [[lambda].sub.3] [member of] R, then [absolute value of (arg([[lambda].sub.3]))] = 0 < [alpha][pi]/2.

(2) Let one of the eigen values be real number, [[lambda].sub.3] [member of] R, and two others eigen values are complex numbers, [[lambda].sub.4] = a + bi, [[lambda].sub.5] = a-bi, where a,b [member of] R. Because [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] > 0 and [[lambda].sub.4] [[lambda].sub.5] > 0, then we obtained [[lambda].sub.3] > 0. As a result, [absolute value of (arg([[lambda].sub.3]))] = 0 < [alpha][pi]/2.

From the two conditions above, we find that the immune control equilibrium E2 of the model (2) is unstable whenever [LAMBDA](c[X.sup.*] -q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] < 1.

Theorem 6. The immune control equilibrium [E.sub.2] of the model (2) is asymptotically stable if [LAMBDA](c[X.sup.*]-q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] > 1 for some fractional order [alpha].

Proof. Let we consider the polynomial

[[lambda].sup.3] + [a.sub.1][[lambda].sup.2] + [a.sub.2][lambda] + [a.sub.3] = 0 (18)

where [a.sub.1] = [LAMBDA][X.sup.*], [a.sub.2] = [[beta].sup.2][X.sup.*][Y.sup.*] + p[Y.sup.*][W.sup.*] (c[X.sup.*]-q), and [a.sub.3] = p[Y.sup.*][W.sup.*]([LAMBDA]c-[LAMBDA]q/[X.sup.*]-c[beta][X.sup.*][Y.sup.*]).

Let [a.sub.3] = p[Y.sup.*][W.sup.*]([LAMBDA]c-[LAMBDA]q/[X.sup.*]-c[beta][X.sup.*][Y.sup.*]) > 0. Hence we obtained [LAMBDA](c[X.sup.*]-q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] > 1. If [a.sub.3] > 0, then based on the third condition (iii) above we have [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] < 0.

(i) Let the polynomial characteristic (18) above have complex number roots. Let one of the eigen values be real number, [[lambda].sub.3] [member of] R, and two other eigen values are complex number, [[lambda].sub.4] = a + bi, [[lambda].sub.5] = a-bi, where a,b [member of] R. Since [[lambda].sub.3] [[lambda].sub.4] [[lambda].sub.5] < 0 and [[lambda].sub.4] [[lambda].sub.5] > 0, then we obtained [[lambda].sub.3] < 0. As a result, [absolute value of (arg([[lambda].sub.3]))] = [pi] > [alpha][pi]/2. In addition, we can found a fractional order value [alpha] [member of] (0,1] such that the arguments of the third and fourth eigen values satisfy [absolute value of (arg([[lambda].sub.4]))] > [alpha][pi]/2 and [absolute value of (arg([[lambda].sub.5]))] > [alpha][pi]/2. As a result, the immune control equilibrium [E.sub.2] of model (2) is asymptotically stable for some fractional order [alpha].

(ii) Let all the eigen values of the polynomial characteristic (18) above [[lambda].sub.3], [[lambda].sub.4], [[lambda].sub.5] be real numbers. Then we will determine the conditions such that all roots of the polynomial (18) above are either negative or complex roots with negative real parts by using Routh-Hurwitz criterion.

Based on Routh-Hurwitz criteria, the polynomial (18) will have real negative or real negative parts roots if it satisfies [a.sub.1], [a.sub.2], [a.sub.3] > 0 and [a.sub.1][a.sub.2]-[a.sub.3] > 0. Hence, we obtained the following conditions:

(a) Since all parameters are assumed to be positive and [X.sup.*] [not equal to] 0, then it is clear that [a.sub.1] > 0.

(b) The coefficient [a.sub.2] will have positive value if [a.sub.3] is positive.

(c) The coefficient [a.sub.3] will have positive value if it satisfies [LAMBDA](c[X.sup.*]-q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] > 1.

(d) It is clear that the coefficients [a.sub.1],[a.sub.2],[a.sub.3] satisfy [a.sub.1] [a.sub.2]-[a.sub.31] > 0.

Based on [20], for [alpha] [member of] (0,1] these conditions are sufficient but not necessary. Therefore, the argument of [[lambda].sub.3], [[lambda].sub.4], [[lambda].sub.5] will satisfy [absolute value of (arg([[lambda].sub.4]))] = [absolute value of (arg([[lambda].sub.4]))] = [absolute value of (arg([[lambda].sub.5]))] = [pi] > [pi]/2 if the condition [LAMBDA](c[X.sup.*]-q)/c[beta][([X.sup.*]).sup.2][Y.sup.*] > 1 is satisfied. As a result, the immune control equilibrium E2 of the model (2) is asymptotically stable. This completes the proof.

4. Sensitivity Analysis

In this section, we will analyze the sensitivity of parameters of the model (2) above. Sensitivity analysis is used to determine the relative importance of model parameters to disease transmission and prevalence based on the sensitivity index of basic reproduction number [R.sub.0]. The calculation of sensitivity index of [R.sub.0] follows the approach in Chitnis [21].

Definition 7 (see [21]). The normalized forward sensitivity index of a variable, m, that depends differentially on a parameter, k, is defined as

[Y.sup.m.sub.k] [??] [partial derivative]m/[partial derivative]k x k/m. (19)

Based on Definition 7, the sensitivity indices of [R.sub.01] with respect to each parameter such as [beta], [LAMBDA], a, and [theta] can be computed in the same way as (19). The sensitivity indices of [R.sub.01] are given in Table 2. Based on Table 2, it can be seen that parameters [beta] and [LAMBDA] positively affect the rate of model changed. Respectively, the parameters a and [theta] negatively affect the rate of model changed. But, we did not know what is the relative importance of model parameters of [R.sub.01]. So, we also analyze the sensitivity parameters of [R.sub.02].

Based on Definition 7, the sensitivity indices of [R.sub.02] with respect to each parameter such as [beta], [LAMBDA], a, [theta], [b.sub.1], q, and c can be computed in the same way as (19). The sensitivity indices of parameter [LAMBDA] are ([partial derivative][R.sub.02]/[partial derivative][LAMBDA])([LAMBDA]/[R.sub.02]) = [LAMBDA](ca-q[beta])/[beta](-a[b.sub.1]-[LAMBDA]q) + ac[LAMBDA]; we next substituted the value of parameters of Table 4 so we obtain ([partial derivative][R.sub.02]/[partial derivative][LAMBDA])([LAMBDA]/[R.sub.02]) = 1. The sensitivity indices of [R.sub.02] are given in Table 3.

Based on Table 3, it can be seen that if production rate of uninfected CD4+ T cells (A) is increased (decreased) about 10%, then [R.sub.02] will increase (decrease) about 10%. Respectively, if death rate of infected CD4+ T cells (a) is increased (decreased) about 10%, then [R.sub.02] will decrease (increase) about 9.94%, and the same for other parameters.

5. Numerical Simulation

In this section, we present numerical simulations to show the dynamics of the model at free virus condition, CTL-exhaustion condition, and control immune condition using MATLAB R2009a. The initial conditions of simulations of free virus condition are [X.sub.0] = 3.6, [Y.sub.0] = 0.2, [W.sub.0] = 0.01, [Z.sub.0] = 0.01, and [H.sub.0] = 0.001. The parameter values are presented in Table 4. Based on Table 4, we have that basic reproduction number [R.sub.0] is [R.sub.01] = [beta][LAMBDA]/a[theta] = (0.15)(0.2)/(0.8)(0.05) = 0.75 < 1 that shows free virus condition or there is no spread of viral infection. This simulation used variations of fractional order a and time interval for 2500 time units.

The result of simulation of free virus condition can be seen in Figure 2. Based on Figure 2, greater the fractional order a leads to the increase of concentrations of uninfected CD4+ T cells but concentrations of infected CD4+ T cells are decreasing that show free virus conditions. Besides, concentrations of precursor CTL cells are decreasing caused bythedecreaseofCD4+Tcells. That also implies the decrease of effector CTL. Meanwhile, concentrations of CD8+ T cells are decreasing because they were differentiated into precursor CTL cells.

Next, we will interpret the simulation of CTL-exhaustion condition. The initial conditions of simulations of CTL-exhaustion condition are [X.sub.0] = 1.5, [Y.sub.0] = 0.8, [W.sub.0] = 0.6, [Z.sub.0] = 0.3, and [H.sub.0] = 0.001. The parameter values are presented in Table 4 except for [beta], a, and [theta], which are [beta] = 0.216, a = 0.4, and [theta] = 0.02. Therefore, we have basic reproduction number [R.sub.01] = (0.216)(0.2)/(0.4)(0.02) = 5.4 > 1 and [R.sub.02] = [beta][LAMBDA]/a[theta] + [b.sub.1][[beta].sup.2]/[theta](q[beta]-ca) = 10.1 > 1 which is condition when there is spread of viral infections in the human body. This simulation used variations of fractional order a and time interval for 2500 time units.

Based on Figure 3, it can be seen that greater the fractional order [alpha] leads to the decrease of concentrations of uninfected CD4+ T cells but concentrations of infected CD4+ T cells are increasing. That implies the conditions of viral infections. But, concentrations of effector CTL cells are decreasing. That implies the CTL-exhaustion conditions where CTL are not activated yet so they cannot against the viral infections. Meanwhile, concentrations of CD8+ T cells were decreased because they differentiated into precursor CTL cells.

Last, we will show the simulation of control immune conditions. The initial conditions of simulations of CTL-exhaustion condition are [X.sub.0] = 5,[Y.sub.0] = 1,[W.sub.0] = 2.5,[Z.sub.0] = 0.8, and [H.sub.0] = 0.001. The parameter values are showed in Table 4 except for [theta], that is, [theta] = 0.02. Therefore, we have basic reproduction number [R.sub.01] = 1.87 > 1 and [R.sub.02] = 1.01 > 1 which is condition when there is spread of viral infections in the human body This simulation used variations of fractional order a and time interval for 2500 time units. Based on Figure 4, it can be seen that greater the fractional order a leads to the increase of concentrations of uninfected CD4+ T cells but concentrations of infected CD4+ T cells are decreasing. Meanwhile, concentrations of precursor and effector CTL cells are increasing, which means the precursor and effector CTL cells are activated so they can fight against the viral infections. That implies control immune conditions. Meanwhile, concentrations of CD8+ T cells were decreased because they differentiated into precursor CTL cells.

6. Conclusions

We were discussed about the stability analysis of fractional order mathematical model of HSCs gene-based therapy. The model has three equilibriums points. The existence and stability of the equilibriums were achieved. Based on the numerical simulations, we conclude that HSCs gene-based therapy has potential to reduce the viral infections because it can decrease the concentration of infection cells and increase the concentration of the immune cells.

https://doi.org/10.1155/2018/6180892

Data Availability

No data were used to support this study.

Conflicts of Interest

The authors state that there are no conflicts of interest concerning the preparation and publication of the manuscript.

References

[1] R. Nairn and M. Helbert, Immunology for Medical Students, Mosby Elsevier, Philadelphia, 2nd edition, 2007.

[2] D. Wodarz, Killer cell dynamics, vol. 32 of Interdisciplinary Applied Mathematics, Springer-Verlag, New York, 2007.

[3] S. G. Kitchen, B. R. Levin, G. Bristol et al., "In vivo suppression of HIV by antigen specific T cells derived from engineered hematopoietic stem cells," PLoS Pathogens, vol. 8, no. 4, Article ID e1002649, 2012.

[4] A. Korpusik, "Hematopoietic stem cell based therapy of immunosuppressive viral infection--Numerical simulations," Biocybernetics and Biomedical Engineering, vol. 34, no. 2, pp. 125-131, 2014.

[5] A. Korpusik and M. Kolev, "Single injection of CD8+ T lymphocytes derived from hematopoietic stem cells--Mathematical and numerical insights," BioSystems, vol. 144, pp. 46-54, 2016.

[6] S. Das and P. K. Gupta, "A mathematical model on fractional Lotka-Volterra equations," Journal of Theoretical Biology, vol. 277, pp. 1-6, 2011.

[7] E. Ahmed and A. S. Elgazzar, "On fractional order differential equations model for nonlocal epidemics," Physica A: Statistical Mechanics and its Applications, vol. 379, no. 2, pp. 607-614, 2007

[8] F. Fatmawati, E. M. Shaiful, and M. I. Utoyo, "A Fractional-Order Model for HIV Dynamics in a Two-Sex Population," International Journal of Mathematics and Mathematical Sciences, vol. 2018, Article ID 6801475, 11 pages, 2018.

[9] J. Huo and H. Zhao, "Dynamical analysis of a fractional SIR model with birth and death on heterogeneous complex networks," Physica A: Statistical Mechanics and its Applications, vol. 448, pp. 41-56, 2016.

[10] J. Huo, H. Zhao, and L. Zhu, "The effect of vaccines on backward bifurcation in a fractional order HIV model," Nonlinear Analysis: Real World Applications, vol. 26, pp. 289-305, 2015.

[11] C. M. A. Pinto and A. R. M. Carvalho, "The HIV/TB coinfection severity in the presence of TB multi-drug resistant strains," Ecological Complexity, vol. 32, pp. 1-20, 2017

[12] M. Saeedian, M. Khalighi, N. Azimi-Tafreshi, G. R. Jafari, and M. Ausloos, "Memory effects on epidemic evolution: the susceptible-infected-recovered epidemic model," Physical Review E: Statistical, Nonlinear, and Soft Matter Physics, vol. 95, no. 2, 022409, 9 pages, 2017

[13] T. Sardar, S. Rana, and J. Chattopadhyay, "A mathematical model of dengue transmission with memory," Communications in Nonlinear Science and Numerical Simulation, vol. 22, no. 1-3, pp. 511-525, 2015.

[14] G. Huang, Y. Takeuchi, and A. Korobeinikov, "HIV evolution and progression of the infection to AIDS," Journal of Theoretical Biology, vol. 307, pp. 149-159, 2012.

[15] S. Z. Rida and A. A. Arafa, "New method for solving linear fractional differential equations," International Journal of Differential Equations, Art. ID 814132, 8 pages, 2011.

[16] O. Diekmann, J. A. P. Heesterbeek, and M. G. Roberts, "The construction of next-generation matrices for compartmental epidemic models," Journal of the Royal Society Interface, vol. 7, no. 47, pp. 873-885, 2010.

[17] F. Brauer and C. Castillo-Chavez, Mathematical Models in Population Biology and Epidemiology, vol. 40 of Texts in applied mathematics, Springer, New York, NY, USA, 2001.

[18] M. S. Tavazoei and M. Haeri, "Chaotic attractors in incommensurate fractional order systems," Physica D: Nonlinear Phenomena, vol. 237, no. 20, pp. 2628-2637, 2008.

[19] P. van den Driessche and J. Watmough, "Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission," Mathematical Biosciences, vol. 180, pp. 29-48, 2002.

[20] E. Ahmed, A. M. A. El-Sayed, and H. A. A. El-Saka, "On some Routh-Hurwitz conditions for fractional order differential equations and their applications in Lorenz, Rossler, Chua and Chen systems," Physics Letters A, vol. 358, no. 1, pp. 1-4, 2006.

[21] N. Chitnis, J. M. Hyman, and J. M. Cushing, "Determining important parameters in the spread of malaria through the sensitivity analysis of a mathematical model," Bulletin of Mathematical Biology, vol. 70, no. 5, pp. 1272-1296, 2008.

Mohammad Imam Utoyo (iD), Windarto (iD), and Aminatus Sa'adah

Department of Mathematics, Faculty of Science and Technology, Universitas Airlangga, Indonesia

Correspondence should be addressed to Mohammad Imam Utoyo; m.i.utoyo@fst.unair.ac.id

Received 18 May 2018; Accepted 7 July 2018; Published 1 August 2018

Academic Editor: Hans Engler

Caption: Figure 1: HSCs gene-based therapy transmission diagram.

Caption: Figure 2: The dynamics of mathematical model of HSCs gene-based therapy for free virus conditions.

Caption: Figure 3: The dynamics of mathematical model of HSCs gene-based therapy for CTL-exhaustion conditions.

Caption: Figure 4: The dynamics of fractional order mathematical model of HSCs gene-based therapy for control immune conditions.

Table 1: The description of parameters in the model. Parameter Description Unit [LAMBDA] Production rate of Concentration cells/time uninfected CD4+ T cells unit [theta] Death rate of uninfected 1/ time unit CD4+ T cells [beta] Infection rate Concentration [cells.sup.-1]/time unit a Death rate of infected 1/ time unit CD4+ T cells p Effector CTL rate of Concentration killing infected CD4+ T [cells.sup.-1]/time unit cells m Fractions of HSCs that -- successfully differentiated into functional CTL cells c Proliferation rate of Concentration precursor CTL [cells.sup.-1]/time unit q Fraction of precursor CTL -- cells that differentiated into effector CTL cells [b.sub.1] Death rate of precursor 1/ time unit CTL [b.sub.2] Death rate of effector 1/ time unit CTL [gamma] Differentiation rate of 1/ time unit CD8+ T lymphocytes Table 2: The sensitivity indices of [R.sub.01]. Parameter Sensitivity indices [beta] 1 [LAMBDA] 1 a -1 [theta] -1 Table 3: The sensitivity indices of [R.sub.02]. Parameter Sensitivity indices [beta] 0.945 [LAMBDA] 1 a -0.994 [theta] -1 [b.sub.1] -0.049 q -0.006 c 0.055 Table 4: Parameters value of the model. Parameter Description Value Source [LAMBDA] Production rate of 0.2 [14] uninfected CD4+ T cells [theta] Death rate of uninfected 0.05 Assumption CD4+ T cells [beta] Infection rate 0.15 [14] a Death rate of infected 0.8 [14] CD4+ T cells p Effector CTL rate of 0.016 [14] killing infected CD4+ T cells m Fractions of HSCs that 0.01 Assumption successfully differentiated into functional CTL cells c Proliferation rate of 0.1 Assumption precursor CTL q Fraction of precursor CTL 0.1 Assumption cells that differentiated into effector CTL cells [b.sub.1] Death rate of precursor 0.05 [14] CTL b.sub.2] Death rate of effector 0.05 [14] CTL [gamma] Differentiation rate of 0.005 Assumption CD8+ T cells

Printer friendly Cite/link Email Feedback | |

Title Annotation: | Research Article |
---|---|

Author: | Utoyo, Mohammad Imam; Windarto; Saadah, Aminatus |

Publication: | International Journal of Mathematics and Mathematical Sciences |

Article Type: | Report |

Date: | Jan 1, 2018 |

Words: | 5823 |

Previous Article: | Annular Bounds for the Zeros of a Polynomial. |

Next Article: | New Extension of Beta Function and Its Applications. |

Topics: |