A calcium model with random absorption: A stochastic approach

J. theor. Biol. (1992) 154, 485-483

A Calcium Model with Random Absorption: A Stochastic Approach PALI SEN, DENIS BELL AND DONNA MOHR

University of North Florida, Department of Mathematics and Statistics, Jacksonuille, Florida 32216, U.S.A. (Received on 6 December 1990, Accepted in revised form on 6 June 1991) Absorption of calcium, or any mineral, by the body is subject to the random fluctuations typical of diffusion through membranes. In this paper we consider the absorption of calcium from the gut as a white noise process added to the deterministic model of Sen & Mohr (1990, J. theor. Biol. 142, 179-188). The first two moments for the amount of calcium in the extracellular fluid (ECF) have been derived using the Ito Calculus. A confidence interval for the total amount of calcium in the ECF is constructed. The equations for the first two moments of the fraction of dose calcium in the ECF are also given. Suggestions are made for the collection of experimental data in a form which should be helpful in investigating the magnitude of the stochastic effect. Introduction In the human body, calcium moves between various c o m p a r t m e n t s - - f r o m the gut to extraceilular fluid (ECF) and from there either to storage or out o f the body (Rubin, 1974). Absorption of dietary calcium plays a crucial role in determining the level of calcium in different sites of the body. The process of calcium transfer from the gut to ECF involves diffusion through a membrane. Such physical diffusion processes are always subject to random fluctuations. Previous models for calcium levels took no account o f these random effects. Marshall & Nordin (1969), Heideger & Ferguson (1985), and Sen & Mohr (1990, hereafter referred to as SM) all explore models where absorption from the gut is a deterministic quantity depending only on time. In this work we aim to rectify this deficiency by modeling calcium absorption as the sum of a constant and a random part assumed to consist o f white noise. This will allow us to model the stochastic, temporal variation in the calcium balance which a person may experience. Our primary goal is to construct confidence intervals for the amount o f calcium in E C F over a moderate time period. Given these intervals, we can infer dose levels which will keep calcium at or above specified levels with a stated probability. Model In SM we introduced a deterministic four-compartment model for calcium kinetics whose pathways are summarized in Fig. 1. This model generalized work by Marshall & Nordin (1969) and Heideger & Ferguson (1985). 485 0022-5193/92/040485 + 09 $03.00/0

~) 1992 Academic Press Limited

486

P. S E N comp. I

I

I

ET

AL, comp. 3

comp. 2

I

-I

I,

I,

. .-'i. . . .

eon'~ 4 FIG. I. Proposed model.

We define the following quantities. Q = amount of calcium present in the system (compartments 2 and 3) at time t = 0 D = anaount of calcium in an oral dose administered at time t = 0

Ri(t) = amount of calcium present in compartment i at time t, for i = 2, 3 .£(t) = fraction of dose present in compartment i at time t, for i = 1, 2, 3, 4. T o model the random nature o f calcium absorption, we suppose that the flow from the gut to ECF is controlled by the following differential equation' .....

a(t)f,(t)+cr~,

dt

(1)

where e, is a white noise process and W, = ~"0es ds is a Wiener process with E(W,) = 0 and Var (HI,) = ~y2t. a(t) represents the fractional absorption rate at time t, and is assumed to be given by the decreasing function a(t) = a ( 1 - t/to) for times 0 < t < to. All absorption is assumed complete by time to. Throughout the rest o f this paper the domains of all functions are the time period 0_< t < r0. The movement of calcium which has been absorbed is described as follows. The parameter fl represents the fractional transfer rate of calcium from the ECF to storage. (Here, storage refers to the pool o f calcium not in the ECF but available for exchange on slow scales of hours to a day.) 7 is the fractional return rate from storage to the ECF. From the blood pool, calcium is lost at a constant rate o f ~b per unit time. While a , / 3 and y are fractional rates, ~b is a constant rate. Flow of calcium once it has been absorbed is governed by the following equations: dR2 dt

dr,

D ~ ' - p R ~ U ) + yR3U)-,~. dt

(2)

STOCHASTIC

MODEL

OF CALCIUM

ABSORPTION

487

The solution of eqn (2) uses the identity: Ra(t) = Q + D - Dfj (t) - Rz(t) - ~bt.

(3)

Replacing R3(t) in (2) with this expression and using (1) yields: dR2 dt

