O Angulo

O Gandrillon

F Crauste
email: crauste@math.univ-lyon1.fr

Investigating the role of the experimental protocol in phenylhydrazine-induced anemia on mice recovery

Keywords: erythropoiesis model, experimental protocol modeling, nonlinear age-structured system

Production of red blood cells involves growth-factor mediated regulation of erythroid progenitor apoptosis and self-renewal. During severe anemia, characterized by a strong fall of the hematocrit followed by a recovery phase, these controls allow a fast recovery of the hematocrit and survival of the organism. Using a mathematical model of stress erythropoiesis and an ad hoc numerical method, we investigate the respective roles of anemia-inducing phenylhydrazine injections and physiological regulation on the organism's recovery. By explicitly modeling the experimental protocol, we show that it mostly characterizes the fall of the hematocrit following the anemia and its severeness, while physiological process regulation mainly controls the recovery. We confront our model and our conclusions to similar experiments inducing anemia and show the model's ability to reproduce several protocols of phenylhydrazine-induced anemia. In particular, we establish a link between phenylhydrazine effect and the severeness of the anemia.

Introduction

Erythropoiesis, the process of production of red blood cells, is performed through complex regulatory processes, as part of the more general process of blood cell production. Red blood cells are produced in the bone marrow, where hematopoietic stem cells, that have abilities of self-renewal and differentiation in all blood cell types, generate immature erythroid cells called progenitors. Throughout successive differentiating divisions, erythroid progenitors acquire maturity to ultimately become mature red blood cells, called erythrocytes, that enter the bloodstream in order to carry oxygen to, and CO 2 from organs and tissues.

The continuous production of erythroid cells is permanently controlled in order to adapt very quickly to changes in, or needs of the organism. One of the main feedback controls, discovered in the early 1990's by [START_REF] Koury | Erythropoietin retards DNA breakdown and prevents programmed death in erythroid progenitor cells[END_REF], deals with erythroid progenitor death by apoptosis, a programmed cell death. [START_REF] Koury | Erythropoietin retards DNA breakdown and prevents programmed death in erythroid progenitor cells[END_REF] showed that, during an anemia (lack of red blood cells), a growth factor named erythropoietin (Epo) was released by the kidneys and inhibited progenitor apoptosis, allowing a fast production of numerous erythrocytes to recover from the anemia. This control of red blood cells upon their production is crucial for erythropoiesis regulation.

Apart from dying by apoptosis, erythroid progenitors experience differentiation in more mature cells or self-renew. Self-renewal [START_REF] Watt | Out of Eden: stem cells and their niches[END_REF] is the ability of a cell to divide and give two daughter cells that have the same maturity than the mother cell, while keeping at all time the ability to engage in a differentiation process (i.e. to give two daughter cells, one of which at least being more mature than the mother cell). Self-renewal of erythroid progenitors has been shown during stress erythropoiesis (an anemia for instance) by Bauer er al. (1999), [START_REF] Gandrillon | TGF-beta cooperates with TGF-alpha to induce the self-renewal of normal erythrocytic progenitors: evidence for an autocrine mechanism[END_REF], and [START_REF] Pain | EGF-R as a hemopoietic growth factor receptor: The c-erbB product is present in normal chicken erythrocytic progenitor cells and controls their self-renewal[END_REF]. It is controlled by glucocorticoids, growth factors whose production negatively depends on the number of erythrocytes.

