PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
Bull Math Biol. Author manuscript; available in PMC 2010 October 1.
Published in final edited form as:
PMCID: PMC2834560
NIHMSID: NIHMS146225

Interpreting CFSE Obtained Division Histories of B Cells In Vitro with Smith-Martin and Cyton Type Models

Abstract

The fluorescent dye carboxyfluorescin diacetate succinimidyl ester (CFSE) classifies proliferating cell populations into groups according to the number of divisions each cell has undergone (i.e., its division class). The pulse labeling of cells with radioactive thymidine provides a means to determine the distribution of times of entry into the first cell division. We derive in analytic form the number of cells in each division class as a function of time using the cyton approach that utilizes independent stochastic distributions for the time to divide and the time to die. We confirm that our analytic form for the number of cells in each division class is consistent with the numerical solution of a set of delay differential equations representing the generalized Smith-Martin model with cell death rates depending on the division class. Choosing the distribution of time to the first division to fit thymidine labeling data for B cells stimulated in vitro with lipopolysaccharide (LPS) and either with or without interleukin-4 (IL-4), we fit CFSE data to determine the dependence of B cell kinetic parameters on the presence of IL-4. We find when IL-4 is present, a greater proportion of cells are recruited into division with a longer average time to first division. The most profound effect of the presence of IL-4 was decreased death rates for smaller division classes, which supports a role of IL-4 in the protection of B cells from apoptosis.

I. INTRODUCTION

The fluorescent dye carboxyfluorescin diacetate succinimidyl ester (CFSE) has been used to measure the number of times a cell has divided since it was labeled [1]. CFSE is a fluorescent dye that binds to intracellular proteins and is distributed equally between each daughter cell at mitosis, resulting in a two-fold dilution. The precise partitioning of this dye allows the clear resolution of 7-8 sequential division cycles, both in vitro and in vivo. Hence at a specific time point, one can measure how many cells are in each division class, where a cell in the nth division class has divided n times since labeling. CFSE labeling has been used to measure the correlation between cell division and surface molecule expression or internal expression of cytokines [2-4], division-linked immunoglobulin class switching in B cells [5-7], and cytokine regulation of T cell kinetics [8, 9], and B lymphocyte proliferation and death rates in vivo [10].

Several approaches have used CFSE data to estimate cell kinetic paramaters such as division rate, death rate, and the proportion of the cell population that has divided (precursor frequency). Several groups have used systems of linear birth-death type ordinary differential equations (ODEs) to model the dynamics of the number of cells in each division class [10, 11]. One problem of such ODE models is that under some circumstances they predict the unrealistically rapid appearance of cells in high division classes. To address this issue, features of the cell cycle were incorporated into these models, most commonly using the Smith-Martin cell cycle model [12-18]. Using this approach, the kinetic parameters of murine bone marrow cells were estimated by dividing the cell cycle into a variable length resting or A phase (or A state), and a fixed length division or B phase [14]. Other adaptations to improve the biological fidelity of ODE models have included considering the time to completion of the first mitosis as a separate parameter within the context of the “generalized Smith-Martin” model [9, 19-21]. This has then been experimentally measured using tritiated thymidine measurement of the time of entry into the first division. Since thymidine is incorporated into the DNA of dividing cells during the S phase, the level of radioactivity is proportional to the total number of cells in the S phase at a specific time point [22]. The beginning of the S phase locates close to the boundary between the A and B phases of the cell cycle. If in the experimental protocol a mitosis inhibitor is used, cells are then able to replicate their DNA, but can not undergo cell division [23]. As a consequence, DNA replication is only measured for cells entering their first division.

Recently introduced cyton model assumes that cell division and cell death are governed by independent stochastic processes [24]. Given a starting cell number and then representing a cyton by using two probability distributions, one for the time to divide and one for the time to die, the number of cells dividing or dying per unit time at time t was calculated [24]. Then the total cell number in each division class as a function of time t was evaluated by the difference between the number of cells that have entered that division class until time t and the number that have left by means of division or death [24]. By assuming lognormal distributions for the time to death and for the time to divide, numerical evaluations of the cyton model successfully described CFSE-labeled B cell proliferation and death data [24].

Here we report a closed analytical solution of the cyton model for the case of the generalized Smith-Martin model with an inhomogeneous death rate that is a function of the division class. The generalized Smith-Martin model uses a time-delayed gamma distribution to model the time to enter the first division and to fit the measured distribution of times to first division obtained from thymidine incorporation experiments. In the subsequent division classes, the standard assumptions of the Smith-Martin model [12] are made: a constant rate of activation from the A phase to the B phase of the cell cycle and a fixed length of the B phase. These assumptions enable us to write the probability distribution of time to cell division as a delayed exponential distribution. The probability distribution for the time to die is assumed to be an exponential distribution with a different death rate in each division class. Using these assumptions, we were able to solve the cyton model analytically. We confirmed that the analytical solutions based on the cyton approach were consistent with numerical solutions of a set of delay differential equations (DDEs) describing the generalized Smith-Martin model.