(fl+z)Rz(t)+D[a(t)-7]f~(t)+[z(a+D-c~t)-c~]-Dae,.

(4)

Equation (!) can be considered as a Stratonovich integral equation driven by IV, (Gard, 1988) for which the corresponding lto form is J](t)= 1 -

a(l-s/ro)fl(s)ds+a

dW,.

(5)

The first two moments for the random variable J](t) can be obtained using the Ito Calculus, and are shown in eqns (6) and (7). Denoting E[fm(t)]=m, and E[f~(t)] =r, m, = e x p [ - a ( t

- t2/2ro)]

O
(6)

t', = tn~ + cr2m2G,

where G, =

exp [2a(s-s2/2ro)] ds=

m-2(s) ds.

(7)

The function m, is the same as the solution of the deterministic equation. As mentioned above, we are interested in the effect of the dose on the amount of calcium in ECF and in storage. We need to find the first two moments of Rz(t), which is now a random variable since df~/dt [as it appears in eqn (2)] is a random variable. E[R2(t)] =p, is derived from eqn (4) using the Ito Calculus; it is the same as the solution of the deterministic equation given in SM

de_

----(fl+r)p,+D(a-at/ro-r)m,+{r(Q+D)-4rt-(~}. dt

(8)

The equation describing the covariance, c,, offl(t) and R2(t) is dc

dt

-

[a(t)+fl+)']c,-DaZ+D[a(t)-~,],',+[),(Q+D)-d~t]m,.

(9)

This covariance function appears in the equation describing the second moment of R2(t), E[R~(t)] = v, dv

- - = -2(fl + ),)v,+ D 2 a ~ + 2D[a(t) - ),]c,+ 2{ y ( Q + D - ~bt)- ~}p,. dt

(]o)

488

v. SEN ET At..

The variance of R2(t) is a linear function of o-2, that is, as o.2 doubles, F[R2(t)] also doubles. The variance is also a linear function of the D 2. The analytical solutions of (9) and (10) are obtained using techniques for differential equations, and are given in the Appendix. A typical example of experimental data for the absorption of calcium is given in Avioli et al. (1965). A radioactive calcium tracer is administered in either an oral or intravenous dose. The concentration of the tracer in blood samples (compartment 2) is measured at several time points. The concentration is usually expressed as a percentage of the original dose, hence, the observation is usually f 2 ( t ) . 100/V2, where V2 is the volume of extracellular fluid in the body. The differential equation forJ~(t) is given by df2_ df, djq dt - - d'-t- fl~( t) + y ~ ( t) - d---t"

(11)

As in SM, we approximate df4/dt as qSf2(t)/Ro, where Ro = R2(0). Denoting E[f2(t)] as h,, we get the differential equation for h, as shown below

~by f0' h, ds. r---Ro

dh, --=dr [ a ( 1 - t / r ° ) - Y ] m ' - ( f l + T + c ~ / R ° ) h ' +

(12)

A general solution for h, yields the same result as the solution off2(t) under the deterministic equation. Non-linear least squares estimates of a, fl, y and ~b can be based on the solution for h, (as described in SM). However, h, does not involve the stochastic variance cr2. Hence, for a data-based investigation of a 2 we must look at the equation for var [f~(t)]. Solving (11) for f2(t), then taking the variance gives var [J~(t)] = exp ( - 2 b t t ) . var [h], where h=

e "r a ( r ) - ~ e -~r

eX~a(s) d s + c r X e -~r

fo

ea~dW~ d r - o r

;o

eUrdW,

;o '

(13)

and/J,/1, are the roots of the quadratic equation x 2+ Bx+ C = 0 where 8 = p + r + e~/no

c = 4,r/no

a(t) = {[a(l - t / r o ) - y]fj(t) + y} dr. Var [~(t)] is independent of dose (D) and is linear in crz. Results

The differential equations for the means and variances off~ and R2 contain the term G, (7) precluding closed form expressions for their solution. Instead, these

STOCHASTIC MODEL OF CALCIUM ABSORPTION

489

equations are solved numerically for given values of a, fl, y, ~b, ~'o and 0.2 using Runge-Kutta algorithms (Press et al., 1986). In this section, means and variances were calculated using parameter values obtained for the deterministic model as described in SM (see Table l for a summary of these values.) The variance o.2= 0.001 was decided by setting the standard deviation of J~(t) at t = 4 hr to be 0-05 (i.e. the amount unabsorbed might typically vary in a range +10% about the mean). TABLE l Parameter oalues used in numerical results a =0.463 fl = 1-457 7 = 0' 336 ~ = 5.083 R0= 1039.35 ro= 4 o.2=0.001 D= 100 Values of a, fl, ?', ~ and Ro are based on fitted values for normal subjects given in SM.

Figure 2 shows the standard deviations o f f , ( t ) and R2(t) as functions of time. Both functions grow most quickly in the first half-hour of time, when absorption occurs at the highest rate. While the standard deviation o f / ] continues to grow rapidly, the standard deviation of R2(/) grows only slowly after time t = 0.5. Figure 3 shows the expected value of R2(t) =Pt. The graph also shows an approximate 95% confidence interval for the value of R2(t) obtained using p,+2[var (R2(t))] °'5. Note that a change in o- will not affect the value of p, or the shape of the confidence interval, but will scale the width of the interval by the same factor as the change in or. Both the expected value of R2 and the variance of R2 depend on the dose D, which was set to 100 mg. Notice that the standard deviation of R2 increases rapidly during the first hour of time, and then stabilizes, as reflected in the width of the confidence interval. The width of the interval is large (about +4 mg) compared to the maximum expected change in R2 (about 14 mg at time 1 hr) and compared to the net change 4 hr after the close was administered (about 6 mg). While the standard deviation for f2 could be computed numerically, it is more practical to use a simulation to estimate it. Two simulated sample paths for ~ ( t ) . 100/V2 are shown in Fig. 4. These paths were generated from a discrete-time simulation with time-step=0.001. One can see that the serum concentration of the dose tends to peak at about time t-- l hr. Increasingly wide swings in the path begin to occur at about this time. The estimated standard deviation for f2 shown in Fig. 5 is based on 1000 repetitions of the simulation, and is scaled to the same units as those reported in the Avioli data (1965). The pattern is very similar to the pattern for the standard deviation of R2. [The simulation also produced estimates of the means and variances for3q(t) and R2(/) which checked closely with the numerical solutions for the equations.]

490

P. S E N

E T A L.

0-05

0.04

0-05

0'02

(3"01

I

I

2

5

,,

4

2-50

/

2.00

1"50

d

1-00

Y

0-50

0

t

i

2 Time (hr)

:5

4

FIG. 2. Standard deviation o f f , U) and R2(t) computed using parameter values in Table 1.

Discussion

Results such as those displayed in Fig. 3 allow the practitioner to appreciate the variability in the effect dietary calcium intake will have on serum calcium levels. This variability is quite large compared to the change in the expected value of R2, even for very modest values of cr2. The interval allows for many of the uncontrollable environmental processes which affect the system. Thus, a stochastic model, like that in the present paper, represents the temporal variation caused by random influences on the absorption process. We remark that the mean of Rz(t) is an affine function of D, while the variance of R~(t) is quadratic in D. This enables us to determine a level for the dose, D, so

STOCHASTIC

MODEL

OF

CALCIUM

491

ABSORPTION

1060

I056

1052

~1o48 I//

"~ -...

1044

/ 1040

1036

I I

I 2

I 3

Time (hr)

FIG. 3. Expected value and confidence intervals for 1. ( ), Mean o f R~ ; ( . . . . ), mean + 2 S.D.

R2(t)

c o m p u t e d using parameter values in Table

2"00

od 2~ o o

0

L.

I

I

2 Time (hr)

I

3

4.

FIG. 4, Two simulated sample paths o f fz(t), computed using parameter values in Table l, ( ~ ~ ~ ~ ) , Sample path 1; (~ ¢ ~ *), sample path 2.

492

P. SEN" E T A L . 0,20

oa 0.15

O"

o OJO

0.05

0

I .... I

I 2 Time

I 3

4

(hr)

FIG. 5. Estimated standard deviation off2(t) using parameter values in Table 1. Estimation is based on 1000 simulated sample paths.