Several mathematical models of erythropoiesis have been proposed over the last 30 years, in order to address the regulatory mechanisms and their role in stress or pathological conditions. [START_REF] Bélair | Age-structured and two-delay models for erythropoiesis[END_REF] proposed a model of erythropoiesis considering the influence of growth factors on stem cell differentiation in erythroid progenitors. This model was then improved by [START_REF] Mahaffy | Hematopoietic model with moving boundary condition and state dependent delay: applications in erythropoiesis[END_REF], and later analyzed in detail by Ackleh et al. (2002Ackleh et al. ( , 2006Ackleh et al. ( , 2008) ) and [START_REF] Banks | Modelling and Optimal Regulation of Erythropoiesis Subject to Benzene Intoxication[END_REF]. Another erythropoiesis model, inspired by the same article, in which Epo is the only growth factor supposed to act during erythropoiesis, was introduced in Adimy et al. (2006).

An important contribution to mathematical modeling of erythropoiesis is also due to Loeffler and his collaborators [START_REF] Loeffler | A comprehensive mathematical model of stem cell proliferation which reproduces most of the published experimental results[END_REF][START_REF] Loeffler | A mathematical model of erythropoiesis in mice and rats. Part 1. Structure of the model[END_REF][START_REF] Roeder | Quantitative stem cell biology: computational studies in the hematopoietic system[END_REF][START_REF] Roeder | A novel dynamic model of hematopoietic stem cell organization based on the concept of within-tissue plasticity[END_REF]Wichmann and Loeffler, 1985;[START_REF] Loeffler | A mathematical model of erythropoiesis in mice and rats. Part 1. Structure of the model[END_REF]. Their models consider feedback controls from progenitors on the stem cell level and from mature cells on progenitors, and are fitted to various experiments that induce severe anemias.

Most of these works were nevertheless performed before the role of Epo was definitely identified and long before erythroid progenitor self-renewal was hypothesized.

We proposed a new model of erythropoiesis [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF], based on previous works by [START_REF] Mackey | Unified hypothesis of the origin of aplastic anaemia and periodic hematopoïesis[END_REF] and co-authors, and Loeffler's group [START_REF] Loeffler | A mathematical model of erythropoiesis in mice and rats. Part 1. Structure of the model[END_REF][START_REF] Roeder | A novel dynamic model of hematopoietic stem cell organization based on the concept of within-tissue plasticity[END_REF][START_REF] Roeder | Quantitative stem cell biology: computational studies in the hematopoietic system[END_REF], to describe stress erythropoiesis in mice. This model is based on the description of erythroid progenitor and erythrocyte dynamics using delay equations. Progenitor apoptosis and self-renewal are regulated by the number of erythrocytes, hence implicitly describing growth factor-mediated regulation. The model's outputs were compared to data on phenylhydrazine-induced anemia in mice and showed first that it was relevant and necessary to account for erythroid progenitor self-renewal in order to explain experimental data. Second, they showed that it was necessary to account for a modification of erythrocyte lifespan after phenylhydrazine injections in order to reproduce the data, even though this had not been assumed in the initial modeling. Phenylhydrazine is a substance that kills circulating red blood cells but also damages surviving red blood cell membranes and yields long term deficiency of the erythropoiesis system [START_REF] Walter | Abnormal membrane surface properties during maturation of rat reticulocytes elicited by bleeding as measured by partition in two-polymer aqueous phases[END_REF]. Notably, phenylhydrazine-induced anemia is described in this model as a perturbation of the initial state of the system and is hence not explicitly modeled. In particular, the biological mechanisms involved in the reaction to phenylhydrazine injections are not considered: it is a model of recovery from an anemia but not of anemia itself.

The way anemia is induced in mice influences the strength of the anemia and the dynamics of the recovery. The two main ways of inducing anemia are either to bleed mice in order to remove red blood cells, or to inject mice with phenylhydrazine. Hematocrit values following bleeding-induced anemia display smooth recovery dynamics, while following phenylhydrazine-induced anemia they show a fast recovery phase that overcomes the initial hematocrit value and then slowly goes back to its initial state (see Figure 1), even Figure 1: Phenylhydrazine versus bleeding-induced anemia. Hematocrits of mice submitted to severe anemia either by phenylhydrazine injections (red squares) or bleeding (blue triangles), mean fold change values are displayed (note that standard deviations are small and not visible using square and triangle symbols as well as normalized scales). The experimental protocol for phenylhydrazine-induced anemia is detailed in Section 2.3, for Group G0. Bleeding-induced anemia have been performed on 3 mice, each mouse was bled about 500 microliters with immediate injection of 500 microliters of intraperitoneal saline on three consecutive days (days 0, 1 and 2). Post-bleeding results were obtained with about 20 microliters of blood drawn on each day from day 3 [START_REF] Eymard | The role of spatial organization of cells in erythropoiesis[END_REF][START_REF] Rhodes | Stress reticulocytes lose transferrin receptors by an extrinsic process involving spleen and macrophages[END_REF]. These data are courtesy of Prof. Mark Koury, Vanderbilt University.

though the bleeding-induced anemia is more severe (minimum hematocrit value is about 30% of the initial hematocrit). This points towards different underlying physiological mechanisms associated to the anemia and its recovery.

We propose in this work a mathematical model of stress erythropoiesis in mice, based on the model of [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF], which explicitly includes an age-dependent description of the experimental protocol used to induce a severe anemia, consisting in injecting mice with two consecutive doses of phenylhydrazine. The objective of the current work is to complete our previous work [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF] by highlighting the respective roles of the experimental protocol (here phenylhydrazine-induced anemia) and physiological processes in the induction of the anemia and its recovery. To do so, we consider an age-structured model of erythropoiesis in which an age-dependent red blood cell mortality rate is introduced to account for phenylhydrazine effects. This model is presented in Section 2.1. We introduce the main features of the experimental data of phenylhydrazine-induced anemia in Section 2.3.

We then use this framework to investigate in Section 3 the respective influences of phenylhydrazine injections and of the physiological processes on the induction of the anemia (its strength and its duration) and the organism recovery.

Material, Methods, and Models

Mathematical Model

We introduce a mathematical model of stress erythropoiesis in mice, based on age-structured nonlinear partial differential equations describing erythroid progenitor and erythrocyte dynamics submitted to phenylhydrazine (PHZ) injections inducing anemia. This model is similar to the model introduced in [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF], except for its description of PHZ injections.

Let consider a population of erythroid progenitors and a population of erythrocytes, characterized by cell age a and the time of the observation t (see Figure 2). Among progenitors, two populations are considered: self-renewing progenitors, whose population is denoted by s(t, a), and differentiating progenitors, denoted by p(t, a). All progenitors are supposed to be in the bone marrow and subject to death by apoptosis. The duration of the differentiating progenitor compartment is denoted by τ p , and the duration of one self-renewing cycle by τ c , with τ c < τ p . Erythroid progenitors are assumed to differentiate in erythrocytes when they reach the age a = τ p . The number of erythrocytes, circulating in blood, is denoted by e(t, a), and E(t) denotes the total number of erythrocytes at time t, defined by

E(t) = +∞ 0 e(t, a)da.
Denote by α the erythroid progenitor apoptosis rate, and by σ the erythroid progenitor self-renewal rate. The rate α is assumed to be an increasing function of E, since the more erythrocytes the less erythropoietin [START_REF] Koury | Erythropoietin retards DNA breakdown and prevents programmed death in erythroid progenitor cells[END_REF], and erythropoietin inhibits erythroid progenitor apoptosis [START_REF] Koury | Erythropoietin retards DNA breakdown and prevents programmed death in erythroid progenitor cells[END_REF]. The rate σ is assumed to be a decreasing function of E: the more erythrocytes, the more glucocorticoids and the more glucocorticoids the less erythroid progenitor self-renewal (Bauer er al., 1999).

Denote by γ the mortality rate of erythrocytes. This rate depends on erythrocyte age a and time t when PHZ induces anemia, and is constant otherwise. In [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF], we noticed that in order to correctly reproduce experimental data, a modification of the mortality rate of erythrocytes had to be accounted for, as a consequence of using PHZ to induce anemia: one effect of this substance is to dramatically alter the lifespan of erythrocytes [START_REF] Berlin | Life span of the red blood cell of the rat following acute hemorrhage[END_REF][START_REF] Nagai | Studies on the short-lived reticulocytes by use of the in vitro labeling method[END_REF][START_REF] Nagai | Studies on the erythroblast dynamics based on the production of fetal hemoglobin[END_REF][START_REF] Shimada | The maturation of reticulocytes. II. Life-span of red cells originating from stress reticulocytes[END_REF][START_REF] Stohlman | Humoral regulation of erythropoiesis. VII. Shortened survival of erythrocytes by erythropoietin or severe anemia[END_REF]. Such a change in lifespan could be due to specific membrane properties of the erythrocytes produced during stress erythropoiesis [START_REF] Walter | Abnormal membrane surface properties during maturation of rat reticulocytes elicited by bleeding as measured by partition in two-polymer aqueous phases[END_REF]. To account for this effect, it is necessary to explicitly describe PHZ injections: we hence define a nonnegative function phz that describes increases in the mortality rate γ following PHZ injections.

The function phz first depends on the age of erythrocytes a and the time of the observation t: erythrocytes already circulating in blood when PHZ is injected are affected by PHZ (they either die or are damaged), while erythrocytes produced after the injection are not affected by the PHZ and consequently their mortality rate is not modified. In addition, the function phz is characterized by the time of the injection, denoted by t i , and by the duration of the PHZ effect in the blood (associated with the PHZ clearance age time

t i t i +t f phz = K+R phz = R t i +t f -t t f K +R R Figure 3: Function (a, t)
→ phz(a, t, t i , t f ), for fixed t i and t f . The grey areas correspond to phz = 0. Along the vertical dashed lines, phz is equal to its maximum (phz = K + R, on the left line) or to its minimum (phz = R, on the right line).

rate), denoted by t f . We define, for a, t ≥ 0,

phz(a, t, t i , t f )=          0, t ≤ t i or a ≤ t -t i , K t i + t f -t t f + R, t i ≤ t ≤ t i + t f and a ≥ t -t i , R, t i + t f ≤ t and a ≥ t -t i .
(1)

The function phz is also characterized by a nonnegative residual effect R, and a positive constant K which determines the maximum effect of the PHZ injection on the mortality rate of erythrocytes. A positive residual effect R means that erythrocytes that did not die following the PHZ injection have a larger mortality rate than cells that never encountered PHZ due to membrane damages, even though the organism is cleared with PHZ [START_REF] Walter | Abnormal membrane surface properties during maturation of rat reticulocytes elicited by bleeding as measured by partition in two-polymer aqueous phases[END_REF].

A schematic representation of the effect of the phz function is illustrated in Figure 3.

In order to compare the model simulations with experimental data (see Sections 2.3 and 2.4), we reproduce the experimental protocol consisting in two injections, following closely one another, with the initial injection occurring at time t i = t 1 = 0 and the second injection at time t i = t 2 > 0.

Since experimentally both injections are similar (same dose, same route of injection), we assume that they both have the same effects and we use the same parameter values for the function phz, except for the value of t 2 .

We use the notation γ(t, a) to stress the time and age-dependency of the PHZ-based experimental protocol, with

γ(t, a) = (1 + phz(a, t, 0, t f ) + phz(a, t, t 2 , t f )) γ, (2)
where γ stands for the average erythrocyte mortality rate.

Then, the quantities p, s and e satisfy the following system, for t > 0,

     ∂ t p(t, a) + ∂ a p(t, a)=-[α(E(t)) + σ(E(t))] p(t, a), a ∈ (0, τ p ), ∂ t s(t, a) + ∂ a s(t, a)=-α(E(t))s(t, a), a ∈ (0, τ c ),
∂ t e(t, a) + ∂ a e(t, a)=-γ(t, a)e(t, a), a > 0.

(3)

Boundary conditions associated with (3) describe cell flux between compartments,

         p(t, 0) = HSC + 2s(t, τ c ), s(t, 0) = τp 0
σ(E(t))p(t, a)da, e(t, 0) = Ap(t, τ p ).

(4)

New erythroid progenitors come from the division of self-renewing progenitors and the differentiation of hematopoietic stem cells: HSC denotes the constant flux of hematopoietic stem cells differentiating in erythroid progenitors. Self-renewing progenitors are produced at a rate σ. Erythrocytes are produced from mature (age a = τ p ) progenitors, and A denotes a constant amplification coefficient accounting for divisions of mature progenitors. For instance, A = 2 d where d is the number of differentiation compartments during the reticulocyte stage [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF].

Existence and uniqueness of solutions of system (1)-( 4) follow from the classical theory of age-structured equations [START_REF] Webb | Theory of nonlinear age-dependent population dynamics, Monographs and textbooks in Pure[END_REF], under conditions of smoothness of the nonlinear feedback functions α and σ. In addition one can show that system (3)-(4) has a unique steady state (a constant solution with respect to time t) under appropriate conditions (see Angulo et al. (2017), for a more general model). However, it must be noted that the complexity of the mathematical model limits the theoretical analysis that can be performed.

Numerical Simulations

The effort made on increasing the realism in the model is achieved at the expense of loss in mathematical tractability. Let's point out that, without additional restrictive assumptions, system (1)-( 4) cannot be solved analytically. Therefore the use of efficient methods providing a numerical approach is the most suitable mathematical tool for studying the model and, indeed, it is often the only one available. Besides, numerical methods have been successfully applied to structured models to replicate available field and/or laboratory data for a variety of different systems (e.g. Angulo and López-Marcos (2012); Angulo et al. (2011aAngulo et al. ( ,b, 2013a,b),b)). We completely describe in Angulo et al. (2017) the explicit second order numerical scheme built "ad hoc", which has been adapted to obtain the solutions of system ( 1)-(4).

When performing numerical simulations, we use the functions

α(E) = C α E nα θ nα α + E nα , and σ(E) = C σ θ nσ σ θ nσ σ + E nσ , (5)
where C α , C σ > 0 are the maximal apoptosis and self-renewal rates, respectively, θ α , θ σ > 0 threshold values, and n α , n σ > 1 sensitivity parameters.

In order to compare the model's outputs with experimental data (see Section 2.3), the numerical solution of system (1)-( 4) is used to compute the simulated hematocrit HCT , using the formula [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF])

HCT (t) = E(t) E(t) + E * (1 -HCT * )/HCT * , (6)
where E(t) is the total number of erythrocytes and HCT * is the steady state value of the hematocrit, obtained through experimental data.

Experimental Data

In order to get some insights into stress erythropoiesis, experiments consisting in inducing an anemia in mice and monitoring the recovery have been performed. Let us remind the reader that "hematocrit" is a test that measures the volume of red blood cells in a blood sample. It gives a percentage of erythrocyte volume found in the whole blood volume. Since platelet and white cell volumes are negligible within a blood sample, we assumed without loss of generality that a blood sample is mainly composed with erythrocytes and plasma.

Figure 4 shows the time evolution of the hematocrit of 5 batches of adult mice (groups G0 to G4), over 17 days, being rendered anemic by two consecutive intraperitoneal injections of PHZ with different doses. Even though the 5 batches display quantitative differences, they all reproduce the main basic features: Initially at its steady value, between 45 and 55 %, the hematocrit shows a strong fall following the anemia (the hematocrit reaches very low values, about 23%±3%, 2 to 3 days after the first injection for G0), then rapidly increases to reach a high value (about 55%±4%, 9 days after the first injection, for G0), and then returns to the equilibrium. Note that data for group G0 come from [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF].

Experiments performed on groups G0 to G4 vary in the dose of PHZ used to induce the anemia: G0 mice were injected with a 60mg/kg body weight dose, while G1, G2, and G3 mice were injected with a 30 mg/kg body weight dose and G4 mice with a 15 mg/kg body weight dose. They also differ with respect to the time of the second injection: G0, G1 and G4 mice were injected on days 0 and 1, G2 mice on days 0 and 0.7 (16h), and G3 mice on days 0 and 2. See Figure 4.C for a summary of these protocols.

In order to characterize the strength of the PHZ-induced anemia, we introduce the following quantities:

-HCT min , the minimum value of the hematocrit (HCT min ≈ 23 ± 3% for group G0), -t min , the time after the initial injection at which the minimum of the hematocrit is reached (t min ≈ 3 days for group G0), -t rec , the recovery time, which is the first time the hematocrit reaches a

given value (we chose the value 40%, which is just below the equilibrium value), after the hematocrit reached a minimum (t rec ≈ 6 days for group G0).

The first quantity directly measures the strength of the anemia, because the smaller the hematocrit the higher the probability that the organism cannot recover (and then eventually dies). The second quantity relates to the initial resistance of the organism to the anemia, and measures the speed at which the hematocrit decreases. The third quantity is related to the ability of the organism to recover from the anemia, and it measures the velocity at which the hematocrit recovers. A "strong" anemia will then be characterized first by a low minimum of the hematocrit, and second by a short time to reach this minimum. Potentially, one can also add to this definition a long time of recovery and a slow speed of recovery, while keeping in mind that the organism has to react quickly otherwise its survival can be questioned. A "weak" anemia on the contrary will be characterized by a high minimum of the hematocrit.

Parameter Estimation

We proceed in two steps. We first use data from group G0 to estimate parameter values and evaluate the contribution of each parameter to the quality of the results. Second, we investigate the robustness of the model by confronting it, and estimated parameter values, to data from groups G1 to G4 (see Section 2.3).

In the first step, we first estimate parameter values so that the model correctly reproduces experimental data from group G0. We simulate the model for fixed parameter values and compute the weighted residual sum of squares using G0 data. The weighted residual sum of squares is given by

RSS := N i=1 (SimulatedHCT i -ExperimentalHCT i ) 2 σ 2 i ,
where N is the size of the experimental sample, SimulatedHCT i the value of the simulated hematocrit at the i-th experimental time, ExperimentalHCT i the corresponding average hematocrit value, and σ i the experimental standard deviation. The same formula is used for groups G0 to G4, using each time corresponding data. Identifiability of estimated parameter values [START_REF] Raue | Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood[END_REF] is deduced from the boundedness of 95% confidence intervals, hence estimated parameter values indeed correspond to a global minimum of the residual sum of squares.

This leads us to identify parameter values that best reproduce G0 data, and situations that favor quantitatively relevant dynamics. We then focus on the contribution of some parameters to the severity of the anemia and the recovery's features. The complexity of the model limits the numerical analysis, so it would not be reasonable to try performing a sensitivity analysis on every parameter of the model. This is why we focus on the influence of specific parameters, which characterize the main feedback functions of the system, on the model's features. To that aim, we vary them individually in a given range while fixing the others. For each parameter, the range used as well as the paragraph presenting the results are indicated in Table 1. In the second step, we use the estimated parameter values to confront our model to data from groups G1 to G4, in order to evaluate its robustness in describing different dynamics obtained with very similar experimental protocols. In order to characterize each group dynamics, we also compute the weighted residual sum of squares.

While most parameter values will be estimated to reproduce experimental data, some values will be fixed throughout the parameter estimation procedure, both in steps 1 and 2 (these parameters are however varied for investigating the potential influence of their variation on the system's behavior).

The average lifespan of an erythrocyte is known to be 40 days in mice [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF], so the average mortality rate would be γ = 1/40 d -1 .

Values of the amplification coefficient of progenitors into erythrocytes A, and of HSC flux into the erythroid lineage HSC, are taken from Crauste et al.

(2008), A = 2 8 and HSC = 10 4 cells.vol -1 respectively. Also, the progenitor compartment duration τ p and the cell cycle duration τ c are fixed to τ p = 4 days and τ c = 1 day [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF].

Additionally parameter K in (1) is fixed to K = 20 d -1 , after multiple tests showed no real influence of variations of K. Finally, the simulation step h is fixed throughout all simulations to h = 1/128 day (see Angulo et al. (2017)).

Results

We use the model presented in Section 2 to simulate a PHZ-induced anemia in mice and the response of the organism, and we mainly focus on two points. First, we investigate the influence of the experimental protocol upon the features (strength, duration) of the anemia. Second, we evaluate the roles of physiological processes regulating erythrocyte production. To do so, we compare results of the simulations of system ( 1)-( 4) with experimental data presented in Section 2.3, group G0. We then show that our model is able to describe other experimental data (groups G1 to G4) and draw conclusions about PHZ-induced anemia.

Data

Hematocrit evolves similarly for groups G1, G2 and G3 (see Figure 4.A).

Quantitative differences are however important between the hematocrits of groups G0, G1, G4, for which the injected doses of PHZ are different (see these high values appear to be caused by a strong reaction of the organism that needs to quickly produce new red blood cells in order to ensure survival associated with a lag (or inertia) in the regulation of progenitor apoptosis that does not compensate fast enough this production of new red blood cells, resulting in an overshoot of the hematocrit.

Influence of PHZ injections

In order to investigate the influence of the injection protocol on the dynamics of stress erythropoiesis and blood recovery, we focus on 3 parameters -t 2 , the time of the second injection (the initial injection occurs at t = 0). System (1)-( 4) is numerically solved and the corresponding simulated hematocrit HCT is computed using (6). The simulation that best reproduces experimental data for group G0 is shown in Figure 5. It has been obtained with R = 5 d -1 , t f = 3 days, t 2 = 1 day, while other parameter values are given in Table 2. This set of parameter values for (R, t f , t 2 ) is characteristic of a "good" situation: the model correctly reproduces experimental data when

t f ≥ 2 days and -either the residual effect R is large (R ≥ 3 d -1 ), -or R is small (1 d -1 < R ≤ 3 d -1
), provided that the smaller R the larger t f , and t 2 is not too large (t 2 ≤ 1.5 days). These situations are a priori biologically realistic. The first one corresponds to a strong effect of PHZ on erythrocytes that survived the first injection, while the second one corresponds to a second injection following quickly the first one in order to sustain the action of PHZ. One may note that, in the experiments, t 2 = 1 day which is in agreement with our estimation, while no information on neither R nor t f is available.

When R = 0 or R is close to zero, the model is unable to correctly describe the fall of the hematocrit that follows the two injections, and varying K in (1)

does not improve the results, since the organism always reacts much faster in this case than what is experimentally observed (data not shown). This indicates that the residual effect plays an important role in inducing severe anemia, following PHZ injections, as suggested in [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF].

We further investigate the roles of R, t f and t 2 by focusing on specific features of the anemia: HCT min , t min and t rec .

Parameters R and t f have a major influence on the strength of the ane- 2.

A 1 2 3 4 5 6 7 8 R HCT min t f =1 t f =2 t f =3 t f =5 B 1 2 3 4 5 t f HCT min R=1 R=3 R=5 R=8
mia. The minimal hematocrit value, HCT min , is mainly controlled by R and t f (Figure 6): the larger R and t f , the lower HCT min and the more severe the anemia. In addition, combined increases of R and t f strengthen the anemia by more strongly decreasing HCT min than each parameter separately. This complements the above conclusion regarding the obtention of good reproduction of data: for large values of both R and t f the simulated solution correctly reproduces experimental data, and in particular it correctly captures the value of the hematocrit at its minimum.

Second, the time at which the second injection occurs, t 2 , is the main parameter controlling the time at which the hematocrit reaches its minimum value, t min (Figure 7): t min increases as t 2 increases. This influence of t 2 on t min indicates that the organism does not react fast enough to the first injection to limit the action of the second injection. The residual effect R

shows very limited influence on the value of t min . The parameter t f shows no influence on t min , except for situations that have been identified as unable to generate a correct response to the anemia (small values of t f , R and t 2 ). In conclusion, t 2 strongly influences t min , whereas R and t f have either limited or no influence on t min .

None of these three parameters is however significantly influencing the recovery (Supp. Fig. 1 to 3), which remains fairly constant for all values of the three parameters (around 4.5-5.5 days). We then hypothesize that the recovery is mainly determined by the response of the organism to the anemia, and the investigation of the roles of physiological processes in the next sections will indeed confirm this hypothesis.

Apoptosis and Differentiation Processes drive the Recovery

The apoptosis rate α is characterized by three parameters: C α , θ α and 2.

Differentiation of erythroid progenitors is characterized by the differentiation time τ p . Figure 8.B shows that long differentiation times first contribute to strong anemias and second, slow down recovery. The differentiation time τ p , although not dependent on a feedback control in this model (no biological evidence), appears to play a role almost as important as the apoptosis rate in controlling the recovery. Shortening the differentiation process might be an efficient way to rapidly recover from a strong anemia. However, to our knowledge, a reduction of the duration of the differentiation process has never been experimentally observed in stress erythropoiesis, contrary to a control of the apoptosis rate by the number of circulating red blood cells.

Comments on the Role of the Self-Renewal Rate

Let us focus on the contribution of the self-renewal rate σ to the organism's recovery from an anemia. The rate σ is characterized by three parameters C σ , θ σ and n σ , First, C σ has no significant influence on HCT min and t min (Supp. Fig. 5), and consequently does not contribute to the severeness of the anemia. Second, for small values of C σ (< 0.5 d -1 ) one observes fast recoveries, characterized by t rec ≤ 5 days, while larger values of C σ are associated with slow recoveries (see Figure 9). Other parameters associated with the self-renewal rate have a very limited influence (θ σ ) or no influence at all (n σ ) on both the anemia and the recovery (Supp. Fig. 6 and7).

Even though self-renewal appears to be mandatory to recover from the anemia (in agreement with [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF]), these results indicate that the regulation of the self-renewal rate (characterized by parameters θ σ and n σ ) does not play a crucial role in the recovery and is potentially dispensable. Consequently, we may hypothesize that the self-renewal rate is indeed constant among the erythroid progenitor population, equal to C σ .

Comments on Stem Cell Regulation

The erythropoiesis process depends upon the production of immature erythroid progenitors by hematopoietic stem cells. This process is highly controlled (see for instance Pujo-Menjouet et al. (2005); [START_REF] Roeder | A novel dynamic model of hematopoietic stem cell organization based on the concept of within-tissue plasticity[END_REF]), even though we did not model it in details since it is not part of the erythropoiesis process. We can however focus on the influence of the parameter associated with a continuous and constant production of erythroid progenitors, HSC. Variations of this parameter show important modifications of HCT min and t rec (Figure 10), and of t min (Supp. Fig. 8), indicating a significant influence of HSC both on the strength of the anemia and on the recovery.

Even though we cannot exclude that the entire recovery process could be controlled by the regulation of HSC dynamics, strong variations of the flux of HSC entering the erythroid lineage are usually associated with hematological disorders or diseases [START_REF] Mackey | Unified hypothesis of the origin of aplastic anaemia and periodic hematopoïesis[END_REF][START_REF] Pujo-Menjouet | Contribution to the study of periodic chronic myelogenous leukemia[END_REF].

Hence, only small variations of the parameter HSC should be considered, and in this case the contribution of the regulation of the HSC compartment would be much less significant and essential to stress erythropoiesis.

Moreover, the regulation of HSC dynamics, leading to a modification of the flux of newly created erythroid progenitors, would not have instantaneous consequences on the erythropoiesis process: it would take several divisions, hence several days, to significantly modify the flux of stem cells differentiating in erythroid progenitors. One can observe on Figure 5 that the organism recovers 2 to 3 days after reaching its lowest level, a priori excluding the HSC regulation as the only mechanism allowing a fast recovery. Group

R (d -1 ) t f (d) G0 5 3 G1 0 2.2 G2 6 1 G3 0 3 G4
4.5 0.2

Reproducing various PHZ-induced anemia protocols

The above-mentioned analysis has been performed with parameter values that were estimated to reproduce experimental data from group G0. One may question the robustness of these values and their relevance in describing stress erythropoiesis dynamics in other PHZ-induced anemia protocols.

Groups G1 to G4 display similar features than group G0: fall of the hematocrit following PHZ injections, then a recovery period with a return to a steady hematocrit by day 17 post-injection ; yet they quantitatively behave differently depending on the dose (G1, G2, G3 vs G4) or the time between two injections (G1, G4 vs G2 vs G3). Since parameters R and t f have been identified in Section 3.2 as responsible for the strength and duration of the anemia, we focus on their contribution in reproducing experimental data.

To investigate the model's behavior when used to reproduce data in groups G1 to G4, we first fix all parameter values associated with physiological processes to values obtained with the G0 group (Table 2), and we only vary parameters associated with PHZ injections: R and t f . Value of the parameter t 2 is set according to the group (t 2 = 1 day for G1 and G4, t 2 = 0.7 day for G2, and t 2 = 2 days for G3). Parameters R and t f are determined to fit the data (see Table 3). Results are presented in Figure 11 by dashed lines.

For each group, estimating only values of R and t f does not provide correct reproduction of experimental data: although the anemia is rather properly described, the recovery is mostly overestimated whatever the group (too fast for all groups, and reaching values higher than experimental ones for groups G1, G2 and G4), and the simulated dynamics of the hematocrit are not in agreement with experimental observations.

In a second step, in addition to parameters R and t f , we allow parameters Figure 11: Estimation of hematocrit dynamics for groups G1 (A), G2 (B), G3 (C) and G4 (D). For each figure, the value of t 2 is fixed according to the experimental protocol. The solid line is the best fit of the data, obtained by varying all parameters (Table 4), while the dash line is the fit of the data obtained by only varying R and t f (Table 3), other parameter values being given by Table 2. Table 4: Estimated values of R, t f , and physiological parameters, for groups G1 to G4. Parameter values for group G0 are recalled for comparison purposes (in gray).

Group

R t f C α θ α n α C σ RSS (d -1 ) (d) (d -1 ) (cells.vol -1 ) (N.U) (d -1 ) (N.U) G0
5.0 3.0 20 10 6.5 11 0.7 6.0 * G1 8.0 0.6 55 10 6.2 24 0.7 2.1 G2 8.0 0.8 21 10 6.0 37 0.9 1.6 G3 0.5 2.0 4.5 10 6.2 2 0.7 0.9 G4 5.0 0.2 68 10 6.1 58 1.1 2.9 * This value may seem much higher than the other ones in the same column, but it has simply been computed using more experimental data (see Figure 4.B).

related to the apoptosis rate (C α , θ α , n α ) and to the self-renewal rate (C σ )

to vary, with the assumption that the self-renewal rate is constant. We then obtain the solid curves in Figure 11 and values given by Table 4. All results show a very good agreement of the simulation and the experimental measurements (see column 'RSS' in Table 4), variations of the hematocrit during the recovery phase (between day 3 and day 15 post-injection) are particularly well captured by the model. This stresses the ability of the model to describe stress erythropoiesis induced by various protocols of PHZ injections.

One may however notice that correct reproduction of the data is obtained for parameter values that are very different from one group to the other (see Table 4).

Regarding physiological rates, values of the self-renewal rate C σ are consistent for all groups, with a doubling time between 7.5h and 12h. On the contrary, the apoptosis rate shows very different values of its three characteristic parameters among all groups. A careful investigation of the values of the apoptosis rate all along the response to PHZ injections first shows that for all groups the maximal value C α is never reached. Second, contrary to what would be expected from the shape of the function α (see ( 5)), the apoptosis rate does not behave as a saturating function for groups G1, G2 and G4, but rather as an exponential function (see Figure 12), with values ranging between 0 (when the number of erythrocytes is high, fully inhibiting progenitor apoptosis) and a maximum value in the interval [3 d -1 ; 12 d -1 ], 4.

by the time of the second injection (2 days after the initial injection), we may hypothesize that the 2-days period between PHZ injections allowed the organism to partly recover and in particular to perform a tighter control of erythroid progenitor death by apoptosis, that resulted in a slower recovery and a less strict regulation of Epo-mediated apoptosis.

Regarding values of parameters R and t f , a trend appears on the way they are influenced by the dose or the time of second injection (Figure 13): the higher the dose the higher t f , and to a lesser extent the sooner the second injection the higher R. Yet in each case a lot of variability is measured for groups with either the same dose or the same time of the second injection.

Since R and t f both are associated to the influence of PHZ on erythrocyte death, one parameter may compensate the influence of the other (R and t f have been identified as key regulators of the value HCT min in Section 3.2 without clearly identifying their respective roles), resulting in the observed variability. In order to check that assumption, we drew the product R × t f (adimensionalized quantity) as a function of the dose and the second time of the injection (Figure 14). It clearly shows an influence of the experimental protocol in the estimated values of R and t f : the higher the dose, the larger the product R × t f , and the sooner the second injection, the larger R × t f . This is particularly true when the dose (respectively, second injection time) effect is observed for a shared second injection time (respectively, dose): this is shown by squared dots in Figure 14, and the trend clearly appears.

Therefore, despite differences in parameter values between all groups, estimated values account for the various experimental protocols used to generate the data (groups G1 to G4) and the associated physiological responses. In addition, the influence of the experimental protocol clearly appears through the value of the product R × t f that sums up the specifics of PHZ injections.

Discussion

We modified a previously published model of stress erythropoiesis in mice to explicitly describe and account for the specific phenylhydrazine-based way of inducing anemia experimentally. Contrary to bleeding-induced anemias for instance, phenylhydrazine-induced anemias not only reduce the quantity of red blood cells but also affect the entire red blood cell production process: indeed, phenylhydrazine either kills red blood cells -hence inducing an anemia -or damages red blood cell membranes and then shortens red blood cell lifespans with potential consequences on the organism behavior on the days following the anemia. Although the model cannot be solved analytically, we implemented an appropriate numerical method (Angulo et al., 2017) in order to numerically identify the respective influences of the experimental protocol and physiological regulatory process on the model's behavior. This gave us insights into the way the organism deals with phenylhydrazine injections inducing anemia and the recovery from the anemia.

First, our model proved its ability to reproduce experimental data consisting in hematocrit measurements at several time points following phenylhydrazine-induced anemia, with relevant estimated parameter values. We showed that the strength of the anemia, induced by phenylhydrazine injections, is mainly controlled by the experimental protocol (the quantity of injected phenylhydrazine, the number of injections, and the time between the two injections). The fall of the hematocrit observed in the days following phenylhydrazine injections is strongly dependent upon three parameters, R and t f on one side, and t 2 on the other side, suggesting different roles and targets of these parameters in the induction of the anemia. Parameter t 2 characterizes the experimentalist's choice (the time at which the second injection occurs) while R and t f are mainly specific of phenylhydrazine's properties (its strength in inducing death and how long its effects are felt by the organism). On the contrary, we showed that the recovery phase is mainly controlled by the regulation of physiological processes. Epo-mediated inhibition of erythroid progenitor apoptosis, identified as the main regulator of erythropoiesis for years [START_REF] Koury | Erythropoietin retards DNA breakdown and prevents programmed death in erythroid progenitor cells[END_REF], is a powerful mean to speed up the recovery, because it can quickly adapt the production of erythrocytes to strong modifications in circulating red blood cell counts.

Second, our model showed its ability to reproduce similar experimental data, obtained with modified protocols: either the dose of phenylhydrazine This is hardly the case, as HSC regulation is performed through several cell cycles and important variations of the HSC population are not observed over a short period of time in healthy individuals. These additional contributions to the organism's response to anemia hence do not appear as relevant as the regulation of erythroid progenitor apoptosis.

In [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF], we concluded on the importance of accounting for erythroid progenitor self-renewal when describing stress erythropoiesis, and described a nonlinear erythrocyte-mediated regulation of self-renewal.

Our results confirm the relevance of accounting for erythroid progenitor selfrenewal, all parameter value estimation procedures leading to positive selfrenewal rates, whatever the experimental protocol. Nevertheless, our results also indicate that regulation of the self-renewal rate is not necessary to account for appropriate dynamics: a constant self-renewal rate, corresponding to a robust doubling time around 10h, allows to reproduce experimental dynamics. This can be explained by the nature of the apoptosis rate regulation: erythroid progenitor apoptosis is highly regulated, and for normal red blood cell counts estimations of the apoptosis rate are rather high (see [START_REF] Crauste | Adding self-renewal in committed erythroid progenitors improves the biological relevance of a mathematical model of erythropoiesis[END_REF] for instance): in this case, a small population of self-renewing cells among the progenitor population can hardly be detected. During stress erythropoiesis, and in particular during anemia, the strong regulation of apoptosis can almost completely inhibit it, then making the self-renewal rate influential, without necessarily an opposite regulation (see Figure 15). Additionally, this can lead to further analysis of the model, in particular system (1)-( 4) can be mathematically analyzed when σ is constant whereas the analysis has been impossible, to this day, for nonlinear self-renewal rates σ.

Regarding the erythroid progenitor apoptosis rate, our analysis showed that its regulation could be slightly different than the one we implemented in our model, that is a saturating function characterized by very low values for small erythrocyte counts and very high values for high erythrocyte counts. As discussed in Section 3.6, values of the apoptosis rate measured through an entire simulation do not correspond to the spectrum that the saturating function could span: only the first increasing, convex part of this function is actually used to perform stress erythropoiesis. From a modeling point of view, this indicates that instead of using a classical saturating Hill function to describe the erythrocyte count-dependent apoptosis rate, an exponential increasing function could be used. This may not seem to represent an important contribution to erythropoiesis modeling, yet it would allow describing the apoptosis rate using two parameters, instead of three in the current version, hence favoring a better, more sensitive control of the apoptosis regulation.

Aknowledgements

OA has been partially supported by project MTM2014-56022-C2-2-P of the Spanish Ministerio de Economía y Competitividad and European FEDER Funds. We thank Quentin Lamy for his technical help in collecting the data.

We thank Prof. Mark Koury for his support and for providing the bleedinginduced anemia data in Figure 1.

Figure 2 :

 2 Figure 2: Schematic representation of the erythropoiesis model. Dashed lines represent feedback controls. See text for notations.

CFigure 4 :

 4 Figure4: Experimentally measured hematocrits. Hematocrits of groups (A) G1 (red), G2 (blue), G3 (green), and (B) G1 (red), G4 (black), G0 (blue). Average and standard deviations over 9 (group G0) or 8 (groups G1 to G4) adult mice are displayed on the 17 days period. (C) Characteristics of PHZ-induced anemia experiments: for each of the 5 groups, the time interval between two consecutive injections of PHZ and the dose that has been injected are presented.

Figure 4 .

 4 Figure 4.B). Both the fall of the hematocrit phase and the recovery phase display different features, which appear to be dose-dependent: distinct values of the minimum value of the hematocrit and different slopes of the recovery phase. Additionally, high values observed for group G0 around day 9 are not observed for the other groups. As mentioned in Crauste et al. (2008),

(Figure 5 :

 5 Figure5: Experimental and simulated hematocrits for group G0. Experimental data are given by the red crosses with bars (mean value ± standard deviation). The simulated hematocrit is given by the black curve. All estimated parameter values are from Table2.

Figure 6 :

 6 Figure 6: Influence of the values of R and t f on the values of HCT min . Values of HCT min are computed for t 2 = 1 day and (A) R ∈ [0, 8] d -1, with t f = 1 day (blue curve), t f = 2 days (red curve), t f = 3 days (black curve), t f = 5 days (green curve), and (B) t f ∈ [0.2, 5] days, with R = 1 d -1 (blue curve), R = 3 d -1 (red curve), R = 5 d -1 (black curve), R = 8 d -1 (green curve). Black dash-dotted lines indicate the range of experimental values of HCT min (from Figure5, group G0). Other parameter values are given by Table2.

Figure 7 :

 7 Figure7: Influence of the value of t 2 on the values of t min . Values of t min are computed for t 2 ∈ [0.1, 1.9] days, with R = 1 d -1 (blue curve), R = 3 d -1 (red curve), R = 5 d -1 (black curve), R = 8 d -1 (green curve), and t f = 3 days. Other parameter values are given by Table2.

Figure 8 :Figure 9 :

 89 Figure 8: Influence of θ α and τ p on HCT min and t rec . (A) Influence of θ α on HCT min and t rec . (B) Influence of τ p on HCT min and t rec . In both cases, values of HCT min are displayed on the left axis and correspond to the plain curve, while values of t rec are displayed on the right axis and correspond to the dashed curve.

Figure 10 :

 10 Figure 10: Influence of HSC on HCT min and t rec . HCT min values are on the left axis and represented by the plain curve, while t rec values are on the right axis and represented by the dashed curve. Note that HSC is represented in logarithmic scale on the x-axis.

Figure 12 :Figure 13 :

 1213 Figure12: Normalized values of α(E) over the anemia and the recovery phases, for groups G1 to G4. On the x-axis, values of E taken during the simulation shown in Figure11have been normalized (0 represents the minimum value, 1 the maximum one). On the y-axis, corresponding values of α(E) have also been normalized.

  or the time of the second injection were modified in 4 different experiments, and qualitatively similar yet quantitatively different dynamics were generated. Model simulations also managed to characterize each experimental protocol by highlighting relationships between the injected dose or the time of the second injection and an adimensionalized variable R × t f that embeds phenylhydrazine properties. This allows to directly link the experimental protocol to the severeness of the anemia via the action of phenylhydrazine.Consequently, respective influences of the experimental protocol and physiological processes on the features of stress erythropoiesis can be separately investigated. Moreover phenylhydrazine-induced anemias can be compared even though the experimental protocol changes (different dose, different time of the second injection), but can also be compared to other experimental protocols (for instance bleeding) thanks to the characterization of the protocol by the combination of two variables, R and t f . Our description of the influence of phenylhydrazine on erythrocyte dy-namics could be improved. Pharmacodynamics and pharmacokinetics of phenylhydrazine could for instance be incorporated to the description of the induction of the anemia, to directly relate features of the induced anemia to phenylhydrazine properties, instead of indirect properties (residual effect, etc.) as we did in this work. It is however difficult to properly define the PK-PD of phenylhydrazine, as few information is accessible and so developing a relevant PK-PD model of phenylhydrazine could represent an entire research work.Our study stresses the potential roles of additional processes: shortening or lengthening of differentiation times (characterized by τ p in our model) and increase or decrease of the stem cell flux from the hematopoietic stem cell compartment could also strongly influence both the effect of phenylhydrazine on the organism and the recovery(Figures 8.B and 10). However, to our knowledge, no experimental evidence ever showed that differentiation times (or cell cycle durations) were modified during stress erythropoiesis.Therefore, such a control of the organism's response to anemia remains hypothetical at this stage. Regarding the influence of the differentiation of HSC into erythroid progenitors, our results show that this could influence the recovery provided that large variations of the HSC flux are allowed, and that these changes can occur on a very short time scale (almost instantaneously).

Figure 15 :

 15 Figure 15: Two scenarii for explaining regulatory process during stress erythropoiesis. (A) The self-renewal rate depends on the erythrocyte count, see (5) (B) The self-renewal rate is constant. In both cases the apoptosis rate depends on the erythrocyte count, as in (5). Normal erythropoiesis corresponds to high erythrocyte counts, while anemia corresponds to low erythrocyte counts. The gray area corresponds to an unclear, a priori transient zone of low erythrocyte counts that does not necessitate a strong reaction from the organism.

  A.S. Ackleh, H.T. Banks, K. Deng, A finite difference approximation for a coupled system of nonlinear size-structured population, Nonlinear Analysis 50 (2002) 727-748, doi:10.1016/S0362-546X(01)00780-5 A.S. Ackleh, K. Deng, K. Ito, J. Thibodeaux, A structured erythropoiesis model with nonlinear cell maturation velocity and hormone decay rate, Math. Bios. 204 (2006) 21-48, doi:10.1016/j.mbs.2006.08.004. A.S. Ackleh, J. Thibodeaux, Parameter Estimation in a Structured Erythropoiesis Model. Math. Bios. Eng. 5 (2008) 601-616. M. Adimy, F. Crauste, S. Ruan, Modelling hematopoiesis mediated by growth factors with applications to periodic hematological diseases. Bull. Math. Biol. 68 (8) (2006) 2321-2351. O. Angulo, F. Crauste, J.C. López-Marcos, Numerical integration of an erythropoiesis model with explicit growth factor dynamics, Journal of Computational and Applied Mathematics (2017) https://doi.org/10.1016/j.cam.2017.01.033. O. Angulo, J.C. López-Marcos, M.A. Bees, Mass Structured Systems with Boundary Delay: Oscillations and the Effect of Selective Predation, J. Nonlinear Sci. 22 (2012) 961-984. O. Angulo, J.C. López-Marcos, M.A. López-Marcos, Numerical approximation of singular asymptotic states for a size-structured population model with a dynamical resource, Math. Comput. Model. 54 (2011a) 1693-1698.

Table 1 :

 1 Parameter sensitivity ranges. For each parameter of interest (column 'Parameter') the prospective interval in which its influence on the dynamics of the model is investigated appears in the column 'Range', and the subsection label in which these results are presented in the column 'Paragraph'. The upper part of the table deals with the 9 parameters whose values are also estimated to fit the data, while the bottom part of the table deals with the 4 parameters whose values are fixed throughout simulations.

	Parameter	Range Paragraph Parameter	Range Paragraph
	C α	[0; 70] not shown C σ	[0; 4]	3.4
	θ α	[10 4 ; 10 8 ]	3.3 θ σ	[10 4 ; 10 8 ] not shown
	n α	[1; 58] not shown n σ	[1; 35] not shown
	R	[0; 8]	3.2 t f	[0; 5]	3.2
	t 2	[0; 2]	3.2
	τ p	[1; 8]	3.3 τ c	[0; τ p ] not shown
	HSC	[10 2 ; 10 6 ]	3.5 A	[100; 500] not shown

Table 2 :

 2 Parameter values corresponding to Figure5. All parameter values listed in the table have been estimated: the best value is indicated in the column 'Value' while its 95% confidence interval is indicated in the column 'Confidence Interval'. Key: vol = volume unit; N.U. = no unit.

	Parameter Value Confidence Interval Unit
	Apoptosis Rate
	C α	20	[19.6; 20.4] d -1
	θ α	10 6.5	[10 6.4 ; 10 6.6 ] cells.vol -1
	n α	10.6	[10.4; 10.8] N.U.
	Self-Renewal Rate
	C σ	0.7	[0.63; 0.77] d -1
	θ σ	10 6.8	[10 6.6 ; 10 7.0 ] cells.vol -1
	n σ	6.2	[5.8; 6.6] N.U.
	PHZ Injections
	R	5	[4.9; 5.1] d -1
	t f	3	[2.95; 3.04] d
	t 2	1	[0.97; 1.03] d

Table 3 :

 3 Estimated values of R and t f for groups G1 to G4. Values of R and t f for group G0 are recalled for comparison purposes (in gray).

  Estimated values of R × t f as function of the dose (A) and t 2 (B). Within-boxes values share the same t 2 value (A) or dose (B). Values taken from Table4.

	A						B
		16						16
		12						12
	R x t f	8					R x t	8
		4						4
		0						0
		10	20	30	40	50	60	15	20	25	30	35	40	45	50
				Dose (mg/kg body weight)					t 2	(hours)

f Figure 14: