Parameter identification of hematopoiesis mathematical model – periodic chronic myelogenous leukemia

Introduction

Leukemia is a kind of white blood cell cancer. It may be categorized based on the characteristics of the disease (acute or chronic). A patient with acute leukemia usually dies within several months, while a patient with chronic leukemia may live for years.

The diagnosis of chronic myelogenous leukemia (CML) is usually based on detection of the Philadelphia chromosome. This abnormality arises as the result of a problem in cell division within the bone marrow. It is a type of translocation, meaning that genetic material from two different chromosomes switches places. In this particular abnormality, chromosomes 9 and 22 are involved. CML arises in a bone marrow stem cell that is the precursor to all the types of blood cells (white blood cells, red blood cells and platelets). This hypothesis is derived from two pieces of evidence. First, in CML the Philadelphia chromosome is found in all of the hematopoietic lineages which come from the hematopoietic stem cell. Second, periodic oscillations are observed in three types of blood cells in one rare sort of CML, which leads to the mentioned hypothesis.

To provide a quantitative basis for prediction of behavior of stem cells and blood cells under hematopoietic diseases, several models have been developed. To our knowledge, Mackey, in 1978 [1] proposed the first model of hematopoietic stem cell dynamics, which has been developed by many authors, including Haurie and Mackey [2], who focused on some hematological diseases which originate from stem cells and specified differentiation and maturation processes in hematopoiesis, Bernard and Mackey [3], who concentrated on the production and oscillations of white blood cells, observed in patients with cyclic neutropenia, and Pujo-Menjouet and Mackey [4], who focused on applications of the model proposed in [3] to chronic myelogenous leukemia. Coljin and Mackey [5, 6] extracted a global model of the hematopoietic system that concentrates on either stem cell dynamics or the leukocyte, erythrocyte and platelet cell lines. Foley and Mackey [7] proposed a delay differential equation model for the regulation of neutrophil production under the effects of G-CSF. Models in the form of partial differential equations have been proposed by Scholz, Engel and Loeffler [8] that simulate the damage caused by cytotoxic drugs on the stem cells and leukocytes. The hematopoietic system has also been formulated as partial differential equations (PDE) by Ostby et al. [9] and ordinary differential equations (ODE) by Panetta et al. [10], Friberg et al. [11], Scholz et al. [8] and Shochat et al. [12].

Recently, a description of dynamics of chronic and acute myeloid leukemias was proposed by Adimy et al. [13, 14]. In [13], a compartmental model with discrete delays was suggested. The conditions for local stability of steady states and for Hopf bifurcation due to delay parameters were also given. In [14], a model with distributed delays was considered.

The goal of this paper is to identify some unknown parameters of the model and improve the discrepancy between the model and clinical data.

The model of the hematopoietic system

Differential equations (partial, ordinary or delay), stochastic processes, Boolean networks, etc. are some mathematical tools used in hematological modeling. The desired level of description of the model often determines the mathematical tools. Since delays arise in biological systems because of properties inherent to the different processes (time lag due to maturation, transmission of an impulse, etc.), delay differential equations are of great importance in describing biological systems such as the hematopoietic system.

Previously, leukocyte [15], erythrocyte [16] and platelet [17] models have been proposed separately. Mackey in [18] connected all these three distinct models. The model describes four distinct compartments: the hematopoietic stem cells (HSC) and the circulating leukocytes, platelets and erythrocytes. The stem cells are pluripotential as well as self-renewing, and are differentiable into the leukocytes, erythrocytes or platelets. Fig. 1 shows the full model of the system represented by Colijn and Macky [18]. The model considers the hematopoietic system as a set of compartments: the pluripotential, non-proliferating stem cells, leukocytes, erythrocytes and platelets are denoted by Q, N, R and P in unit of cells/kg, respectively. KN(N), KR(R) and KP(P) represent the differentiation rate of HSC into these three lines in units of days–1, respectively. These differentiation rates are a monotone decreasing function of the numbers of circulating cells of the relevant type (negative feedback). More discussion of the forms of KN(N), KR(R) and KP(P) will follow [18].