Finally, we tested the usefulness of the generalized Smith-Martin model with inhomogeneous death rates to analyze thymidine incorporation and CFSE data obtained from B cell cultures stimulated with lipopolysaccharide (LPS) in the presence or absence of the immunomodulatory cytokine interleukin-4 (IL-4). IL-4 is a cytokine produced by multiple cell types that regulates survival and proliferation of B cells. Numerous studies reported that IL-4 enhances the survival of B cells through up-regulation of antiapoptotic proteins and through regulating glucose energy metabolism [23, 27-29]. Enhancement of survival of B cells in in vitro culture by adding IL-4 was directly observed [24]. The cytokine IL-4 modulates human B cell expansion after B cell receptor crosslinking, and thus would be expected to alter either time to first cell division after B cell receptor crosslinking or the cell cycle time in subsequent divisions [25, 26]. In addition, IL-4 fosters isotype switching [30, 31]. Using the developed cyton model, we sought to discern how the presence or absence of IL-4 after B cell activation influenced cell division and cell death of the B cells.

We observe that IL-4 has more effect on the death rate of B cells than on the rate of proliferation. The protection of B cells from apoptosis was significant by the addition of IL-4 especially in the lower division classes. While IL-4 also affects the proportion of cells recruited into the first division, the average cell cycle time after the first division is hardly affected by IL-4.

II. EXPERIMENTAL METHODS

A. B cell isolation and cell culture

Single-cell suspensions were prepared from spleens and lymph nodes of mice, and red cells were lysed and run on discontinuous Percoll density gradients as described [32]. Small dense cells were harvested from the 65%/80% interface, and B cells were purified by means of negative selection by using magnetic beads (Miltenyi Biotec, Bergisch Gladbach, Germany). Cells were labeled with CFSE (Molecular Probes, Eugene, OR) according to the originally published method [33]. B cells were > 95% B220+, CD19+, IgM+, IgD+ as determined by flow cytometry. B cells were cultured in B cell medium consisting of RPMI 1640 medium (Gibco BRL, Grand Island, NY), 2 mM L-glutamine, 0.1 mM nonessential amino acids, 10 mM HEPES, pH 7.4, 100 g/ml of streptomycin, 100 U/ml of penicillin, 50 M beta-mercaptoethanol (all supplements were purchased from Sigma, St. Louis, MO) and 10% heat-inactivated FCS (CSL, Parkville, Australia).

B. Colcemid analysis

For direct analysis of time of entry into first division, cells were stimulated in the presence of 25 ng/ml colcemid (Roche, Indianapolis, IN). Cultures were pulsed with 1 Ci/well [methyl-3H] thymidine (Amersham Biosciences, UK) for 1 hour before harvesting. Incorporation of radioactivity was measured using a scintillation counter (Topcount NXT, Packard, Meridan, CT).

III. MODEL DEVELOPMENT

A. Delay differential equations of Smith-Martin model with inhomogeneous death rates

Assume that there are N0 cells at t = 0 and that they are recruited into division at rate, R(t). We call R(t) the recruitment function. Then according to the Smith-Martin model [12], the number of cells at time t in the A phase, equation M1 and the B phase, equation M2, of division class 1 satisfy the following differential equations [20, 21],

equation M3
(1)

equation M4
(2)

where λ is a constant activation rate from phase A to phase B, d1 is the death rate of the cells in division class 1, and Δ is the amount of time a cell spends in the B phase, i.e., the length of the B phase. For division classes, n ≥ 2, the dynamics of the number of cells in each phase is given by

equation M5
(3)

equation M6
(4)

where dn is the death rate of cells in division class n. The first term on the right hand side of Eq. (3) is the rate of entrance of cells into the A phase of the nth division class. Each cell that exits the B phase of the (n − 1)st division class contributes two progeny. The number of cells that exit this B phase per unit time is the number that enter per unit time, equation M7, and survive for time Δ, edn−1Δ, assuming a constant death rate dn−1. The second term is the rate at which cells exit the A state either by activation or death. In Eq. (4), cells that leave the A state enter the B phase, die at constant rate dn, and assuming they survive for a time Δ then leave the B phase.

B. Time to first division distribution and the recruitment function

The total population of B cells can be classified into two groups: one group of cells eventually divides and the cells in the second group die at rate dnd without ever completing a division. The proportion of cells completing a first division and entering the A phase of division class 1 is called the precursor frequency, ϕ [9, 13, 20]. In the Smith-Martin model, the A phase is generally considered as corresponding to G1 or G0 and G1, and the B-phase as corresponding to the S, G2, and M phases of the cell cycle. While the length of B phase is assumed to be fixed, the duration of A phase may differ considerably between cells. Since a cell is not committed to division until it enters the S phase of the cell cycle, following Lee and Perelson [20], we choose a probability distribution for the time to enter the first division, TFD(τ), as the delayed gamma distribution,

equation M8
(5)

and TFD(τ) = 0 for τ < τd. Note that this distribution in the appropriate parameter range appears similar to the log-normal distribution used by Deenick et al. [9] used to fit T cell thymidine incorporation data. Hawkins et al. [24] used a similar construct in the cyton model to describe the time to B cell division. Assuming a constant B phase length in division class 0, Δ, the recruitment function R(t) in Eq. (1) is

equation M9
(6)