that the lower bound of the resulting confidence interval for the serum calcium level at any given time exceeds a preassigned amount. Hurwitz et al. (1987) showed that absorption, loss and turnover rates are functions of the calcium level, R2. Here, we dealt with only one source of variation, that which occurs in the absorption process. More complex models which consider all sources of variability rapidly become intractable. This model is simple yet informative. It gives projected calcium levels in the serum at any time after absorption has stopped. It can be used for further study, if the absorption rate is a function of diet, age, or other personal factors. The model is rigorous to statistical moment theory, incorporating the equations for the stochastic integrals and their solutions. The interval estimates for R2(t) require some knowledge of an appropriate value for 0-2. In a research setting, it should be possible to administer an oral dose of a radioactive tracer under controlled conditions, and observef2(t) at discrete intervals (see e.g. Avioli, 1965). This process should be repeated on several occasions for a single individual. The mean f2(t) can be used to estimate this person's values for a, fl, 7' and ~b. The variation in J~(t) could conceptually result in an estimate for 0 -2, using a numerical approximation to the solution of (13). While this function is very complicated, we need to compute or simulate it only once using the already estimated a,/3, T and ~, and a convenient initial choice of 0-2, say 0-~ =0" I. Taking advantage of the fact that Var [J~(t)] is linear in 0-2 we then choose a new 0-2 to rescale this function to the observed variances. Alternatively, the variation in R2(t) around its predicted level could result in an estimate of tr 2 using the results given in the Appendix. However, though estimates of R2(t) are available from 5 cc of blood, a dose has only a small percentage effect on R2, so that measurement error may make this approach impractical. In the past, pharmacokinetic methods dealt primarily with deterministic models. In recent years, however, researchers have begun using stochastic models, as these

STOCHASTIC

MODEL

OF CALCIUM

ABSORPTION

493

provide a b r o a d e r theoretical f r a m e w o r k . As d e m o n s t r a t e d in this paper, such models allow the investigator to examine the m a g n i t u d e o f variabilities generated by u n o b servable factors affecting the kinetics o f the system, and should therefore be increasingly i m p o r t a n t in the study o f biological systems. REFERENCES AV1OLI,L. V., McDONALD,J. E., SINGER, R. A. & HENNEMAN, P. H. (1965). A new oral isotopic test of calcium absorption. J. clin. Invest. 44, 128-139. GARD, T. C. (1988). Introduction to Stochastic Differential Equations. New York: Marcel Dekker. HEIDEGER, W. J. & FERGUSON, M. E. (1985). A theoretical model for calcium absorption from the intestinal lumen. J. theor. BioL 114, 657-664. HURWITZ, S . FISHMAN, S. & TALPAZ,H. (1987). Calcium dynamics: a model system approach. J. Nutr. 117, 791-796. MARSHALL, D. H. & NORDm, B. E. C. (t969). Kinetic analysis of plasma radioactivity after oral ingestion of radiocalcium. Nature, Lond. 222, 797. PRESS, W. H., FLANNERY,B. P., TEUKOLSKY,S. A. & VE'I-rERLWNG,W. Z. (1986). Numerical Recipes: The Art of Scientific Computing. Cambridge: Cambridge University Press. RualN, R. P. (1974). Calchtm and the Secretory Process. New York: Plenum Press. SEN, P. & MOHR, D. (1990). A kinetic model for calcium distribution. J. theor. Biol. 142, 179-188. APPENDIX Let k = f l + y , p,=E[R2(t)], c , = E [ f l ( t ) . R2(t)] and o~=E[R~(t)]. Using the Ito Calculus, we derive

P,= [ fl(Q k D )+d:

k2

) ~ ] e -k' -t ~ ( Q + kD ) - q b R3(0)_._£7" Din, + DflZ( t) e -r°[t*-")2 - k t l / 2 a

where Z ( t ) = S'o exp [ a ( u + r o ( k - a)/a)2/2ro] du

(a - a s ~ t o - 9,) e 1.'- ~ + ~':/2~01 ds

c , = D e -I~ +k~'+ ~':/2r° t-~ 0

_Dcr 2

+ Oct z 2

:o

et~ +kJ,--,2/2~ol ds +

f'

[7(Q+D-dps)-d~]ek'ds

~0

fo

(a - a s ~ t o - y) e lk'- ~"+ ~'v2~°la(s) ds + R2(0)

[-e2k'

v,=O tr2[Tk

1-1

fr

lj+2De-'k'

( a - a s / r o - ? , ) e 2 " S c , ds ~0

+ 2 e -u'

f, ~0

[y(Q+O-4~s)-4~leZ~'p, ds+[R2(0)l~e -2k'.