A portion of non-proliferating stem cells differentiates into leukocyte, erythrocyte and platelet cell lines at the rate of KN(N), KR(R) and KP(P), respectively. The remaining portion proliferates at the rate of  (Q) (in units of days–1). In the proliferation compartment the cells divide taking a time of s (units of days). So, in the equation describing Q, the number of cells passing the proliferation phase is twice the number of those entering this phase a time s ago. But the cells traversing the proliferation compartment are not quite twice those entering, because some random losses occur at a rate of s (in units of days–1) over the proliferating time. So, the number of non-proliferating cells entering per unit time from the proliferation compartment is 2e–ss (Q(t – s). In this paper the notation X is used for X (t – ) hereinafter. Thus, the delay differential equation describing Q is (1).

In each leukocyte, erythrocyte and platelet cell line there is a dimensionless amplification parameter AN, AR and AP, respectively, which represents several stages of cell divisions.

We use the equations related to red blood cells and platelets without any changes from [18], but for the leukocyte line we take some ideas from the model represented by Foley and Mackey [19], which is mainly about modeling the effects of G-CSF on leukocytes and stem cells. We consider the proliferative and non-proliferative population of precursor cells of leukocytes separately. This separation leads to some parameters being added to the model. We assume that all these new parameters are fixed and find their best values so that the model has a better fit of the quantitative data on leukocytes and platelets. Fig. 2 provides an overview of the modified model structure.

In the leukocyte line, the upper box and lower box are the compartment of proliferating precursor cells and non-proliferating precursor cells, respectively. The leukocytes which enter the circulating leukocyte compartment after amplification and maturation are denoted by N. They are randomly lost at a rate N due to apoptosis. NN is the rate of random loss of non-proliferating precursor cells while NP represents the rate of random loss of proliferating precursor cells. We assumed that the parameters N, AN, KN, NN and NP are independent of G-CSF. NN is non-proliferative neutrophil precursors maturation time and NP is proliferative neutrophil precursors maturation time. The equation governing the leukocytes is then

dN

dt = – NN + ANKN (NNP) QNP e–NP NP –NN –n [2]

with –n satisfying NN = Vn–n. Vn is the velocity that cells in the non-proliferating precursor compartment age with. An increase in Vn leads to decreased transit time NN.

Simulation and fitting

We extracted clinical data for platelets and leukocytes from [20] and tried to find the parameters of the leukocyte line such that the model fits real data more precisely. We have used B. Ermentrout’s implementation of the Marquar-

dt-Levenberg method, described by Press et al. [21]. Our

N

method to do so is minimizing the function x2 = (ys – yi)2

i = 1

where yi are the observed data points and ys are the simu-

lated data points. We use MATLAB software for all of our calculations. In this paper, hereinafter LS represents the least

N

square value (ys – yi)2. Considering proliferating and non-pro-

i = 1

liferating precursor compartments in the leukocyte line and

estimating the parameters nn, nn, nP and Vn displayed

in Table 1 improves the LS value for both platelets and

leukocytes data in comparison with Mackey's (2005) work (Table 2 and Table 3).

Results

• Although the parameters nn, nP, nn and nP are just involved in leukocyte dynamics, a change in their value affects the count of both leukocytes and platelets.

• High values of nP (proliferative neutrophil precursors maturation time) cause a decrease in the mean value of the leukocytes.

• Variation in the parameter Vn has the least effect on the model in comparison with the other parameters.

• The dynamics of the system are more sensitive to variations in nP in comparison with nn. Looking at equ. (2) more precisely supports this fact. As can be observed the parameter KN depends on nP in the mentioned equation.

The lower the value of the parameter nP, the higher the mean value of the leukocytes is.

• The higher the value of the parameter 2 is, the more the amount of the peak and the mean value of N and P are.

• The parameters nP and 2 affect the shape of the model at its first and last oscillates.

In conclusions, this article identifies some parameters of the hematopoietic system so that the model has a better fit of the quantitative data on leukocytes and platelets in comparison with fitness done by the others. The effects of mentioned parameters on behavior of leukocytes and platelets are also investigated.

Discussion

The major intention of this paper is to construct a physiological model which can map major identified mechanisms of leukopoiesis to provide a deeper insight into this complex biological process. A new model of leukopoiesis has been proposed. In fact, precursor cells of the leukocyte line have been considered as two separate groups: proliferating precursor cells and non-proliferating precursor cells. As a result, more parameters have appeared in the model and identifying the new parameters has resulted in a better fit of clinical data and the data extracted from the model for both platelets and leukocytes. Using the new model can lead to better prediction of the behavior of leukocytes and platelets for PCML disease.

References

 1. Mackey MC. A unified hypothesis for the origin of aplastic anemia and periodic haematopoiesis. Blood 1978; 51: 941-56.

 2. Haurie C, Mackey MC, Dale DC. Cyclical neutropenia and other periodic hematological diseases: a review of mechanisms and mathematical models. Blood 1998; 92: 2629-40.

 3. Bernard S, Belair J, Mackey M. Oscillations in cyclical neutropenia: new evidence based on mathematical modeling. J Theor Biol 2003; 223.

 4. Pujo-Menjouet L, Bernard S, Mackey MC. Long period oscillations in a G0 model of hematopoietic stem cells. SIAM J Appl Dynam Syst 2001; 4: 312-32.

 5. Colijn C, Mackey MC. A mathematical model of hematopoiesis. I. Periodic chronic myelogenous leukemia. J Theoret Biol 2006; 237: 117-132.

 6. Colijn C, Mackey MC. A mathematical model of hematopoiesis. II. Cyclical neutropenia. J Theoret Biol 2006; 237: 133-46.

 7. Foley C, Mackey MC. Dynamic hematological disease: a review. J Math Biol 2008; DOI 10.1007/s00285-008-0165-3.

 8. Scholz M, Engel C, Loeffler M. Modelling human granulopoiesis under polychemotherapy with G-CSF support. J Math Biol 2005; 50: 397-439.

 9. Ostby I, Rusten LS, Kvalheim G, Grottum P. A mathematical model for reconstitution of granulopoiesis after high dose chemotherapy with autologous stem cell transplantation. J Math Biol 2003; 47: 101-36.

10. Panetta JC, Kirstein MN, Gajjar A, Nair G, Fouladi M, Stewart CF. A mechanistic mathematical model of temozolomide myelosuppression in children with high-grade gliomas. Math Biosc 2003; 186: 29-41.

11. Friberg LE, Henningsson A, Maas H, Nguyen L, Karlsson MO. Model of chemotherapy induced myelosuppression with parameter consistency across drugs. J Clin Oncol 2002; 20: 4713-21.

12. Shochat E, Rom-Kedar V, Segel L. G-CSF control of neutrophil dynamics in the blood. Bull Math Biol 2007; 69: 2299-338.

13. Adimy M, Crauste F, El Abdllaoui A. Asymptotic behavior of a discrete maturity structured system of hematopoietic stem cells dynamics with several delays. Math Model Nat Phenom 2006; 1: 1-22.

14. Adimy M, Crauste F, El Abdllaoui A. Discrete maturity structured model of cell differentiation with applications to acute myelogenous leukemia. J Biol Syst 2008; 16: 395-424.

15. Hearn T, Haurie C, Mackey M. Cyclical neutropenia and the peripherial control of white blood cell production. J Theor Biol 1998; 192: 167-81.

16. Mackey MC. Dynamic haematological disorders of stem cell origin. In: Vassileva-Popova JG, Jensen EV (eds.). Biophysical and Biochemical Information Transfer in Recognition. Plenum Publishing Corp, New York 1979; 373-409.

17. Be´lair J, Mackey M, Mahaffy J. Age-structured and twodelay models for erythropoiesis. Math Biosci 1995; 128: 317-46.

18. Colijn C, Mackey M. A mathematical model of hematopoiesis: I. Periodic chronic myelogenous leukemia. J Theor Biol 2005; 237: 117-32.

19. Foley C, Mackey MC, 2008. Dynamic hematological disease: a review. J Math Biol 2009; 58: 285-322.

20. Iizuka Y, Horikoshi A, Sekiya S, Sawada U, Ohshima T, Amaki I. Periodic fluctuation of leukocytes, platelets and reticulocytes in a case of chronic myelocytic leukemia: The relation between leukocyte counts, CFU-C colony formation, CSA and CIA. Nihon Ketsueki Gakkai Zasshi 1984; 47: 71-9.

21. Press W, Teukolsky S, Vetterling W, Flannery B. Numerical Recipes in C, second ed. Cambridge University, Cambridge 1992.



Address for correspondence



Anahita Safarishahrbijari

K.N. Toosi University of Technology

#19 Pardiss Alley, Sadra Street, Vanak Square

Tehran, Iran, Zip: 19991 43344

PO Box: 19395-1999

Tel: + 98 21 84063235

Fax: + 98 21 88674748

e-mail:ghaffari@kntu.ac.ir, anahita.safari@gmail.com



Submitted: 7.05.2012

Accepted: 23.05.2012