and R(τ) = 0 for τ < τd + Δ. As we show below, using Eqs. (5) and (6) allows us to fit the measured thymidine incorporation data and to derive analytical solutions to the cell kinetics problem.

C. Cyton model

Rather than directly solving the differential equations (1)-(4) we use the cyton approach introduced by Hawkins et al. [24]. The cyton model assumes independent distributions of times to divide and times to die in each division class. In each division class, cell death is described with a probability distribution for the time to death, Kn(τ), and time to cell division is described with a probability distribution, ψn(τ). Thus, the probability of a cell surviving until time τ since entering division class n, Ln(τ), is given by

equation M10
(7)

The probability that a division does not occur until time τ since entering the division class n, is given by

equation M11
(8)

The probability of cells entering nth division class at time τ, Sn(t), is the probability of entering the (n − 1)th division class at time t′, Sn−1(t′), multiplied by the probability of survival in division class n − 1 for the duration tt′, Ln−1(tt′), multiplied by the probability of dividing at time tt′, ψn−1(tt′), yielding

equation M12
(9)

Likewise the probability of cells dying in the nth division class at time t, Dn(t), is given as the probability of entering the nth division class at time t′, Sn(t′), multiplied by the probability that the division to the next division class does not occur for the duration tt′, ηn(tt′), multiplied by the probability of dying after being in the nth division class for time tt′, Kn(tt′), yielding

equation M13
(10)

The number of cells in each division class is given as the difference between the total number of cells that entered the division class and the number of cells that are lost either by dying or exiting the division class by further division,

equation M14
(11)

where N0 is the initial number of cells at time 0 and the factor 2n accounts for the clonal expansion of cells that have divided n times starting from a single cell [15]. Since Sn(t) is the probability of cells entering the nth division class at time t, which depends on all the probability distributions in previous division classes as in Eq. (9), Eq. (11) holds for any arbitrary distributions of cell death and cell division depending on the division class.

D. Cyton approaches for differential equations of the Smith-Martin model with inhomogeneous death rates depending on the division class

We next formulate Eq. (9) for the Smith-Martin model with the delayed gamma distribution recruitment function in Eq. (6). The death rate of the precursor ϕ N0 cells in division class 0 is zero since these cells all eventually divide. Hence we have K0(τ) = 0 and the distribution of the time to the first division,

equation M15
(12)

and ψ0(τ) = 0 for τ < τd + Δ. Then from Eqs. (7) and (8), we have L0(τ) = 1 and η0(τ) = ϕ(1 + (ττd − Δ)/s)e−(ττd−Δ)/s. For the subsequent division classes, we assume a constant death rate which is a function of division class, dn. Thus, the distribution for the time to death in the division class n is Kn(τ) = dne−dnτ and the distribution for surviving the duration of τ since entering the division class n is Ln(τ) = ednτ from Eq. (7). Assuming a fixed length of the B phase, Δ, and a constant activation rate from the A to the B phase, λ, the probability of dividing (i.e., exiting the B phase) is ψn(τ) = λe−λ(τ−Δ) for τ > Δ and ψn(τ) = 0 for τ < Δ. From Eq. (8), the probability of not dividing for the duration of τ since entering the division class n is ηn(τ) = e−λ(τ−Δ) for τ > Δ and 1 for τ < Δ. We illustrate the distributions for the time to division and the time to death in Fig. 1.

FIGURE 1
Cyton representations for the generalized Smith-Martin model. (A) Distribution of time to first division, ψ0(τ) in Eq. (12) with τd = 10 h, s = 10 h, Δ = 10 h, and ϕ = 1 and distribution of time to divide in subsequent ...

The probability of entering the 0th division class is S0(t) = δ(t − 0), i.e., all cells enter the division class 0 at time 0. The probability of entering division class 1 is calculated as

equation M16
(13)

where t0 = tτd − Δ and H(t) is the unit step function, i.e., H(t) = 1 when t ≥ 0 and H(t) = 0 when t < 0. The probability of entering division class 2 is calculated as

equation M17
(14)

with m1 = 1/s − λ − d1 and t1 = tτd − 2Δ. We generalize this and find

equation M18
(15)

with mi = 1/s − λ − di and tn−1 = tτdnΔ.

The probability of cells dying in the 0th division class at time t, D0(t), is calculated as

equation M19
(16)

The probability of cells dying in the first division class at time t, D1(t), is given as

equation M20
(17)

where m1 = 1/s − λ − d1, r1 = 1/sd1, and t1 = tτd − 2Δ. For n ≥ 2, we have

equation M21
(18)

where mi = 1/s − λ − di, rn = 1/sdn, and tn = tτd − (n + 1)Δ.

The integral of Sn(t) and Dn(t) that are needed to calculate Nn(t) in Eq. (11) are straightforward. The results are

equation M22
(19)

equation M23
(20)

Non-dividing cells are assumed to have a different death rate, dnd, than the cells that are capable of dividing. With these generalizations, the number of cells in each division class is given as

equation M24
(21)

IV. RESULTS

A. Direct comparison between numerical solutions of the delay differential equation model and the cyton model

