|Home | About | Journals | Submit | Contact Us | Français|
Most pharmacodynamic (PD) models of cellular response assume a time-invariant (i.e., constant) cellular disposition despite known changes in the disposition with time, such as the reticulocyte residence time in the systemic circulation during stress erythropoiesis. To account for changes in cellular disposition, a comprehensive PD model that involves endogenous erythropoietin (Epo), reticulocytes, and hemoglobin responses was developed in phlebotomized sheep that considers a time-variant reticulocyte residence time and allows for the simultaneous determination of changes in the cellular disposition and cellular production. Five sheep were phlebotomized to hemoglobin concentrations of approximately 4 g/dl. Epo concentrations, reticulocytes, and hemoglobin concentrations were frequently sampled for 5–7 days prior to and 25–30 days following the phlebotomy. Initial steady-state conditions were assumed and the time-variant reticulocyte residence time in the systemic circulation was semiparametrically represented using a constrained spline function. Hemoglobin production was modeled using a Hill function via an effect site compartment. The initial steady state reticulocyte residence time in the systemic circulation was estimated as 0.477 (0.100) (mean (SD)) days, which maximally increased 2.01- to 2.64-fold higher than the initial steady-state residence time 5.95 (0.899) days post-phlebotomy (P < 0.01). On average, the residence time returned to steady-state values 15.4 (2.36) days post-phlebotomy, which was not significantly different from the initial steady-state value (P > 0.05). The baseline hemoglobin production rate was estimated at 0.0929 (0.0472) g/kg/day and the maximum production rate under stress phlebotomy was estimated at 0.504 (0.0422) g/kg/day. These data indicate that endogenously released Epo under acute anemic conditions can increase hemoglobin production approximately 5-fold. The determined increase in reticulocyte residence time produced under stress erythropoiesis is similar to the commonly reported 2- to 3-fold increase observed in human patients.
Reticulocytes are immature red blood cells (RBCs) that initially reside in the bone marrow and subsequently are released into the systemic circulation where they develop into mature RBCs . Reticulocytes are produced from erythroid progenitor cells (BFU-E and CFU-E) located primarily in the bone marrow. They are commonly used to monitor the “real-time” functional state of the erythroid bone marrow and play an important role in diagnoses and management of anemia [2–7]. Measurement of the reticulocyte count is valuable for determining if a patient has a functionally normal response to either endogenously produced or exogenously administered erythropoietin (Epo) [5,8]. Reticulocyte counts are also used to monitor bone marrow suppression by chemotherapy  and bone marrow engraftment following bone marrow transplantation [10–13]. The maturation of erythroid progenitor cells into reticulocytes and ultimately RBCs is primarily controlled by Epo, a-35kD glycoprotein hormone produced by the pertibular cells of the kidney in response to oxygen need . During the development from erythroid progenitor cells their hemoglobin content increases until it develops into a reticulocyte upon nucleus extrusion, where further maturation primarily involves the removal of ribosomal RNA, remodeling of the plasma membrane, and a progressive decrease in cell size, increasing the concentration of hemoglobin on a mass/volume basis, but not increasing the overall amount of hemoglobin within the cell (approximately 30 picograms per cell in adult humans) [1,3,4]. Therefore, a reticulocyte cell has approximately the same amount of hemoglobin as a RBC .
Generally, under normal erythropoietic conditions in humans (i.e., non-anemic or non-erythropoietically stimulated), reticulocytes reside in the systemic circulation approximately 24 h before maturation into RBCs [3,4]. However, during stress erythropoiesis (i.e., stimulated erythropoietic conditions), the reticulocyte residence time in the circulation increases to an estimated 2–3 days . Reticulocytes produced under stress erythropoiesis also contain more residual ribosomal RNA, are larger, and have less flexible plasma membranes than those produced under normal conditions, and therefore are thought to be immature reticulocytes that under “normal” physiological conditions reside in the bone marrow until being released as more mature reticulocytes [2,4,14,15]. Accordingly, the reticulocyte counts increase under stress erythropoiesis not only due to increased reticulocyte production in response to Epo stimulation, but also due to a longer residence time in the systemic circulation. To account for the effect of this increased reticulocyte residence time in the systemic circulation on the reticulocyte count in clinical practice, one “rule of thumb” approach has been to divide the observed reticulocyte count or percentage by two to obtain a “more accurate” index of reticulocyte production when they are produced under stress erythropoiesis . More sophisticated quantitative approaches utilize the maturity of the produced reticulocytes, referred to as the reticulocyte maturity index (RMI). The RMI is based on of individual reticulocyte RNA content commonly determined by fluorescence intensity measurements. However, the problem with existing methods utilizing RMI is that there is no standard method that corrects for the observed reticulocyte counts based on the RMI. Furthermore, these methods suffer from the difficulty in standardizing and cross-correlating across different instruments and laboratories .
Although different terminology is often used to describe reticulocyte or cellular residence times, in this manuscript it is defined as the time a cell resides in the sampling pool (the systemic blood circulation). Analogous terminology used by other authors includes “lifespan” and “maturation time” [16–23]. Previous Epo-reticulocyte pharmacokinetic/pharmacodynamic (PK/PD) models have commonly assumed that all cells have a time invariant and common lifespan or maturation time (i.e., “point distribution”) [16–21]. Recently, Epo-reticulocyte PK/PD models have been presented that allow for a fixed (time invariant) distribution of cell lifespans . The time invariant lifespan assumption has been partially addressed by allowing for two different reticulocyte maturation times, a baseline and a stress-erythropoiesis maturation . Other models that partially address the stress erythropoiesis effects on increasing reticulocyte lifespan incorporate an erythropoietin dose dependent lifespan . However, the dynamic changes in reticulocyte residence time in the systemic circulation throughout the time course of stress erythropoiesis were not determined in any of the above studies. Given the known change in reticulocyte residence time in the systemic circulation following stress erythropoiesis, the clinical importance of the reticulocyte count/residence time change, and the limited attempts to quantify these changes; the primary objective of the present analysis was to dynamically determine the changes in reticulocyte residence time due to phlebotomy induced stress erythropoiesis. To achieve this aim, a time-variant cellular disposition model was developed and integrated into a comprehensive pharmacodynamic analysis of Epo, reticulocyte, and hemoglobin response in phlebotomized sheep.
The unit response (UR) function of a linear cellular disposition, (Fig. 1, panel A), for hematopoietic cells may be described as:
where ‘a’ is the time it takes for an individual cell to appear in the systemic circulation (the sampling compartment) following progenitor stimulation, ‘b’ is the time the subsequent cell exits, and ‘t’ is measured, as ‘a’ and ‘b’, relative to the time of stimulation (time ‘s’ in Fig. 1, panel A). ‘U’ is the unit step function defined by:
Therefore, UR(t) takes on the value of 1 when the cell is present, and a value of 0 when the cell is not present in the sampling compartment.
Consider hematopoietic progenitor cells, which in the following will be referred to simply as “cells”, that are to be stimulated over time by a stimulation function, fstim(t), which typically depends on endogenous growth factor(s) and/or exogenous drug. The dependence of fstim on these factors will be considered subsequently. The units for fstim are numbers of cells stimulated per unit time. Let it be assumed that these cells act fairly independently of each other, then one can apply the superposition principle. In the case where all cells have the same time invariant disposition parameters (a and b, Fig. 1, panel A) it follows that the number of cells present in the systemic circulation, N(t), is given by:
where * denotes the convolution operator. The integration limit of −∞ in Eq. (3) is to be interpreted to consider “all prior history” of the system. The above equation assumes a point distribution of cell residence times (i.e., zero variance). In order to readily fit Eq. (3) to observed data, it is convenient to transform it into a differential equation form. Equation (3) can be separated into two components as:
where t0 is the initial observation time. If we define:
Then Eq. (4) becomes:
If steady-state (SS) conditions are assumed to exist prior to t0 so that:
where is the steady-state stimulation rate, then:
Differentiating Eq. (9) with respect to time gives:
Next consider the case where the cellular disposition parameters are time variant. That is, the residence time in the systemic circulation, ‘b − a’, is a function of time. For simplicity we consider the time variance of the residence time to be due only to the ‘a’ parameter (Fig. 1, panel B). In this case the UR function becomes:
Here ‘s’ is the time that the cell is stimulated (Fig. 1, panel B). The value of ‘a’ is considered variable and determined at the time the cell is stimulated and remains fixed at that value thereafter. Again, a point distribution of cell residence times, at a given stimulation time, is assumed.
For this time-variant case, Eq. (4) becomes:
If similar to before we assume steady-state conditions for t ≤ t0, then:
where aSS is the steady-state ‘a’ value. Due to the steady-state assumption the first term, Ig(t), becomes nearly identical to Eq. (5), and similarly to before we get:
However, this time we write:
For convenience in derivation, we define:
Due to the unit step function, U(t), the integrand of Eq. (20) will be zero in the integration region u = t0 to u = t for ‘u’ values for which t − u − a (u) < 0. Thus, the switch between zero and non-zero integrand values will occur at u = x value(s) where ‘x’, that is a function of time, is the solution(s) to the equation:
Potentially, Eq. (22) may have multiple solutions for ‘x’ at the current time ‘t’. However, as evident from the graphical representation of the solution(s) to Eq. (22) given in Fig. 2, if the a(t) never has a slope, a′(t), less than or equal to −1 then Eq. (22) only has one solution (Fig. 2, time ‘t1’). For a more detailed proof, see Appendix B. Accordingly, Eq. (20) becomes:
It is recognized from Eq. (22) that ‘x’ is simply the time for the progenitor cell stimulation of those “new” cells that are currently entering the systemic circulation at time ‘t’.
Let FR denote the fraction of cells removed momentarily after the time of the phlebotomy, denoted ‘T’, then the equation to be integrated changes to a new “initial” condition at time ‘T’ according to:
where ε denotes an infinitesimally small time increment. However, just a change in a boundary condition according to Eq. (34) is insufficient to correct for the cells mechanically removed. The differential equation, Eq. (29), being integrated across the phlebotomy boundary point (t = T) must also be changed. If the original equation (Eq. (29)) is not changed then this equation will for some time period beginning immediately after the phlebotomy underestimate the number of cells because it contains an elimination component (i.e., ) that includes those cells mechanically removed at time ‘T’, which should not be included in the elimination. To correct for this, one needs to add a term that accounts for the elimination rate of those cells removed at time ‘T’ identical to their elimination had they not been removed by phlebotomy. Accordingly, the correction term (Appendix C) to be added is:
In deriving Eq. (35), it is recognized that the oldest cells in the removed cell population would be eliminated in vivo at t = T with a (projected) rate of FR · fstim (T − b), while the youngest cells would be eliminated at t = T + b − a(x(T)) with a rate of FR · fstim (T − a(x(T))). Thus, the correction term is active from the time period t = T to t = T + b − a(x(T)), consistent with the above considerations.
The above Eq. (36) describes the general case considering any arbitrary variation over time in the a(t) “parameter”. The following section is aimed at specifically modeling the a(t) variation over time.
The time-variant cellular disposition parameters are not readily determined by modeling cellular count data from the cells of interest alone. The observed response is a function of both the stimulation rate and the disposition “parameters”, which are both changing as functions of time, causing difficulty in the estimation procedure. One possible solution is to utilize shared disposition parameters or functions from “coupled” cell data (i.e., cells that transform into each other) to assist in the simultaneous determination of the stimulation rate and the disposition “parameters”.
To determine the time-variant reticulocyte residence time, the present analysis utilizes the shared input and disposition functions of reticulocytes and hemoglobin to overcome the above problem, since a reticulocyte has approximately the same amount of hemoglobin as a mature RBC . Accordingly, both the reticulocyte count and hemoglobin response share the same lag time, a(t), for appearance in the systemic circulation. The release of the hemoglobin into the systemic circulation is assumed to enter the circulation within a reticulocyte before it matures into an RBC. The ‘b’ parameter of the reticulocyte disposition was assumed to be time invariant. Furthermore, both the reticulocytes and hemoglobin share the stimulation function, fstim, which is multiplied by a conversion constant, KR, to change from units of hemoglobin to number of reticulocytes (i.e., has units of 103 cells/g of Hb). Additionally, since the same fractional removal in cells will occur for RBCs (i.e., hemoglobin) and reticulocytes, both share the FR parameter as well. Thus, reticulocytes and hemoglobin are “coupled” through several shared parameters and functions which generally assist in the analysis.
Measurement of hematopoietic cells in vivo is done in terms of concentrations, such as number of erythroid cells or grams of hemoglobin per volume of blood, therefore the number of cells (or grams of hemoglobin), ‘N’, is divided by the volume, ‘V’. However, in the analysis, due to the measurements units, the ‘V’ fuses with the parameters of the stimulation function, fstim, so that it is not estimated directly in the regression analysis. Erythropoietin was considered the stimulator for the stimulation of the progenitor cells. Accordingly, Ce is the concentration of Epo at the biophase (effect site). The analysis indicated that Ce is kinetically distinct from the Epo plasma concentration, Cp, and the following convolution relationship:
was used to model the plasma/biophase relationship resulting in the following differential equation:
where and denote the steady-state effect site and plasma Epo concentrations, respectively, and was set equal to the initial (first observation) fitted plasma Epo concentration. The stimulation of the erythroid progenitor cells (quantified in terms of hemoglobin) was related to Ce by the Hill equation:
where qstim is the redefined stimulation function that depends on Ce, Emax is the maximal hemoglobin stimulation rate in g/dl/day, EC50 is the effect site Epo concentration that results in 50% of maximal hemoglobin stimulation rate (Emax), and qstim is the redefined steady-state stimulation rate.
The analysis involved simultaneously fitting to both the hemoglobin and reticulocyte count concentration versus time profiles, denoted Hb(t) and R(t), respectively. Let ‘d’ denote the time the RBCs are removed from the systemic circulation relative to progenitor stimulation (i.e., the RBC equivalent ‘b’ parameter from Eq. (14) and Fig. 1) and ‘b’ denote the time the reticulocytes mature into RBCs relative to the stimulation. Both the ‘b’ and ‘d’ parameters are considered time invariant. As before, let a(t) denote the time that reticulocytes enter into the systemic circulation relative to progenitor stimulation. The time variance of a(t) is shared by the hemoglobin kinetics because the release of hemoglobin into the systemic circulation is assumed to initially enter within a reticulocyte. To ensure that b > a(t) ≥ 0 for all ‘t’ and for computational convenience, it is helpful to parameterize a(t) as:
Therefore, τ(t) is the (time-variant) reticulocyte residence time in the systemic circulation. The parameterization of Eq. (36) using Eq. (42) (Appendix D) with the appropriate substitutions (i.e., Eqs. (40) and (41)) from above results in the following specific key equations:
The final model is schematically detailed in Fig. 3.
All animal care and experimental procedures were approved by the University of Iowa Institutional Animal Care and Use Committee. Five healthy young adult sheep ranging in age from 2- to 4-months old and weighing 21.2 (3.60) kg (mean (SD)) were selected. Animals were housed in an indoor, light- and temperature-controlled environment, with ad lib access to feed and water. Prior to study initiation, jugular venous catheters were aseptically placed under pentobarbital anesthesia. Intravenous ampicillin (1 g) was administered daily for 3 days following catheter placement.
Blood samples (~0.5 ml/sample) for plasma Epo, reticulocyte counts, and hemoglobin were collected for 5–7 days to determine baseline values prior to conducting a single controlled phlebotomy over several hours to induce acute anemia. Animals were phlebotomized to hemoglobin concentration of approximately 4g/dl. To maintain a constant blood volume during the procedure, equal volumes 0.9% NaCl solution were infused for each volume of blood removed. Blood samples were subsequently collected at least daily for 25–30 days post-phlebotomy. No iron supplementation other than that in the animal’s feed was given. Plasma iron and total iron binding capacity (TIBC) were monitored to ensure animals did not become iron deficient. To minimize hemoglobin and red cell loss due to frequent blood sampling, blood was centrifuged, the plasma removed, and the red cells re-infused.
Plasma Epo concentrations were measured in triplicate using a double antibody radioimmunoassay (RIA) procedure as previously described (lower limit of quantitation 1mU/ml) . All samples from the same animal were measured in the same assay to reduce variability. The reticulocyte count was determined by flow cytometry (FACScan, Bection–Dickinson, San Jose, CA) as previously described . Hemoglobin concentrations were measured spectrophotometrically using an IL482 CO-oximeter (Instrumentation Laboratories, Watham, MA). In total, 30–45 samples were analyzed per variable per subject.
The proposed PD model is schematically detailed in Fig. 3. The pharmacokinetics of Epo (i.e., Epo plasma concentrations) were non-parametrically represented using a generalized cross-validated cubic spline function, with nodes at every data point . The hemoglobin and reticulocyte-plasma Epo concentration relationship was modeled by simultaneously fitting Eqs. (38)–(40) and Eqs. (43)–(50) to the modeled plasma Epo, observed blood hemoglobin, and observed reticulocyte count concentration–time data of each subject. The lag time parameter, ‘d’, between progenitor cell stimulation and RBC removal from the systemic circulation was fixed at ‘aSS + red cell (RBC + reticulocyte) lifespan’. The normal red cell lifespan has previously been determined in sheep using [14C] cyanate label and found to be 114 days . A fifth order constrained spline was used to semiparametrically estimate the residence time of the reticulocytes, τ(t), as a function of time. The residence time by the spline function was constrained to remain at the pre-phlebotomy steady-state value (τ0 = τSS) until some time (T1) after the phlebotomy (τ′(t) = 0 for t ≤ T1 < T2) and then return to a steady-state value (τend) at some later point (T2, τ′(t) = 0 for t ≥ T2), as detailed in Appendix E.
All modeling was conducted using WINFUNFIT, a Windows (Microsoft) version evolved from the general non-linear regression program FUNFIT , using weighted least squares. Data points were weighted by , where yobs is the observed hemoglobin or reticulocyte concentration. The analysis required the numerical solution to delay differential equations (DDE) because of the lag times involved in the mathematical model. This was done by WINFUNFIT using a DDE solver module that is based on a adaptation of the code “RETARD”  (see electronic supplementary material). The fraction of hemoglobin and reticulocytes removed, FR, was estimated by WINFUNFIT using an generalized events processing module that allows parameterized events to change during the iterative parameter estimation procedure. Accordingly, the magnitude and time for the events and delays can be estimated simultaneously with other unknown parameters of the mathematical model. Akaike’s Information Criterion (AIC) was used to compare the time variant model with 11 estimated parameters (four cellular production parameters (Emax, EC50, ke, and FR) and seven residence time parameters/coefficients (‘b’ and the spline parameters/coefficients)) to the identical time invariant model with only six estimated parameters (four cellular production and two residence time parameters (‘a’ and ‘b’)) .
To summarize the uncertainty in the individual subject parameter estimates, the mean percent standard error (MSE%) of the estimate was calculated for each parameter as:
where SEi and Pi are the standard error of the parameter and the estimate of the parameter for the ith subject, respectively, and ‘n’ is the number of subjects.
Mean reticulocyte residence times in the systemic circulation were compared using paired two tailed Student’s t-test at selected points in time using Microsoft® Excel 2002 SP3 (Microsoft Corporation, Redmond, WA). Statistically significant differences were determined at the α = 0.05 type I error rate.
The simultaneous fit to the Epo, reticulocyte, and hemoglobin data for two representative subjects is shown in Fig. 4. Good agreement between the observed and predicted values is evident from the fits (overall r = 0.978 (0.0084) (mean (SD)), hemoglobin r = 0.981 (0.0203), reticulocyte r = 0.966 (0.0114)). The PD model parameters are summarized in Table 1. The mean determined maximal rate of hemoglobin production, Emax, was estimated to be 1.29 (0.471) g/dl/day. The mean lag time between progenitor stimulation and reticulocyte disappearance from the systemic circulation, ‘b’, was estimated at 1.24 (0.167) days and the mean pre-phlebotomy steady-state reticulocyte residence time in the systemic circulation, τ0, was estimated to be 0.447 (0.100) days. Accordingly, the mean pre-phlebotomy steady-state lag time between progenitor stimulation and reticulocyte (hemoglobin) appearance in the systemic circulation, ‘aSS’, was estimated at 0.797 (0.193) days.
The reticulocyte systemic circulation residence time determined for each subject rapidly increased shortly after the phlebotomy, reaching a maximal value of 0.761–1.24 days at 5.95 (0.899) days post-phlebotomy (Fig. 5). Figure 5 displays the residence time for cells currently being stimulated, not the residence of reticulocytes currently entering the bloodstream or circulating in the systemic blood. The time post-phlebotomy that the residence time began to deviate from steady-state conditions was estimated at 0.730 (0.834) days. The maximal residence time was 2.01–2.64-fold higher than the initial steady-state residence time (P < 0.01). Following establishment of the maximal residence time, the residence time rapidly dropped in all subjects below its pre-phlebotomy value, to a minimum value of 0.471–0.906-fold of the initial steady-state residence time (P < 0.05). It then remained at that value or increased to a new steady-state value. On average, the residence time returned to steady-state conditions 15.4 (2.36) days following the phlebotomy. In all but one subject, the ending steady-state reticulocyte residence time, τend, was near the initial value (P > 0.25, including and excluding the single subject with the ending residence time well above the initial value). The reason for this different behavior in the single subject is not clear, but it is evident from the subject’s data (Fig. 4, panel B) that the reticulocyte count was substantially elevated compared to the baseline count well past the phlebotomy, supporting that the reticulocyte residence time likely remained elevated.
The MSE% for the primary parameter estimates are also displayed in Table I. The other residence time parameters (see Appendix E), T1 and T2, had MSE% of 1.2% and 1.4%, respectively. The two unconstrained coefficients of the residence time fifth order spline function, α4 and α5, had MSE% of 0.7% and greater than 1000%, respectively. The high MSE% for α5 was due to two subjects, as the other three subjects had a MSE% of 0.6%. Modeling the data instead with a fourth-order residence time spline function still resulted in a good model fit to the data for the two subjects with the high SE for the α5 coefficient, however, it resulted in an unacceptably poor fit to the data from the other three subjects. Therefore, a fifth-order spline function was used to model the time-variant residence time as the minimal model that resulted in an acceptable fit to the all the subject data.
In all subjects, the time variant reticulocyte residence model was superior to the identical time invariant reticulocyte residence model, where a single residence time in the systemic circulation was estimated, based on AIC. The difference in the fits between these two models for a representative subject is displayed in Fig. 6.
The determined steady-state lag time between effect-site Epo stimulation and reticulocyte or hemoglobin appearance in the systemic circulation, the ‘aSS’ value, of 0.797 (0.193) days is similar to the corresponding parameter we have determined in previous sheep phlebotomy studies of 0.72 days . The current model includes an effect site compartment, while the previous analysis did not. In contrast to the previous analysis , this analysis also includes hemoglobin data in the estimation procedure. The determined ‘aSS’ value (0.797 days or 19.1 h) is approximately the same as a similar physiological parameter determined in cynomolgus monkeys of 14.95 h and in humans of 10.76 h [19, 20].
The determined value and qstim function for each subject estimates a steady-state hemoglobin production of 0.0929 (0.0472) g/kg/day, assuming a blood volume of 0.75 dl per kg bodyweight. Furthermore, the maximum production rate of hemoglobin following the observed Epo stimulation (i.e., maximum Ce) was estimated at 0.504 (0.0422) g/kg/day, indicating that endogenous Epo released under maximally stimulated acute anemic conditions can increase hemoglobin (RBC) production approximately 5-fold. The determined 5-fold increase in hemoglobin production following acute anemia in sheep is similar to the 3- to 5-fold increase reported in human patients .
As would be predicted from the current understanding in the scientific community of reticulocyte maturation under stress erythropoiesis, the reticulocyte residence time in the systemic circulation increased post-phlebotomy (Figs. 5 and and6,6, panel B). The mean pre-phlebotomy steady-state reticulocyte residence time, τ0, was estimated at 0.447 (0.100) days (10.7 h) (Table 1). The estimated τ0 is similar to other estimates in sheep of 0.61 days and the commonly reported value of 1 day in humans [3, 23, 32]. From the steady-state value, the reticulocyte residence time increased 2.01- to 2.64-fold during stress erythropoiesis, before dropping below pre-phlebotomy values. The determined maximally increased reticulocyte residence time following a phlebotomy is approximately the same as that previously estimated in vivo in sheep and humans of 4.5-fold and 2–3-fold, respectively [23, 32]. It is also very close to the commonly recognized increase in reticulocyte residence time in the systemic in humans due to stress erythropoiesis of 2–3 days . The drop of the mean reticulocyte residence time below the pre-phlebotomy steady-state value is not generally noted in the literature, and may be due to an overshoot/rebound phenomenon often observed in biological systems.
It is important to note that the determined reticulocyte residence time is not the reticulocyte lifespan. The latter is not possible to determine accurately in vivo by only sampling in the systemic circulation. A reticulocyte initially normally exists in the bone marrow of sheep and humans [1, 3, 32], and therefore cannot be observed until it is released from the marrow into the systemic circulation. Accordingly, the determined drop below baseline reticulocyte residence time post-phlebotomy most likely represents the release of more mature reticulocytes from the bone marrow into the circulation through some rebound/overshoot phenomenon.
Utilizing a time-variant reticulocyte residence time model was superior to a time invariant reticulocyte residence time model for all subjects based on the AIC (Fig. 6). As can also be seen in Fig. 6, an adequate simultaneous fit to both the reticulocyte and hemoglobin data is not possible with a time invariant reticulocyte residence time model. The current model is also an improvement over previous Epo, reticulocyte, and hemoglobin integrated PK/PD models published by our laboratory in accounting for hemoglobin production , since previous models did not include the hemoglobin in reticulocytes as part of the total measured hemoglobin. The hemoglobin contained in the reticulocyte, as previously modeled, only added to the measured hemoglobin upon maturation of the reticulocyte into a RBC. Furthermore, previous models used the percent reticulocytes instead of the absolute reticulocyte count, the former which is not as easily interpreted physiologically. The current analysis also provides an exact correction for the effect of removing erythroid cells during the phlebotomy process.
Because the current model shared functions through the improved hemoglobin accounting, a time-variant reticulocyte residence time could be modeled, which is widely known to occur physiologically, but apparently has not been previously modeled. It may be physiologically relevant to include time variance in the ‘b’ parameter in addition to the ‘a’ parameter. However, the time-variant residence could not be modeled as such because a variable ‘b’ value in addition to a variable ‘a’ value would not be possible to determine in a reliable way. Since reticulocytes are released from the bone marrow at a younger age following stress erythropoiesis, modeling the time variance due to changes in the ‘a’ value accounts for the primary physiological features of reticulocytes produced under stress erythropoiesis. The ability to estimate a time-variant ‘a’ value is primarily due to the simultaneous fitting to the “coupled” hemoglobin and reticulocyte responses. Two key features of the “coupled” responses are the longer residence time in the systemic circulation of hemoglobin compared to reticulocytes and the larger pool of hemoglobin relative to reticulocytes. Accordingly, a substantial change in the ‘a’ value, and thus the hemoglobin and reticulocyte input, has a relatively large and rapid effect on the observed reticulocyte response, but only a relatively minor effect on the observed hemoglobin response. Alternatively, a model with a time invariant ‘a’ parameter but with a time variant ‘b’ parameter could be formulated and fitted to the data to estimate the reticulocyte residence time in the systemic circulation. However, a model formulated as such would not acknowledge the primary features of erythropoiesis physiology. The model also determined a time-variant hemoglobin residence time in the systemic circulation. However, we only report and discuss the results of the model in terms of the change in the reticulocyte residence time, τ, and the observed reticulocyte response due to the small magnitude of the percent change in the hemoglobin residence time (<1%), the fixing of the ‘d’ parameter, and the limited observation time relative to the residence time of hemoglobin.
Another important model assumption is the fixed ‘d’ value of ‘aSS + red cell lifespan’ for all cells, since it has been experimentally demonstrated in mice and rats that red cells produced under stress erythropoietic conditions have a shortened lifespan or residence time in the systemic circulation [33, 34]. However, the lifespan of the red cells would need to be reduced by greater than 70% (i.e., to less than 30 days) to have an effect, since this was the longest period of time of observation post-phlebotomy. A final simplifying assumption is that a single, calculated hemoglobin to reticulocyte count conversion constant, KR, was used for all subjects. The KR parameter acts only as a scalar though, so any changes in the magnitude of the reticulocyte residence time would still be observed even if the value for KR varied from subject-to-subject. Additionally, a constant KR would be inappropriate if the amount of hemoglobin per erythroid cell varied with time, as may occur if subjects became iron deficient. However, monitoring of iron status did not indicate that any subjects became iron deficient throughout the time course of the experiments.
The model, as currently described, does not account for a distribution(s) of cell residences, since it was assumed that all cells stimulated at a given time have the same “pre-determined” residence time. However, previous work has demonstrated little advantage in moving to a more complicated residence or lifespan distribution model compared to a point distribution model due to problems estimating the standard deviation of the cell lifespans . The simultaneous determination of a cellular distribution(s) and changes in the individual ‘a’ and/or ‘b’ values would be mutually confounding and would create tremendous estimation problems.
The present analysis did not indicate that a′(t) at any time fell at or below a value of −1 (nor did τ′(t) fall at or above a corresponding value of +1). Thus, the problem associated with accounting for multiple solutions, possible according to Eq. (22), fortunately did not arise in this analysis (see Fig. 2, showing multiple solutions when a′(t) becomes < −1). Physiologically, what occurs when multiple solutions to Eq. (22) arise is that cells that were stimulated at a later point in time are entering the sampling compartment prior to cells that were stimulated earlier (i.e., “younger” cells are entering prior to “older” cells), resulting in multiple groups (solutions) of cells that were stimulated at different times entering the systemic circulation at exactly the same time. Additionally, mathematically Eq. (30) becomes undefined when a′(t) = −1, because there is an infinite number of solutions along the finite range where a′(t) continuously equals −1. If the analysis had indicated that a′(t) at any time fell at or below −1, we would have had to resort to much more complex numerical solutions to the already complex delay-type differential equations describing the proposed model.
The time-variant cellular disposition model presented in the theoretical section is applicable to modeling any cellular response where changes in cellular residences is expected and may be extended to incorporate unique features of the specific physiological system being modeled. Furthermore, the specific Epo, reticulocyte, and hemoglobin PK/PD model that was presented could readily be applied to humans, both adult and neonate, and other species. Future modeling work and more “complete” data may also allow for the relaxation of certain model assumptions (i.e., a point distribution of cell residences), developing an even more general PD model of cellular responses.
In summary, by sharing parameters and functions between the measured hemoglobin and reticulocyte count response, the reticulocyte residence time in the systemic circulation was accurately determined as it changed with time. The utilized methodology allows for the dynamic determination of time-variant cellular disposition or residence time without directly measuring the cellular residence times. The determined increase in reticulocyte residence time from reticulocytes produced under stress erythropoiesis and the determined increase in hemoglobin production is consistent with other experimental data and with the scientific community’s current understanding of erythropoiesis physiology. Furthermore, the determined reduction below baseline residence time post-phlebotomy may serve to generate future experiments for determining why or how physiologically the reduction in reticulocyte residence times occurs in vivo in sheep, and ultimately in humans.
The authors thank Dr. Robert T. Cook and the personnel of the Iowa City VAMC Pathology and Laboratory Medicine Services for assistance in performing the flow cytometric measurements of the reticulocytes. This work is supported by the United States Public Health Service, National Institute of Health grants PO1 HL46925, and R21 GM57367, and by the Veterans Administration Medical Center, Iowa City, Iowa.
Derivation of Eq. (11)
By evaluation of Eq. (6):
Proof of unique solution to Eq. 22 when a′(t) > −1.
Equation 22 may be written as:
with y x(t). Treating ‘t’ as fixed the Jacobian determinant, |JG|, of Eq. (B1) is:
Thus by the Implicit Function Theorem , Eq. 22 has a unique solution if |JG| ≠ 0, i.e., if 1 + a′(y) ≠ 0, which is satisfied when a′(y) > −1 or when a′(y) < −1 for all ‘y’. The condition a′(t) > −1 for all ‘t’ values ensures that a′(y) > −1 for all ‘y’ values, thus ensuring a unique solution to Eq. 22.
Derivation of Eq. (35)
It is recognized that the phlebotomy at time t = T can only remove cells already present in the circulation. If we separate Eq. (15) such as:
It is apparent that the Ij(t) term gives all the cells that have already entered into the circulation at time, t = T, and thus all the cells that can be removed by the phlebotomy.
By evaluation of Eq. (C2):
Since the Ij(t) term represents the cells that can be removed by the phlebotomy, it is multiplied by (1 − FR) of Eq. (34) at time t = T + ε, resulting in:
The first term on the right hand side of Eq. (C7), Ij(t) + Ik(t), is recognized as N(t) of Eq. (15), which if differentiated will result in Eq. (29). Therefore, −FR · Ij(t) represents the correction term. Since t − T > 0 in Eq. (C7), Ij(t) only has non-zero values from time t = T to t = T + b − a (x (T)) as given by Eq. (C5). Differentiation of Eq. (C5) results in:
Thus, the correction term generalizes to:
which proves Eq. (35).
Alternative parameterization of Eq. (36)
The exponential transformation ensures that a (x (t)) > 0 for any x(t) value. Rearrangement and differentiation give:
Residence time spline function
The reticulocyte systemic circulation residence time, τ(t), was modeled as follows:
In order to satisfy the above constraints, the fifth-order polynomial that represents τ(t) when t1 ≤ t ≤ t2 was parameterized in terms of the parametric terms, τ0, τend, T1 and T2, and the non-parametric coefficients, α4 and α5.
Kevin J. Freise, Division of Pharmaceutics, College of Pharmacy, The University of Iowa, Iowa City, IA 52242, USA.
John A. Widness, Department of Pediatrics, College of Medicine, The University of Iowa, Iowa City, IA 52242, USA.
Robert L. Schmidt, Department of Pediatrics, College of Medicine, The University of Iowa, Iowa City, IA 52242, USA.
Peter Veng-Pedersen, Division of Pharmaceutics, College of Pharmacy, The University of Iowa, Iowa City, IA 52242, USA.