The cyton model assumes independent control of cell survival and cell division [24]. This assumption would be valid if the cellular machinery governing the time to divide is different and independent from the cell machinery controlling the time to die. When the distribution of the time to first division is assumed to be a delayed gamma distribution and the distribution of the time to subsequent division is assumed to be given by a Smith-Martin model with an exponential distribution of times in the A phase with a constant length of the B phase (Fig. 1), we calculated from Eq. (21), the number of cells in each division class as a function of time. In each division class the time to death was also assumed to be governed by an exponential distribution (Fig. 1). To confirm the validity of Eq. (21), we compared the cyton solutions in Eq. (21) with the numerical solutions of the DDEs of Eqs. (1)-(4). Figure 2 shows that the number of cells in each division class calculated by the cyton model is exactly consistent with that calculated by numerical integration of the DDEs.

FIGURE 2
Cell numbers as a function of time for each division class from cyton solutions in Eq. (21) and the numerical solutions of delay differential equations for the generalized Smith-Martin model in Eqs. (1)-(4). Here the rate of death in division class 1 ...

B. Thymidine incorporation data and the time to first division distribution

The level of incorporation of tritiated thymidine into DNA was measured at various times after a pulse of tritiated thymidine was added to the B cell LPS culture with and without IL-4 stimulation. The level of radioactivity is assumed proportional to the number of cells in the S phase since thymidine is incorporated into the DNA of dividing cells during the S phase. Assuming the number of cells entering the first S phase is the number of cells entering the first B phase of the Smith-Martin model, the distribution of time to enter the first division, TFD(τ) in Eq. (5), describes the proportion of cells labeled with thymidine as a function of time in this experiment.

Figure 3 displays the best-fit of Eq. (5) to the thymidine incorporation data. The distribution of time to enter the first division shows a stretched tail, increasing the mean (τmean = τd + 2s) and the standard deviation equation M25 of the delayed gamma distribution of Eq. (5) with the addition of IL-4 (1000 Unit/ml) [Fig. 3 (E) and (F)]. Greater mean and standard deviations of the distribution of time to enter the first division was also observed for T cells stimulated in vitro under greater level of IL-2 concentrations [20].

FIGURE 3
B cells were cultured in the presence of colcemid and a 1h thymidine pulse was given. Cell were harvested at the times indicated in the figure and scintillation counted from cultures without (A) and with (B) IL-4. The fitted time to enter the first division, ...

The best fit of the lognormal distribution, equation M26, to the thymidine incorporation data is presented as a dotted line in Fig. 3. The lognormal distribution provides a better fit than the gamma distribution. The sum of squared errors using the lognormal distribution (gamma distribution) are 5.02×10−6 (1.43×10−5) and 4.49×10−6 (1.19×10−5) without IL-4 and with IL-4, respectively. Although the lognormal distribution provides a slightly better fit to the data, the use of the gamma distribution allowed us to analytically solve the model.

C. CFSE data and the generalized Smith-Martin model with inhomogeneous death rate

Having fitted the experimental thymidine incorporation data with T F D(τ) in Eq. (5), we next attempted to estimate the cell kinetic parameters from CFSE data obtained from B cells under the same experimental conditions of the thymidine labeling using the cyton solutions of the generalized Smith-Martin model with inhomogeneous death rates. A division-dependent death rate in the generalized Smith-Marin model has previously been introduced to study CD4 T cell proliferation and death kinetics by Ganusov et al. [21]. For the CFSE data without IL-4 stimulation, we assumed that the death rate increases linearly as a function of the division class, n,

equation M27
(22)

In comparison with the model with a constant death rate, using a linear increase of the death rate, Eq. (22), gave a decrease in the sum of squared errors, which is not statistically significant (F-test, F = 3.05, numerator degrees of freedom=1, denominator degrees of freedom=50, p = 0.09). However, when we assumed that death rate increases linearly in two steps,

equation M28
(23)

we were able to improve the fit of the model to the CFSE data (F-test, F = 27.4, numerator degrees of freedom=2, denominator degrees of freedom=49, p=1.02×10−8). Here the sum of squared errors divided by the number of degree of freedom is 4.5 ×106 as compared with the value of 9.2 ×106 for a constant death rate. When the slope of the death rate was changed at n = nx = 2, the improvement of the fit using the two-step increase of the death rate was the most significant.

For the CFSE data obtained from B cells stimulated with LPS plus IL-4, assuming a linear increase of the death rate, Eq. (22), resulted in a decrease in the sum of squared errors as compared with constant death rate. However, this improvement was not statistically significant (F-test, F = 0.74, numerator degree of freedom=1, denominator degree of freedom=50, p = 0.39). When we assumed that death rate increases linearly in two steps as in Eq. (23), the fit of the model to CFSE data was improved (F-test, F = 183, numerator degrees of freedom=1, denominator degrees of freedom=49, p = 3.7 × 10−18). Here the sum of squared errors divided by the number of degree of freedom is reduced from 7.51 ×107 to 1.62 ×107 by assuming a division dependent death rate rather than a constant death rate. When the slope of the death rate was changed at n = nx = 4, the improvement of the fit by the two-step increase of the death rate was the most significant.

Figure 4 compares the number of cells in each division class from the CFSE data with the fit of the model for Nn(t) in Eq. (21) with the bi-linear increase of death rates in Eq. (23). We plot the comparison in two different ways: the dynamics of cell cohorts vs. time as proposed in Ref. [19] [Fig. 4 (a) and (b)] and the number of cells as a function of division class [Fig. 4 (c) and (d)]. To show how the fit is improved by the bi-linear increase of death rates in comparison with a constant death rate, we also plot the result of the fit to the CFSE data with constant death rates in Fig. 5. The bi-linear increase of death rates improves the fit visually as well as statistically.

FIUGRE 4
Fit of the cyton solution (solid line), given by Eq. (21) with a bi-linear increase of the death rate with division class, Eq. (23), to CFSE data (symbols) for B cells stimulated with LPS without (a)/(c) and with (b)/(d) IL-4. In (a) and (b), the data ...
FIGURE 5
Fit of the cyton solution (solid line), given by Eq. (21) with constant death rates after the first division to CFSE data (symbols) for B cells stimulated with LPS without (a)/(c) and with (b)/(d) IL-4. In (a) and (b), the data (each symbol) depicts the ...

Using the parameters determined from fitting to the distribution of the time to first division, τd, s, and then fitting the CFSE data allows us to estimate the average cell cycle time after the first division, CCT = 1/λ + Δ, the precursor frequency, ϕ, the profile of the death rate of cells having completed one division, dn in Eq. (23), and the death rate of non-dividing cells, dnd, with and without IL-4 stimulations, respectively.

In the Table I, we summarize the fitted parameters with 95% confidence intervals calculated by bootstrapping the residuals with 103 simulations [35]. The average cell cycle time of the 0th division class (CCT0 = τd + 2s + Δ) is estimated to be lower in the absence of IL-4 stimulation. The addition of IL-4 stimulation increases the average cell cycle time at the subsequence division classes (CCTn) from 6.1 h to 7.7 h. We also estimated that more cells are recruited into division with the addition of IL-4: the precursor frequency, ϕ is increased around 2.4-fold. We expect that the area under the curve representing the level of incorporated thymidine as a function of time would be proportional to the precursor frequency. The area under the thymidine incorporation level was observed to increase with the addition of IL-4 around 2-fold [Fig. 3 A and B], which is comparable to the estimated ratio of the precursor frequencies with and without IL-4 from our model.

TABLE I
Estimated parameters with 95% confidence intervals. CCT0 is average cell cycle time of division class 0. CCTn is average cell cycle time at the subsequence division classes. s and τd are the parameters in TFD(τ). Δ is the length ...

We estimated that d, the death rate of B cells averaged over division classes after the first division, increases with IL-4 stimulation (Table I). The death rate profile with division class is plotted in Fig. 8. In the presence of IL-4, the death rates are estimated be lower than those of B cells in culture without IL-4 when the division class is lower than 4. At higher division classes, we estimate greater death rates in the culture with IL-4. The death rate of non-dividing cells is estimated to decrease with the addition of IL-4. In the division classes less than 4, our estimate of a smaller death rate of cells with IL-4 is consistent with the observation that IL-4 is an anti-apoptotic cytokine for B cells [23, 27-29]. We plot the total number of cells as a function of time with and without IL-4 in the inset of Fig. 8. The number of cells at 87 hours in the presence of IL-4 is approximately double the numbers in culture without IL-4. On average, five divisions are completed at 87 hours with IL-4 stimulation. The greater number of cells with the IL-4 stimulation might be responsible for the enhanced cell death at higher division classes due to more rapid consumption of cell survival factors.

FIGURE 8
The predicted profile of the death rate as a function of division class. All the parameters for the fitted profile of the death rate in Eq. (23) are listed in the table I. The addition of IL-4 (1000 Unit/ml) results in smaller death rates when division ...

As well as reduced death rates, we found a greater fraction of cells ultimately divide with the addition of IL-4 stimulation, which also contributes to the greater total number of cells with IL-4 stimulation (Fig. 8 inset). Figure 4 shows that there is a higher proportion of cells in higher division classes in the absence of IL-4 stimulation. We can interpret this observation in light of our estimate of smaller death rates in division classes less than 4 and larger death rates in higher division classes greater than 5 with the addition of IL-4 stimulation. Addition of IL-4 stimulation recruits more cells to division, however, the increase of the death rates in the higher division classes leads to a reduced proportion of cells in higher division classes in comparison to the situation without IL-4 stimulation. A slightly lengthened average cell cycle time also contributes to a reduced proportion of cells in higher division classes in the presence of IL-4.

When the cyton solution, given by Eq. (21), with a bi-linear increase of death rate with division class was fit to the CFSE data, the estimated death rates in division class 1 were small, 0.01 h−1 without IL-4 and 6.58 × 10−4 h−1 with IL-4 (Table I). This suggested that d1 might be zero. To test this, we redid the fit fixing d1 = 0. The results of the fits and the estimated parameters with 95% CIs are shown in Figure 6. The estimated parameters when we fixed d1 as 0 are quantatitively similar to those when we allowed d1 to be a free parameter. Furthermore, even though we reduced the number of parameters by one, the sum of squared errors of the fit to the CFSE data decreases to 7.69 × 108 when we fix d1 = 0 from 7.98 × 108 in the presence of IL-4. However, in the absence of IL-4, having d1 as a free parameter improves the fit in a statistically significant way (F-test, F=23.8, numerator degress of freedom=1, denominator degress of freedom=49, p = 1.16 × 10−5).

FIGURE 6
Fit of the cyton solution (solid line), given by Eq. (21) with a bi-linear incerase of death rate with division class by fixing d1 = 0 to CFSE data (symbols) for B cells stimulated with LPS without (a)/(c) and with (b)/(d) IL-4. In (a) and (b), the data ...

D. CFSE data and the generalized Smith-Martin model with lengthening of the cell cycle time

As an alternative for inhomogeneous death rates, we examined the possiblity of describing the CFSE data with the generalized Smith-Martin model with the length of the B phase depending on the division class. When the length of the B phase changes, Eqs. (1)-(4) are modified as,

equation M29
(24)

equation M30
(25)

where Π1 is the length of B phase in the division class 1

equation M31
(26)

equation M32
(27)

where Πn is the length of B phase in the division class n and d is the constant death rate of cells. The recruitment function R(t) in Eq. (6) is modified as

equation M33
(28)

and R(τ) = 0 for τ < τd + Δ0. Here Δ0 denotes the length of B phase in division class 0. Lee and Perelson [20] calculated the analytical solutions of Eqs. (24)-(28) when the length of the B phase increases linearly with division class, Πn = Δ0 + nΔ1 for n ≥ 1 [20]. We fit the solutions to the CFSE data and the results of the fit and the estimated parameters are shown in Fig. 7. Without IL-4 stimulation, the sum of squared errors divided by the number of degrees of freedom is smaller for the model with a bi-linear increase of the death rate (4.52 ×106) compared to that for the model with a linear increase of the B phase (1.02 ×107). With IL-4 stimulation, the sum of squared errors divided by the number of degrees of freedom with the model of bi-linear increase of the death rate is 1.63 ×107 and that with the model of linear increase of the cell cycle length is 5.08 ×107. We can conclude that the model with the increase of death rate provides a better fit than the model with the increase of cell cycle length depending on the division class.

FIGURE 7
Fit of the generalized Smith-Martin model with linear increase of the length of the B phase with the division class (lines), given by Eq. (16)-(17) in Ref. [20] to CFSE data (symbols) for B cells stimulated with LPS without (a)/(c) and with (b)/(d) IL-4. ...

V. DISCUSSION

The cyton model has been formulated as a general tool for examining regulation of growth and survival of cells [24]. It assumes that cells have independent stochastic processes that determine the time to divide and to die. The molecular machinery governing cell division and death and the nature of signaling paths that modulated them have been investigated [36, 37]. The regulation of each process requires sets of molecules that interact with each other through chains of enzyme-substrate interactions. As an example of the cyton model, the lognormal distribution was assumed for the time to divide and the time to die since the lognormal distribution is typically generated by a large number of sequential multiplicative events similar to sequential enzyme-products [24].

We solved the cyton model analytically by employing the delayed gamma distribution and the delayed exponential distribution for the time to divide and the exponential distribution for the time to die. The delayed gamma distribution provides a reasonable fit to the time to enter the first division, comparable to the fit with the lognormal distribution. Gamma distributions arise from a series of consecutive stochastic events. Two random transitions in the cell cycle were previously suggested to describe the lag between stimulation and entry into S phase [38]. Delayed exponential distribution for time to divide at subsequent division classes originates from the Smith-Martin cell cycle model which successfully described the frequency distributions of intermitotic times for various cells obtained by time-lapse cinematography [12]. The generalized Smith-Martin model assumes an exponential distribution for the time to die, which enables obtaining the number of cells in each division class as a function of time in an analytical form. An exponential distribution for the time to die fits a rapid loss of small resting B cells at the initial stage, which was attributed to the manner of preparing cells [24]. The survival curve was obtained by propidium iodide update. While a gradual loss of cells over time by changing the method of preparation of cells, was described with a lognomal time to die curve [24], an exponential distribution for the time to die seems to be a reasonable approximation.

The set of distributions assumed here consists of the generalized Smith-Martin model which can be written as a set of delay differential equations. The generalized Smith-Martin model with explicit consideration of the finite length of the B phase of the cell cycle distinguishes the first division from later divisions by assuming separate distributions for the time to divide [20]. The distinct feature of the first division class with delayed gamma distribution was treated with a recruitment function in the set of delay differential equations. We verified that the numerical solutions of the set of delay differential equations for the generalized Smith-Martin model are consistent with the analytical solutions of the cyton model.

Our analysis of CFSE and thymidine labeling data of B cells in the LPS culture with the generalized Smith-Martin model indicates that the addition of IL-4 results in greater precursor frequency. The observation of greater recruitment of cells into the division by IL-4 can be understood by the measurement that active growth was associated with IL-4 mediated elimination of the cell cycle inhibitor p27 [25]. When we used a linear increase of the death rate as a function of division class, we were unable to improve the fit to the experimental data significantly in comparison with the constant death rate over the division classes. We obtained a statistically better fit to the CFSE data when we used a two-step linear increase of the death rate. The best fit death profiles had smaller death rates in the presence of IL-4 when the division class was smaller than 4. Protecting B cells from apoptosis is one of the many suggested roles of IL-4 [23, 27-29]. We may understand elevated death rates for cells that divided more than 5 times in the presence of IL-4 by the possibility that cell survival factors were consumed in the medium. The number of cells at around 90 hours after the IL-4 stimulation was doubled by adding IL-4. Dependencies of the cell kinetics parameters on the changing concentration of cytokines have previously been addressed in detail [21]. In conclusion, IL-4 has a profound effect on the death rate of B cells. On average, around 80% reduction in the average death rate up to division class 4 was estimated by adding IL-4. This observation supports the role of IL-4 as a survival factor during the proliferation of B cells.

VI. ACKNOWLEDGEMENT

This reserach was supported by NIAID/NIH contract N01-AI-50020. This work was supported in part by the University Rochester Developmental Center for AIDS Research (NIH) P30AI078498. The portion of this work was done under the auspices of the U.S. Department of Energy under contract DE-AC52-06NA25396 and supported by NIH grants P30AI078498, AI28433, P01-AI071195, and RR06555, and the Human Frontiers Science Program grant RGP0010/2004.

VII. REFERENCES

1. Lyons AB. Analyzing cell division in vivo and in vitro using flow cytometric measurement of CFSE dye dilution. J. Immunol. Methods. 2000;243:147–154. [PubMed]
2. Fazekas de St Groth B, Smith AL, Koh W-P, Girgis L, Cook MC, Bertolino P. Carboxyfluorescein diacetate succinimidyl ester and the virgin lymphocyte: A marriage made in heaven. Immnol. and Cell Biol. 1999;77:530–0538. [PubMed]
3. Gett AV, Hodgkin PD. Cell division regulates the T cell cytokine repertoire, revealing a mechanism underlying immune class regulation. Proc. Natl. Acad. Sci. USA. 1998;95:9488–9493. [PubMed]
4. Bird JJ, Brown DR, Mullen AC, et al. Helper T cell differentiation is controlled by the cell cycle. Immunity. 1998;9:229–237. [PubMed]
5. Hodgkin PD, Lee JH, Lyons AB. B cell differentiation and isotype switching is related to division cycle number. J. Exp. Med. 1996;184:277–281. [PMC free article] [PubMed]
6. Deenick EK, Hasbold J, Hodgkin PD. Switching to IgG3, IgG2b, and IgA is division linked and independent, revealing a stochastic framework for describing differentiation. J. Immunol. 1999;163:4707–4714. [PubMed]
7. Hasbold JA, Lyons AB, Kehry MR, Hodgkin PD. Cell division number regulates IgG1 and IgE switching of B cells following stimulation by CD40 ligand and IL-4. Eur. J. Immunol. 1998;28:1040–1051. [PubMed]
8. Gett AV, Hodgkin PD. A cellular calculus for signal integration by T cells. Nat. Immunol. 2000;1:239–244. [PubMed]
9. Deenick EK, Gett AV, Hodgkin PD. Stochastic model of T cell proliferation: a calculus revealing IL-2 regulation of precursor frequencies, cell cycle time, and survival. J. Immunol. 2003;170:4963–2972. [PubMed]
10. Asquith B, Debacq C, Florins A, Gillet N, Sanchez-Alcaraz T, Mosley A, Willems L. Quantifying lymphocyte kinetics in vivo using carboxyfluorescein diacetate succinimidyl ester (CFSE) Proc. Roy. Soc. B. 2006;273:1165–1171. [PMC free article] [PubMed]
11. Revy P, Sospedra M, Barbour B, Trautmann A. Functional antigen-independent synapses formed between T cells and dendritic cells. Nat. Immunol. 2001;2:925–931. [PubMed]
12. Smith JA, Martin L. Do cells cycle? Proc. Natl. Acad. Sci. USA. 1973;70:1263–1267. [PubMed]
13. De Boer RJ, Perelson AS. Estimating division and death rates from CFSE data. J. Comp. Appl. Math. 2005;184:140–164.
14. Bernard S, Pujo-Menjouret L, Mackey MC. Analysis of cell kinetics using a cell division marker: mathematical modeling of experimental data. Biophys. Journal. 2003;84:3414–2424. [PubMed]
15. Leon K, Faro J, Carneiro J. A general mathematical framework to model generation structure in a population of asynchronously dividing cells. J. Theoret. Biol. 2004;229:455–476. [PubMed]
16. Ganusov VV, Pilyugin SS, De Boer RJ, Murali-Krishna K, Ahmed R, Antia R. Quantifying cell turnover using CFSE data. J. Immunol. Methods. 2005;298:183–200. [PubMed]
17. Nordon RE, Nakamura M, Ramirez C, Odell R. Analysis of growth kinetics by division tracking. Immunol. Cell. Biol. 77:523–529. [PubMed]
18. Pilyugin SS, Ganusov VV, Murali-Krishna K, Ahmed R, Antia R. The rescaling method for quantifying the turnover of cell population. J. Theoret. Biol. 2003;225:275–283. [PubMed]
19. De Boer RJ, Ganusov VV, Milutinovic D, Hodgkin PD, Perelson AS. Estimating lymphocyte division and death rates from CFSE data. Bull. Math. Biol. 2006;68:1011–1031. [PubMed]
20. Lee HY, Perelson AS. Modeling T cell proliferation and death in vitro based on labeling data: generalizations of the Smith-Martin cell cycle model. Bull. Math. Biol. 2008;70:21–44. [PubMed]
21. Ganusov VV, Milutinovic D, De Boer RJ. IL-2 regulates expansion of CD4+ T cell populations by affecting cell death: insights from modeling CFSE data. J. Immunol. 2007;179:950–957. [PubMed]
22. Tangye SG, Avery DT, Deenick EK, Hodgkin PD. Intrinsic differences in the proliferation of naive and memory human B cells as a mechanism for enhanced secondary immune response. J. Immunol. 2003;170:686–694. [PubMed]
23. Hodgkin PD, Go NF, Cupp JE, Howard M. Interleukin-4 enhances anti-IgM stimulation of B cells by improving cell viability and by increasing the sensitivity of B cells by improving cell viability and by increasing the sensitivity of B cells to the anti-IgM signal. Cell. Immunol. 1991;134:14–30. [PubMed]
24. Hawkins ED, Turner ML, Dowling MR, van Gend C, Hodgkin PD. A model of immune regulation as a consequence of randomized lymphocyte division and death times. Proc. Natl. Acad. Sci. USA. 2007;104:5032–5037. [PubMed]
25. Wagner EF, Hleb M, Hanna N, Sharma S. A pivotal role of cyclin D3 and cyclin-dependent kinase inhibitor p27 in the regulation of IL-2-, IL-4-, or IL-10- mediated human B cell proliferation. J. Immunol. 1998;161:1123–1131. [PubMed]
26. Guo B, Blair D, Chiles TC, Lowell CA, Rothstein TL. B cell receptor (BCR) cross-talk: the IL-4-induced alternate pathway for BCR signaling operates in parallel with the classical pathway, is sensitive to Rottlerin, and depends on Lyn. J. Immunol. 2007;178:4726–4730. [PubMed]
27. Illera VA, Perandones CE, Stunz LL, Mower DA, Jr, Ashman RF. Apoptosis in splenic B lymphocytes. Regulation by protein kinase C and IL-4. J. Immunol. 1993;151:2965–2973. [PubMed]
28. Wagner EF, Hanna N, Fast LD, Kouttab N, Shank PR, Vazquez A, Sharma S. Novel diversity in IL-4 mediated responses in resting human naive B cells versus germinal center/memory B cells. J. Immunol. 2000;165:5573–5579. [PubMed]
29. Dufort FJ, Beiman BF, Gumina MR, Blair D, Wagner DJ, Roberts MF, Abu-Amer Y, Chiles TC. IL-4-mediated protection of primary B lymphocytes from apoptosis via Stat6-dependent regulation of glycolytic metabolism. J. Immunol. 2007;179:4953–4957. [PubMed]
30. Pene J, Rousset F, Briere F, Chretien I, Bonnefoy J-Y, Spits H, Yokota T, Arai N, Arai K, Banchereau J, de Vries JE. IgE production by normal human lymphocytes is induced by interleukin 4 and suppressed by interferons γ and α and prostaglandin E2. Proc. Natl. Acad. Sci. 1988;85:6880–6884. [PubMed]
31. Gascan H, Gauchat JF, Roncarolo MG, Yssel H, Spits H, de Vries JE. Human B cell clones can be induced to proliferate and to switch to IgE and IgG4 synthesis by interleukin 4 and a signal provided by activated CD4+ T cell clones. J. Exp. Med. 1991;173:747–750. [PMC free article] [PubMed]
32. Deenick EK, Hasbold J, Hodgkin PD. Switching to IgG3, IgG2b, and IgA is division linked and independent, revealing a stochastic framework for describing differentiation. J. Immunol. 1999;163:4707–4714. [PubMed]
33. Lyons AB, Parish CR. Determination of lymphocyte division by flow cytometry. J. Immunol. Methods. 1994;171:131–137. [PubMed]
34. Burnham KP, Anderson DR. Model selection and multimodel inference: A Practical-Theoretic Approach (Springer) 2002
35. Efron B, Tibshirani R. Bootstrap methods for standard errors, confidence intervals, and other measures of statistical accuracy. Stat. Sci. 1986;1:54–77.
36. Cory S, Huang DC, Adams JM. The Bcl-2 family: roles in cell survival and oncogenesis. Oncogene. 2003;22:8590–8607. [PubMed]
37. Nurse P. A long twentieth century of the cell cycle and beyond. Cell. 2000;100:71–78. [PubMed]
38. Brooks RF, Bennett DC, Smith JA. Mammalian cell cycles need two random transitions. Cell. 1980;19:493–504. [PubMed]
39. Takens F. Detecting strange attractors in turbulence. In: Rand DA, Young LS, editors. Proceedings of the Symposium on Dynamical Systems and Turbulence. Springer-Verlag; Berlin: 1981.