- Email: [email protected]
Cell survival probability under ionizing radiation
Cell Survival Probability Under Ionizing Radiation GRACE
L. YANG
Department of Mathematics, Uniuersity of Maryland, College Park, Maryland 20742 Received
7 January 1992; revised 24 August 1992
ABSTRACT The survival probability of a living cell exposed to ionizing radiation in an experimental setup is derived. The survival of a cell depends on the severity of the radiation damage and efficiency of the cellular repair. The formula of the survival probability is expressed as a function of dose, nonlinear rate of lesion induction, nonlinear rate of cellular repair, and a key experimental parameter-the holding time. The result is an extension of the Markovian dose-response model developed by Yang and Swenberg.
1.
INTRODUCTION
This presentation is concerned with the calculation of the survival probability of a living cell exposed to ionizing radiation in an experimental setup. Most cells have the capability of repairing the induced radiation damage. But the repair is not 100% efficient. A cell can survive the radiation damage if lesions induced are repaired completely and correctly. Incorrect repair (misrepair) will result in either cell transformation or cell inactivation, by which we mean that the cell is unable to divide. It is a well-known experimental fact that the survival probability S(D) of cells exposed to sparsely ionizing radiation, such as X or y rays, exhibits a distinct change of slope from zero to exponential as dose D increases. In other words, log S(D) is a quadratic function of dose D. This is known as the shoulder effect (see Figure 1B). This nonlinear phenomenon is also present when either the dose rate decreases or the allowable cell repair time increases up to a certain level (see, e.g., Elkins [61, Frankenberg-Schwager et al. [9, 101, Frankenberg et al. [S], and Metting et al. [151). With the shoulder effect used as a guideline, different biophysics theories and mathematical models have been proposed for studying the dose-response relationship. The repair and MATHEMATICAL
BIOSCIENCES
112:305-317
(1992)
OElsevier Science Publishing Co., Inc., 1992 655 Avenue of the Americas, New York, NY 10010
305 0025-5564/92/$5.00
306
GRACE L. YANG
A
-Single
Strand
10’ r
10-33 0
400
800
1200
Dose (rads) FIG. 1. (A) Illustrates a number of different types of radiation-induced lesions in DNA. (B) Typical survival curves for mammalian cells exposed to radiation. At high LET (densely ionizing) radiation the logarithmic survival fraction versus dose is linear. At low LET (sparsely ionizing) radiation the shoulder effect is present in the logarithmic survival curve.
misrepair model of Tobias and coworkers [17, 211 attributes the shoulder to the cell’s particular repair mechanism. The dual radiation theory of Kellerer and Rossi [12] and the linear-quadratic model of Chadwick and Leenhouts [3] use the track geometry of the radiation particles and proximity of lesions to explain the shoulder effect. Iliakis and Pohlit [ill use nonlinear rates of lesion initiation as a cause for the shoulder. The literature is extensive on the subject. References on mathematical
CELL SURVIVAL
PROBABILITY
UNDER RADIATION
307
modeling include the works of Albright [l, 21, Le Cam and Neyman [14]; Neyman and Puri [16], Tobias [20], Curtis [5], Sachs and coworkers [18, 191, and Yang and Swenberg [23, 241 and references quoted therein. In Yang and Swenberg [24], a stochastic model is constructed that embodies the lesion initiation events, the cell repair mechanism as well as a key experimental parameter-the holding time. A detailed discussion of the model and comparisons with the others are given in [24]. The cell survival probability S(D) in [24] is calculated under the assumption of nonlinear lesion initiation rate and linear cellular repair rates. The model was shown to fit several sets of experimental data reasonably well. However, neglecting the possibility of cells having a nonlinear repair rate clearly limits the applicability of the model. An example of lack of fit of the model is provided in this volume by Le Cam [14], who shows that a modified version of the model does not fit the Frankenberg-Schwager experiments [9]. A plausible cause for this is the omission of the nonlinear repair rate. The exclusion of nonlinear repair rates in [241 was due to mathematical intractability of a nonlinear Kolmogorov equation. In this paper we are able to use an alternative, probabilistic approach to derive the cell survival probability that assumes both nonlinear lesion initiation and repair rates. A formal solution, though complicated, is presented. The paper is organized as follows. Section 2 contains a brief description of the radiation experiments. Section 3 presents the mathematical assumptions and the models. The derivation of the cell survival probability is given in Section 4. Concluding remarks are in Section 5. 2.
RADIATION
EXPERIMENT
To calculate cell survival probability, it is necessary to describe how cell survival is determined in a radiation experiment. The survival of a cell is not directly observable after irradiation as one cannot distinguish a living cell from a dead cell. A practical criterion for determining cell survival is that it must have the proliferative ability to form a colony of a given (observable) size within a specified time after irradiation. Thus to determine survivability, after irradiation, the cell has to be plated onto a petri dish filled with nutrients to allow for its proliferation. In some experiments there is an intermediate step, in which case the cells are held in a nongrowth environment (e.g., no nutrients) for a period of time before plating. This period is called the holding time. Experimental data have demonstrated that a substantial amount of cellular repair occurs during the holding time, which is therefore a key parameter in studying cell repair mechanisms. As such, a typical radiation experiment is composed of three consecutive time periods: the irradiation period (0, T), the holding period (T, T,), and the cloning period (7”, T,). Cells of
308
GRACE L. YANG
a given type are irradiated during (0, T) at a specified dose rate corresponding to a given total dose, D. Upon termination of the radiation, the cells are transferred to an environment that allows for repair of radiation-induced lesions but not for cell division. Some induced lesions are lethal and irreparable. The reparable lesions are referred to as “potentially lethal lesions” for the reason that incorrect repair may still result in cell death. The correct repair as well as the incorrect repair of lesions that causes cell transformation will result in cell survival. This will be discussed in Section 3 under Assumption A2 in the construction of stochastic models. At time T,, cells are plated for cell division. The clones are counted at time T,. The counting time T, is generally selected such that colonies have more than 50 cells and are observable by the naked eye. Each colony is taken as having originated from a single cell. The ratio of the number of colonies and initial number of irradiated cells, corrected for plating efficiency, measures a cell’s survival fraction for a given radiation dose. 3.
ASSUMPTIONS
AND
STOCHASTIC
MODELS
The following mathematical assumptions are made with respect to the occurrences of events in the experiment. Assumptions Al-A3 are the same as those given in [24]; we set the plating efficiency parameter 6 = 0 for simplicity. But Assumption A4 is different. In A4 a repair-misrepair model with a nonlinear repair rate K is introduced. For computing survival probability, alternative formulations different from A4 are needed; they are given in A5 and A6. Assumption Al. Poisson am’val of primary radiation particles. The primary radiation particles arrive at the cell nucleus according to a Poisson process with rate A(t) per unit time. A(t) is a positive continuous function of t. (It is known that the passage of particles in the cytoplasm has little effect on survival; see, e.g., Von Borstel and Rogers [22].) Then the expected number of primary particles arriving in [O, t 1is l,‘A(x) du. Assumption A2. Formation of lesions. Each primary generates a random number, M, of “spurs” with a probability-generating function (pgf) g(s) = Es”. For low-LET radiation, all spurs act independently of one another (Chatterjee et al. [41X Each spur has a probability 7~~ of generating a potentially lethal lesion and a probability rrz of generating an irreparable lesion that inactivates the cell (a lethal event). We have 0 < 7~~+ rTT2< 1, where 1 - rr, - nz is the probability that the spur has no effect on the cell. The values of 7~~ and z-Z may depend on the time at which the spur is created through the total accumulated dose received by the nucleus up to time t. Through this dependence a nonlin-
CELL SURVIVAL
PROBABILITY
UNDER RADIATION
309
ear rate of lesion initiation is introduced into the model. The effect of these M spurs is described by random variables U, and U,, where U1 is the number of potentially lethal lesions and U, is the number of lethal lesions induced by the M spurs. Under Assumption Al, the pgf of Ut and U, is given by E( I&@
) = E( 1- rr, - 7rz + %-IU, + rr2U*)M = g( 1- 9-r, - T2 + 7T,u1 + 7r*u2)
(1)
for lull G 1, luzI < 1.
Assumption A3. Colony-forming ability. We assume that a single lesion will inactivate the cell unless it is correctly repaired or is misrepaired to a mutant (transformed cell). We also assume that repair processes stop at the time of plating TP. Thus for a cell to divide and grow into a visible colony it is necessary that by time TP all potentially lethal lesions are either correctly repaired or misrepaired in transformation. Assumption A4. Repair - misrepair mechanisms. Radiation produces a variety of cellular lesions that immediately activate the cell’s repair mechanisms. Figure 1A illustrates a variety of radiation-induced lesions. Some of these lesions are reparable, and some are not. A stochastic model that describes the outcomes of lesion repair is as follows. The evolution of a single cell during and after radiation is a vector-valued continuous-time Markov process (X,Y, 2) = {(X,,Y,, Z,>; t 2 01, where X, is the number of potentially lethal lesions that the cell has at time t, Y, is the number of transformed lesions that the cell has at time t, and Z, is the number of lethal lesions (inactivation events) the cell has experienced up to time t. At the onset of the experiment, t = 0, no lesions exist; that is, X0 = Y,, = Z,, = 0. During the time interval (t, t + h], either nothing happens or one and only one of the following transitions takes place:
CX,,y,,Z,)
+a,
-l,Y,,Z,)
with probability
[ax,
+ .X,(X,
-
I)1 h + o(h) (correct
ar(X, - l,Y, + l,Z,>
with probability
-+(X,-l,Y,,Z,+l)
withprobabilityyX,h+o(h)
PX,h + o(h) (misrepair in transformation), (misrepair
-+ (Xt+u,,Yt,
Z,+u,)
with probability for uI,u2
repair),
hhP[U,=u,,
in inactivation),
U2=u21+dh),
= O,l;... (initiation
of new lesions),
310
GRACE L. YANG
where (Y, p, y, and K are nonnegative parameters that may depend on the time t and dose rate r. This is the model used in [24] except for the introduction of a nonlinear term K. Insofar as cell survival is concerned, it is not necessary to consider Y explicitly because, according to Assumption A3, the cell survival probability is equal to S = P[X, = 0, Z, = 01, which does not involve Y,. In other words, Y, does not contribute to cell death. This means that we can add the transition rates (Y+ p in computing the survival probability. In view of this, it is sufficient to consider the process (X, Z) = ((X,,Z,>; t a 0) where the transition of (X,,Z,) to (X, - l,Z,) is modified to (a + p>X, + KX,
W,,Z,) with probability
+(X,-1,-Z,>
(correct -+(X,-1,z,+11
repair
[ax,
+ KX,
and misrepair
with probability
yX,h + o(h)
with probability
hhP[U,
in transformation),
(misrepair -(X,
+ u,,ZI
+ u,)
for u,,uz
= ul, lJ, =
in inactivation),
u21= + o(h),
=O,l;... (initiation of new lesions),
where (Y, y, and K are nonnegative parameters that may depend on the time t and dose rate r. As explained in A4, (Y should be replaced by (Y+ p. As it is only a symbol for the parameter, we shall for simplicity continue using (Y even though it represents “(Y + p.” It is convenient to formulate the cell survival probability, S, in terms of the process (X,Z). But, in the presence of the nonlinear repair rate (K > 01, it is difficult to compute S from the (X, Z) process. We shall use a different probabilistic approach to compute S. This requires an alternative formulation of the process (X,, Z,>, t a 0, which is given in A6. It is easy to show that the Markov process (X,,Z,) describes the same repair process as that given in Assumption A6.
CELL SURVIVAL
PROBABILITY
UNDER
311
RADIATION
Assumption A6. Let N, be the Poisson number of primary particles that arrived at the nucleus in the time interval (O,tl, and let r;, for j=l ,.. . , N,, be their arrival times. At time rj, U[ potentially lethal lesions and Ui lethal lesions are created according to A2. Fix j. Assume that rj and U{ are given, for j = 1,. . . , N,. The cellular repair of the l.I/ lesions starts at time 5. The repair process is modeled by a two-dimensional Markov process {Rj(t - rj), y ~~1,where R,(t - TV>is the number of potentially lethal lesions that remain at time t and y-ti equal to zero if t < -rj; likewise for W. For simplicity, we set
Rj( t - TV)= R, The transition
and
qt
process
are
rates of the repair
(R,,U:)+(R,-LY+l)
-Tj)=w,.
with rate yR, with rate aR, + KR,( R, - 1).
-(R,-I,%)
Given N, and rj, for j = 1,. . . , N,, the repair ~~1,E$t - ~~1); t a T~) are stochastically independent.
processes
((Rj(t
-
In summary, under Assumptions Al-A4 (or with A4 replaced by A5 or A6), we can visualize several stochastic processes (not necessarily observable) that are evolving simultaneously or sequentially over time t. They are (1) the Poisson arrivals of the primary radiation particles, (2) the induction of lesions in the nucleus, (3) the Markovian repairmisrepair process during the irradiation and holding periods, and (4) the cell proliferations or cloning. These assumptions are clearly an oversimplification of the actual process. For instance, we have made no distinction of different kinds of radiation damage illustrated in Figure 1A. The reader is referred to [241 for further discussion and to a survey article by Sachs et al. in this volume for discussions of DNA damage and Markov models. 4.
COMPUTING Under
CELL
SURVIVAL
A3, the cell survival
probability
PROBABILITY is defined
s(D)=P[xTp=o,zT~=o]. Note that it is written through h(t) in Al.
as a function
by
(2)
of dose D. The dose is determined
312
GRACE L. YANG
Consider first the survival probability of a cell that was hit by one primary particle, say at time 0. Suppose that the induced number of potentially lethal lesions is ur and the number of lethal hits U, = 0. This initiates a repair process CR, W) with initial conditions R(O) = u1 and W(O) = 0. Under Assumption A6, the state space of CR, W) is a triangle with vertices (O,O), (u,,O), and (0, u,>. Let q,=
KT(T--l)+yr
CZr+
and
(3) for r = l,...,u,. Let p, be the waiting time of the process CR, W) in state (r, w>. Then p, is an exponential random variable with density
fW=s,exp(-q,xL
for x > 0 and Y > 0.
Note that the distribution of p,. does not depend on W. Thus the total cell repair time (correct repair and misrepair combined) is the sum EFr 1 p, of u, independent random variables p,. According to definition (21, the survival probability in the present case is S(D) = P[ R, = 0, W, = 01. It is the probability that the process (R, W) completes [he path {(u,,O) -+ (ul - 1,O) + (ul -2,0) + ... + (O,O)} by the time of plating, Tp, By the Markovian assumption A.5, P[(u,,O)
-+(u,
-130)
+(u*
-210)
+ ... + (0,O) by time T,lR(O) =P
Consequently,
g [ r=l
pr”Tp
= ul, W(0) = O]
&(l,u,).
(4)
I
the survival
probability
(5)
= where -
Remark. According created
= with A2,
=
U, to
initiation and the
0] and
the is
may accumulated
P(u,,O) on
in
derivative with time received
g(l-
-
= on and which spur the
CELL SURVIVAL
PROBABILITY
313
UNDER RADIATION
up to time t. In the present case, 7, is set to zero. This corresponds to the experimental situation of acute irradiation of cells for a very brief time period. We can set T = 0 and consider that lesions are produced instantaneously at time 0. In this case, the rr’s will be constant independent of time. For instance, rTTIstudied in [24] has the form 5 + r,[l exp( - at)], which reduces to 5 at t = 0. In the derivation of (5) we have separated the sample path from the time it requires to complete the path and thus avoided having to solve a differential equation. Formula (5) is valid for any transition rates q, and the repair rates, not necessarily quadratic functions of r. In the case of (31, q, have positive distinct values. The sum CF; 1 p, in the middle has density qlq2...qu,[$,,n
exp( -
91x) + -** + k,n exP( -
4u,x)l, (6)
+zl =(41-qk)...(qk~l-qk)(qk+l-qk)...(qu,-qk) 3
I
fork=1
,..., u,.
See, for example, [7]. If all the q, are equal, then CAL, p, has a gamma density. For other cases of q, the density of the sum can be obtained by inverting Laplace transforms. However, the formulas are messy. The repair process considered so far deals with potentially lethal lesions created under acute irradiation. This is similar to the repair-misrepair (RMR) model 1211. A stochastic version of RMR model is given by Albright [l, 21. The solution for the nonlinear case in finite time studied here is not given in Albright. We now to return to the general situation and investigate the repair processes generated by all the lesions (Uj’,Ui), for j = 1,2,. . . , created at times O
K4Tl. It is convenient to consider the cell’s repair activities separately in each time interval CT,,rj+ ,> between the arrivals of newly radiationinduced lesions. The repair begins at time r1 upon the arrival of ZJ, potentially lethal lesions. The repair is described by the process CR,, W,) with the initial condition R,(O) = u!, W,(O) = 0. Recall that a cell cannot survive any lethal hit or with any misrepaired lethal lesion. The cell survival probability at time t E CT,, TV) is
P[R,(t-71)
=x1,
W,(~-T~)=OIR~(~)=~~,W~(O)=~].
(7)
GRACE
314
L. YANG
At t = TV, there are new arrivals of lesions, Uf. Thus, the number potentially lethal lesions that the cell has at time r2 is
of
&(O)=&(r2-q)+q5
(8)
where the summands are independent. In general, the number of potentially at time 7j is
lethal
Rj(O)=Rj~,(7/-Tj_,)+U:”
lesions
j=1,2
that the cell has
>.‘.>
(9)
which together with II$(O) = 0 form the initial condition of the repair process (Rj(t - ~~1, Fi$t -TV)), for t E (7,,rj+ 1). Note that (Rj(t - ~~1, y(t - 7,)) has transition rates similar to that of (R, W) but with possibly different initial conditions, The probability distribution of the initial conditions Rj(0) and II$O) depend on the lesion initiation probabilities 7~~ and rTT2,which by Assumption A2 depend on TV. we are interested m the probability For each j and t E(7j,Tj+l), P[ Rj(t Denote
- TV) = rj, H$(t -TV) = OIlSi
the jth initiation
probability
P(z&O)
= rj_, + 241, y(O)
= 01.
by
= P[ U:’ = Ui, u: = 01.
By Assumption A2, this can be derived from (l), where the rr’s in this case depend on TV.We can immediately obtain the following conditional cell survival probability at time t > T: P[X(t)=O,Z(t)=OIN,=n,~
=c
all r,
,,...,
Tn]
n-1
i
,‘II, (+0)x
xp?,(t
p
-~)=rj,~(t-T,)=OIR,(0)=rj_,+u:,~(O)=O])
x P( u;t,O) x P x[R,(t-7,)=O,W,(t-7,)
=OlR,(O)
=rnpl +u;,W,(o)
=o]
) 1 (10)
where
the summation
extends
to all 0 < rj < rj_ 1 + u{, for j = 1,. 1s7n.
CELL SURVIVAL
PROBABILITY
UNDER RADIATION
The last term in (10) can be obtained from (3). It remains the middle term. We suppress j for simplicity. Let
Si-.u,
=
315
to compute
c
pl.
i=r
Denote (u,-l,O) = rl=
the event by B(r,u,) that the process moves from to (r,O). Using the fact that on the set B(r,u,), ,...,
Applying
(u,,Ol to [R(t -7)
Eq. (6) yields the result
P[S,,,,
- fig, [ 2 I-. ,=r
1
~i~~~‘{l-exP[-qi(f-7)l} 3
i=l
where
*F’ = iii(qj-qi) j=r jti
and P[B(r, u,)] = [(r, u,>, as given in (3). The unconditional survival probability is obtained by integrating the distributions of the T’S and the N(t) in (10). Formally,
S(D)=
c
Tl>O
mv)
=A
h(T)”
T... /0
out
7. /0
XP[X,=Z,=OIN,=~,T,,...,T,]
ifIh(~~)d~~
"'dT,,,
(11)
survival
prob-
i=l
where R(T) = /$A(y)dy. ability.
Setting
t = Tp yields the desired
316 5.
GRACE L. YANG
CONCLUDING
REMARKS
We have derived an explicit expression for the cell survival probability under ionizing radiation under the assumption of nonlinear rates of repair and lesion initiation. In fact, the nonlinear repair rate can be made quite arbitrary, not necessarily a quadratic function as is used in A5. This offers some flexibility in selecting appropriate repair rates to explain the cell repair mechanism. The result is an extension of what is obtained in [24]. Our formula makes explicit use of the holding time TP. TP is a key experimental parameter for studying the cell repair mechanism. It is, however, mostly neglected in the mathematical modeling literature. The survival probability derived here is a complicated analytical expression that requires heavy numerical calculation. It remains to be investigated whether the model will fit the experimental data. The author wishes to thank L. Le Cam for many helpful discussions. Thanks are also due to a referee for valuable suggestions that led to a substantial improuement of the paper. This paper was written while the author was visiting MSRI in Berkeley. Research at MSRI was supported in part by National Science Foundation grant DMS 8505550. REFERENCES N. W. Albright, A Markov model of the repair-misrepair process of cell survival, Lawrence Berkeley Laboratory Rep. LBL-11660, Univ. California, 1981. N. W. Albright, A Markov formulation of the repair-misrepair model of cell survival, &z&t. Res. 118: l-20 (1989). K. H. Chadwick, and H. P. Leenhouts, The Molecular Theory of Radiation Biology, Springer-Verlag, New York, 1981, p. 102. A. Chatterjee, P. Koehl, and J. L. Magee, The theoretical consideration of the chemical pathways for radiation-induced strand breaks, A&. Space Res., 6( I1):97-105 (1986). S. B. Curtis, Lethal and potentially lethal lesions induced by radiation-a unified repair model, Radiut. Res. 106:252-270 (1986). M. M. Elkind, Repair processes in radiation biology, Radiut. Res. 100:425-449 (1984). W. Feller, An Introduction to Probability Theory and Its Applications, Vol. 2, Wiley, New York, 1966. D. Frankenberg, M. Frankenberg-Schwager, and R. Harbich, Interpretation of the shapes of survival curves in terms of induction and repair-misrepair of DNA double strand breaks, Br. J. Cancer 49cSuppl. 6):223-238 (1984). M. Frankenberg-Schwager, D. Frankenberg, D. Blocher, and C. Adamczyk, Repair of DNA double strand breaks in irradiated yeast cells under new growth conditions, Rudiat. Res. 82:498-510 (1980).
CELL SURVIVAL 10
11 12 13 14 15
16 17
18 19 20 21
22 23
24
PROBABILITY
UNDER
RADIATION
317
M. Frankenberg-Schwager, D. Frankenberg, and R. Harbich, Exponential or shouldered survival curves result from repair of DNA double strand breaks depending on postirradiation conditions, Radiat. Res. 11454-63 (1988). G. Iliakis and W. Pohlit, Quantitative aspects of repair of potentially lethal damage in mammalian cells, Int. J. Radiat. Biol. 36(6):649-658 (1979). A. M. Kellerer and H. H. Rossi, A generalized formulation of dual radiation action, Radiat. Res. 75:471-488 (1978). L. Le Cam, Stochastic models of lesion induction and repair in yeast, Math. Biosci., this volume. L. Le Cam and J. Neyman, Eds., Probability Models and Cancer, North-Holland, New York, 1982. N. F. Metting, L. A. Braby, W. C. Roesch, and J. M. Nelson, Dose rate evidence for two kinds of radiation damage in stationary-phase mammalian cells, Radiat. Res. 103:204-218 (1985). J. Neyman and P. S. Puri, A hypothetical stochastic mechanism of radiation effects in single cells, hoc. Roy. Sot. Land. B 213:139-160 (1981). F. Q. H. Ngo, E. A. Blakely, C. A. Tobias, P. Y. Chang, and L. Lommel, Sequential exposures of mammalian cells to low- and high-LET radiations. II. As a function of cell-cycle stages, Radiut. Res. 115:54-69 (1998). R. K. Sachs and L. R. Hlatky, Stochastic dose-rates in radiation cell survival models, Radiat. Emiron. Biophys. 29:169-184 (1990). R. K. Sachs, P.-L. Chen, P. Hahnfeldt, and L. R. Hlatky, DNA damage caused by ionizing radiation, Math. Biosci. (1992). C. A. Tobias, The repair-misrepair model in radiobiology: comparison to other models, Radiat. Res. 104:S-77-S-95 (1985). C. A. Tobias, E. A. Blakely, F. Q. H. Ngo, and T. C. H. Yang, The repair-misrepair model of cell survival radiation biology and cancer research, A. Meyn and R. Withers, Eds., Raven Press, New York, 1980, pp. 195230. R. C. Von Borstel and R. W. Rogers, Alpha-particle bombardment for the Habrobrazon egg. II. Response of the cytoplasm, Radiat. Res. 8:248-253 (1958). G. L. Yang and C. E. Swenberg, Stochastic models for cells exposed to ionizing radiation, in Modeling of Biomedical Systems, J. Eisenberg and M. Witten, Eds., Elsevier, Amsterdam, 1986, pp. 85-89. G. L. Yang and C. E. Swenberg, Stochastic modeling of dose-response for single cells in radiation experiments, Math. Sci. 16:16-65 (1991).
L. YANG
Department of Mathematics, Uniuersity of Maryland, College Park, Maryland 20742 Received
7 January 1992; revised 24 August 1992
ABSTRACT The survival probability of a living cell exposed to ionizing radiation in an experimental setup is derived. The survival of a cell depends on the severity of the radiation damage and efficiency of the cellular repair. The formula of the survival probability is expressed as a function of dose, nonlinear rate of lesion induction, nonlinear rate of cellular repair, and a key experimental parameter-the holding time. The result is an extension of the Markovian dose-response model developed by Yang and Swenberg.
1.
INTRODUCTION
This presentation is concerned with the calculation of the survival probability of a living cell exposed to ionizing radiation in an experimental setup. Most cells have the capability of repairing the induced radiation damage. But the repair is not 100% efficient. A cell can survive the radiation damage if lesions induced are repaired completely and correctly. Incorrect repair (misrepair) will result in either cell transformation or cell inactivation, by which we mean that the cell is unable to divide. It is a well-known experimental fact that the survival probability S(D) of cells exposed to sparsely ionizing radiation, such as X or y rays, exhibits a distinct change of slope from zero to exponential as dose D increases. In other words, log S(D) is a quadratic function of dose D. This is known as the shoulder effect (see Figure 1B). This nonlinear phenomenon is also present when either the dose rate decreases or the allowable cell repair time increases up to a certain level (see, e.g., Elkins [61, Frankenberg-Schwager et al. [9, 101, Frankenberg et al. [S], and Metting et al. [151). With the shoulder effect used as a guideline, different biophysics theories and mathematical models have been proposed for studying the dose-response relationship. The repair and MATHEMATICAL
BIOSCIENCES
112:305-317
(1992)
OElsevier Science Publishing Co., Inc., 1992 655 Avenue of the Americas, New York, NY 10010
305 0025-5564/92/$5.00
306
GRACE L. YANG
A
-Single
Strand
10’ r
10-33 0
400
800
1200
Dose (rads) FIG. 1. (A) Illustrates a number of different types of radiation-induced lesions in DNA. (B) Typical survival curves for mammalian cells exposed to radiation. At high LET (densely ionizing) radiation the logarithmic survival fraction versus dose is linear. At low LET (sparsely ionizing) radiation the shoulder effect is present in the logarithmic survival curve.
misrepair model of Tobias and coworkers [17, 211 attributes the shoulder to the cell’s particular repair mechanism. The dual radiation theory of Kellerer and Rossi [12] and the linear-quadratic model of Chadwick and Leenhouts [3] use the track geometry of the radiation particles and proximity of lesions to explain the shoulder effect. Iliakis and Pohlit [ill use nonlinear rates of lesion initiation as a cause for the shoulder. The literature is extensive on the subject. References on mathematical
CELL SURVIVAL
PROBABILITY
UNDER RADIATION
307
modeling include the works of Albright [l, 21, Le Cam and Neyman [14]; Neyman and Puri [16], Tobias [20], Curtis [5], Sachs and coworkers [18, 191, and Yang and Swenberg [23, 241 and references quoted therein. In Yang and Swenberg [24], a stochastic model is constructed that embodies the lesion initiation events, the cell repair mechanism as well as a key experimental parameter-the holding time. A detailed discussion of the model and comparisons with the others are given in [24]. The cell survival probability S(D) in [24] is calculated under the assumption of nonlinear lesion initiation rate and linear cellular repair rates. The model was shown to fit several sets of experimental data reasonably well. However, neglecting the possibility of cells having a nonlinear repair rate clearly limits the applicability of the model. An example of lack of fit of the model is provided in this volume by Le Cam [14], who shows that a modified version of the model does not fit the Frankenberg-Schwager experiments [9]. A plausible cause for this is the omission of the nonlinear repair rate. The exclusion of nonlinear repair rates in [241 was due to mathematical intractability of a nonlinear Kolmogorov equation. In this paper we are able to use an alternative, probabilistic approach to derive the cell survival probability that assumes both nonlinear lesion initiation and repair rates. A formal solution, though complicated, is presented. The paper is organized as follows. Section 2 contains a brief description of the radiation experiments. Section 3 presents the mathematical assumptions and the models. The derivation of the cell survival probability is given in Section 4. Concluding remarks are in Section 5. 2.
RADIATION
EXPERIMENT
To calculate cell survival probability, it is necessary to describe how cell survival is determined in a radiation experiment. The survival of a cell is not directly observable after irradiation as one cannot distinguish a living cell from a dead cell. A practical criterion for determining cell survival is that it must have the proliferative ability to form a colony of a given (observable) size within a specified time after irradiation. Thus to determine survivability, after irradiation, the cell has to be plated onto a petri dish filled with nutrients to allow for its proliferation. In some experiments there is an intermediate step, in which case the cells are held in a nongrowth environment (e.g., no nutrients) for a period of time before plating. This period is called the holding time. Experimental data have demonstrated that a substantial amount of cellular repair occurs during the holding time, which is therefore a key parameter in studying cell repair mechanisms. As such, a typical radiation experiment is composed of three consecutive time periods: the irradiation period (0, T), the holding period (T, T,), and the cloning period (7”, T,). Cells of
308
GRACE L. YANG
a given type are irradiated during (0, T) at a specified dose rate corresponding to a given total dose, D. Upon termination of the radiation, the cells are transferred to an environment that allows for repair of radiation-induced lesions but not for cell division. Some induced lesions are lethal and irreparable. The reparable lesions are referred to as “potentially lethal lesions” for the reason that incorrect repair may still result in cell death. The correct repair as well as the incorrect repair of lesions that causes cell transformation will result in cell survival. This will be discussed in Section 3 under Assumption A2 in the construction of stochastic models. At time T,, cells are plated for cell division. The clones are counted at time T,. The counting time T, is generally selected such that colonies have more than 50 cells and are observable by the naked eye. Each colony is taken as having originated from a single cell. The ratio of the number of colonies and initial number of irradiated cells, corrected for plating efficiency, measures a cell’s survival fraction for a given radiation dose. 3.
ASSUMPTIONS
AND
STOCHASTIC
MODELS
The following mathematical assumptions are made with respect to the occurrences of events in the experiment. Assumptions Al-A3 are the same as those given in [24]; we set the plating efficiency parameter 6 = 0 for simplicity. But Assumption A4 is different. In A4 a repair-misrepair model with a nonlinear repair rate K is introduced. For computing survival probability, alternative formulations different from A4 are needed; they are given in A5 and A6. Assumption Al. Poisson am’val of primary radiation particles. The primary radiation particles arrive at the cell nucleus according to a Poisson process with rate A(t) per unit time. A(t) is a positive continuous function of t. (It is known that the passage of particles in the cytoplasm has little effect on survival; see, e.g., Von Borstel and Rogers [22].) Then the expected number of primary particles arriving in [O, t 1is l,‘A(x) du. Assumption A2. Formation of lesions. Each primary generates a random number, M, of “spurs” with a probability-generating function (pgf) g(s) = Es”. For low-LET radiation, all spurs act independently of one another (Chatterjee et al. [41X Each spur has a probability 7~~ of generating a potentially lethal lesion and a probability rrz of generating an irreparable lesion that inactivates the cell (a lethal event). We have 0 < 7~~+ rTT2< 1, where 1 - rr, - nz is the probability that the spur has no effect on the cell. The values of 7~~ and z-Z may depend on the time at which the spur is created through the total accumulated dose received by the nucleus up to time t. Through this dependence a nonlin-
CELL SURVIVAL
PROBABILITY
UNDER RADIATION
309
ear rate of lesion initiation is introduced into the model. The effect of these M spurs is described by random variables U, and U,, where U1 is the number of potentially lethal lesions and U, is the number of lethal lesions induced by the M spurs. Under Assumption Al, the pgf of Ut and U, is given by E( I&@
) = E( 1- rr, - 7rz + %-IU, + rr2U*)M = g( 1- 9-r, - T2 + 7T,u1 + 7r*u2)
(1)
for lull G 1, luzI < 1.
Assumption A3. Colony-forming ability. We assume that a single lesion will inactivate the cell unless it is correctly repaired or is misrepaired to a mutant (transformed cell). We also assume that repair processes stop at the time of plating TP. Thus for a cell to divide and grow into a visible colony it is necessary that by time TP all potentially lethal lesions are either correctly repaired or misrepaired in transformation. Assumption A4. Repair - misrepair mechanisms. Radiation produces a variety of cellular lesions that immediately activate the cell’s repair mechanisms. Figure 1A illustrates a variety of radiation-induced lesions. Some of these lesions are reparable, and some are not. A stochastic model that describes the outcomes of lesion repair is as follows. The evolution of a single cell during and after radiation is a vector-valued continuous-time Markov process (X,Y, 2) = {(X,,Y,, Z,>; t 2 01, where X, is the number of potentially lethal lesions that the cell has at time t, Y, is the number of transformed lesions that the cell has at time t, and Z, is the number of lethal lesions (inactivation events) the cell has experienced up to time t. At the onset of the experiment, t = 0, no lesions exist; that is, X0 = Y,, = Z,, = 0. During the time interval (t, t + h], either nothing happens or one and only one of the following transitions takes place:
CX,,y,,Z,)
+a,
-l,Y,,Z,)
with probability
[ax,
+ .X,(X,
-
I)1 h + o(h) (correct
ar(X, - l,Y, + l,Z,>
with probability
-+(X,-l,Y,,Z,+l)
withprobabilityyX,h+o(h)
PX,h + o(h) (misrepair in transformation), (misrepair
-+ (Xt+u,,Yt,
Z,+u,)
with probability for uI,u2
repair),
hhP[U,=u,,
in inactivation),
U2=u21+dh),
= O,l;... (initiation
of new lesions),
310
GRACE L. YANG
where (Y, p, y, and K are nonnegative parameters that may depend on the time t and dose rate r. This is the model used in [24] except for the introduction of a nonlinear term K. Insofar as cell survival is concerned, it is not necessary to consider Y explicitly because, according to Assumption A3, the cell survival probability is equal to S = P[X, = 0, Z, = 01, which does not involve Y,. In other words, Y, does not contribute to cell death. This means that we can add the transition rates (Y+ p in computing the survival probability. In view of this, it is sufficient to consider the process (X, Z) = ((X,,Z,>; t a 0) where the transition of (X,,Z,) to (X, - l,Z,) is modified to (a + p>X, + KX,
W,,Z,) with probability
+(X,-1,-Z,>
(correct -+(X,-1,z,+11
repair
[ax,
+ KX,
and misrepair
with probability
yX,h + o(h)
with probability
hhP[U,
in transformation),
(misrepair -(X,
+ u,,ZI
+ u,)
for u,,uz
= ul, lJ, =
in inactivation),
u21= + o(h),
=O,l;... (initiation of new lesions),
where (Y, y, and K are nonnegative parameters that may depend on the time t and dose rate r. As explained in A4, (Y should be replaced by (Y+ p. As it is only a symbol for the parameter, we shall for simplicity continue using (Y even though it represents “(Y + p.” It is convenient to formulate the cell survival probability, S, in terms of the process (X,Z). But, in the presence of the nonlinear repair rate (K > 01, it is difficult to compute S from the (X, Z) process. We shall use a different probabilistic approach to compute S. This requires an alternative formulation of the process (X,, Z,>, t a 0, which is given in A6. It is easy to show that the Markov process (X,,Z,) describes the same repair process as that given in Assumption A6.
CELL SURVIVAL
PROBABILITY
UNDER
311
RADIATION
Assumption A6. Let N, be the Poisson number of primary particles that arrived at the nucleus in the time interval (O,tl, and let r;, for j=l ,.. . , N,, be their arrival times. At time rj, U[ potentially lethal lesions and Ui lethal lesions are created according to A2. Fix j. Assume that rj and U{ are given, for j = 1,. . . , N,. The cellular repair of the l.I/ lesions starts at time 5. The repair process is modeled by a two-dimensional Markov process {Rj(t - rj), y
Rj( t - TV)= R, The transition
and
qt
process
are
rates of the repair
(R,,U:)+(R,-LY+l)
-Tj)=w,.
with rate yR, with rate aR, + KR,( R, - 1).
-(R,-I,%)
Given N, and rj, for j = 1,. . . , N,, the repair ~~1,E$t - ~~1); t a T~) are stochastically independent.
processes
((Rj(t
-
In summary, under Assumptions Al-A4 (or with A4 replaced by A5 or A6), we can visualize several stochastic processes (not necessarily observable) that are evolving simultaneously or sequentially over time t. They are (1) the Poisson arrivals of the primary radiation particles, (2) the induction of lesions in the nucleus, (3) the Markovian repairmisrepair process during the irradiation and holding periods, and (4) the cell proliferations or cloning. These assumptions are clearly an oversimplification of the actual process. For instance, we have made no distinction of different kinds of radiation damage illustrated in Figure 1A. The reader is referred to [241 for further discussion and to a survey article by Sachs et al. in this volume for discussions of DNA damage and Markov models. 4.
COMPUTING Under
CELL
SURVIVAL
A3, the cell survival
probability
PROBABILITY is defined
s(D)=P[xTp=o,zT~=o]. Note that it is written through h(t) in Al.
as a function
by
(2)
of dose D. The dose is determined
312
GRACE L. YANG
Consider first the survival probability of a cell that was hit by one primary particle, say at time 0. Suppose that the induced number of potentially lethal lesions is ur and the number of lethal hits U, = 0. This initiates a repair process CR, W) with initial conditions R(O) = u1 and W(O) = 0. Under Assumption A6, the state space of CR, W) is a triangle with vertices (O,O), (u,,O), and (0, u,>. Let q,=
KT(T--l)+yr
CZr+
and
(3) for r = l,...,u,. Let p, be the waiting time of the process CR, W) in state (r, w>. Then p, is an exponential random variable with density
fW=s,exp(-q,xL
for x > 0 and Y > 0.
Note that the distribution of p,. does not depend on W. Thus the total cell repair time (correct repair and misrepair combined) is the sum EFr 1 p, of u, independent random variables p,. According to definition (21, the survival probability in the present case is S(D) = P[ R, = 0, W, = 01. It is the probability that the process (R, W) completes [he path {(u,,O) -+ (ul - 1,O) + (ul -2,0) + ... + (O,O)} by the time of plating, Tp, By the Markovian assumption A.5, P[(u,,O)
-+(u,
-130)
+(u*
-210)
+ ... + (0,O) by time T,lR(O) =P
Consequently,
g [ r=l
pr”Tp
= ul, W(0) = O]
&(l,u,).
(4)
I
the survival
probability
(5)
= where -
Remark. According created
= with A2,
=
U, to
initiation and the
0] and
the is
may accumulated
P(u,,O) on
in
derivative with time received
g(l-
-
= on and which spur the
CELL SURVIVAL
PROBABILITY
313
UNDER RADIATION
up to time t. In the present case, 7, is set to zero. This corresponds to the experimental situation of acute irradiation of cells for a very brief time period. We can set T = 0 and consider that lesions are produced instantaneously at time 0. In this case, the rr’s will be constant independent of time. For instance, rTTIstudied in [24] has the form 5 + r,[l exp( - at)], which reduces to 5 at t = 0. In the derivation of (5) we have separated the sample path from the time it requires to complete the path and thus avoided having to solve a differential equation. Formula (5) is valid for any transition rates q, and the repair rates, not necessarily quadratic functions of r. In the case of (31, q, have positive distinct values. The sum CF; 1 p, in the middle has density qlq2...qu,[$,,n
exp( -
91x) + -** + k,n exP( -
4u,x)l, (6)
+zl =(41-qk)...(qk~l-qk)(qk+l-qk)...(qu,-qk) 3
I
fork=1
,..., u,.
See, for example, [7]. If all the q, are equal, then CAL, p, has a gamma density. For other cases of q, the density of the sum can be obtained by inverting Laplace transforms. However, the formulas are messy. The repair process considered so far deals with potentially lethal lesions created under acute irradiation. This is similar to the repair-misrepair (RMR) model 1211. A stochastic version of RMR model is given by Albright [l, 21. The solution for the nonlinear case in finite time studied here is not given in Albright. We now to return to the general situation and investigate the repair processes generated by all the lesions (Uj’,Ui), for j = 1,2,. . . , created at times O
K4Tl. It is convenient to consider the cell’s repair activities separately in each time interval CT,,rj+ ,> between the arrivals of newly radiationinduced lesions. The repair begins at time r1 upon the arrival of ZJ, potentially lethal lesions. The repair is described by the process CR,, W,) with the initial condition R,(O) = u!, W,(O) = 0. Recall that a cell cannot survive any lethal hit or with any misrepaired lethal lesion. The cell survival probability at time t E CT,, TV) is
P[R,(t-71)
=x1,
W,(~-T~)=OIR~(~)=~~,W~(O)=~].
(7)
GRACE
314
L. YANG
At t = TV, there are new arrivals of lesions, Uf. Thus, the number potentially lethal lesions that the cell has at time r2 is
of
&(O)=&(r2-q)+q5
(8)
where the summands are independent. In general, the number of potentially at time 7j is
lethal
Rj(O)=Rj~,(7/-Tj_,)+U:”
lesions
j=1,2
that the cell has
>.‘.>
(9)
which together with II$(O) = 0 form the initial condition of the repair process (Rj(t - ~~1, Fi$t -TV)), for t E (7,,rj+ 1). Note that (Rj(t - ~~1, y(t - 7,)) has transition rates similar to that of (R, W) but with possibly different initial conditions, The probability distribution of the initial conditions Rj(0) and II$O) depend on the lesion initiation probabilities 7~~ and rTT2,which by Assumption A2 depend on TV. we are interested m the probability For each j and t E(7j,Tj+l), P[ Rj(t Denote
- TV) = rj, H$(t -TV) = OIlSi
the jth initiation
probability
P(z&O)
= rj_, + 241, y(O)
= 01.
by
= P[ U:’ = Ui, u: = 01.
By Assumption A2, this can be derived from (l), where the rr’s in this case depend on TV.We can immediately obtain the following conditional cell survival probability at time t > T: P[X(t)=O,Z(t)=OIN,=n,~
=c
all r,
,,...,
Tn]
n-1
i
,‘II, (+0)x
xp?,(t
p
-~)=rj,~(t-T,)=OIR,(0)=rj_,+u:,~(O)=O])
x P( u;t,O) x P x[R,(t-7,)=O,W,(t-7,)
=OlR,(O)
=rnpl +u;,W,(o)
=o]
) 1 (10)
where
the summation
extends
to all 0 < rj < rj_ 1 + u{, for j = 1,. 1s7n.
CELL SURVIVAL
PROBABILITY
UNDER RADIATION
The last term in (10) can be obtained from (3). It remains the middle term. We suppress j for simplicity. Let
Si-.u,
=
315
to compute
c
pl.
i=r
Denote (u,-l,O) = rl=
the event by B(r,u,) that the process moves from to (r,O). Using the fact that on the set B(r,u,), ,...,
Applying
(u,,Ol to [R(t -7)
Eq. (6) yields the result
P[S,,,,
- fig, [ 2 I-. ,=r
1
~i~~~‘{l-exP[-qi(f-7)l} 3
i=l
where
*F’ = iii(qj-qi) j=r jti
and P[B(r, u,)] = [(r, u,>, as given in (3). The unconditional survival probability is obtained by integrating the distributions of the T’S and the N(t) in (10). Formally,
S(D)=
c
Tl>O
mv)
=A
h(T)”
T... /0
out
7. /0
XP[X,=Z,=OIN,=~,T,,...,T,]
ifIh(~~)d~~
"'dT,,,
(11)
survival
prob-
i=l
where R(T) = /$A(y)dy. ability.
Setting
t = Tp yields the desired
316 5.
GRACE L. YANG
CONCLUDING
REMARKS
We have derived an explicit expression for the cell survival probability under ionizing radiation under the assumption of nonlinear rates of repair and lesion initiation. In fact, the nonlinear repair rate can be made quite arbitrary, not necessarily a quadratic function as is used in A5. This offers some flexibility in selecting appropriate repair rates to explain the cell repair mechanism. The result is an extension of what is obtained in [24]. Our formula makes explicit use of the holding time TP. TP is a key experimental parameter for studying the cell repair mechanism. It is, however, mostly neglected in the mathematical modeling literature. The survival probability derived here is a complicated analytical expression that requires heavy numerical calculation. It remains to be investigated whether the model will fit the experimental data. The author wishes to thank L. Le Cam for many helpful discussions. Thanks are also due to a referee for valuable suggestions that led to a substantial improuement of the paper. This paper was written while the author was visiting MSRI in Berkeley. Research at MSRI was supported in part by National Science Foundation grant DMS 8505550. REFERENCES N. W. Albright, A Markov model of the repair-misrepair process of cell survival, Lawrence Berkeley Laboratory Rep. LBL-11660, Univ. California, 1981. N. W. Albright, A Markov formulation of the repair-misrepair model of cell survival, &z&t. Res. 118: l-20 (1989). K. H. Chadwick, and H. P. Leenhouts, The Molecular Theory of Radiation Biology, Springer-Verlag, New York, 1981, p. 102. A. Chatterjee, P. Koehl, and J. L. Magee, The theoretical consideration of the chemical pathways for radiation-induced strand breaks, A&. Space Res., 6( I1):97-105 (1986). S. B. Curtis, Lethal and potentially lethal lesions induced by radiation-a unified repair model, Radiut. Res. 106:252-270 (1986). M. M. Elkind, Repair processes in radiation biology, Radiut. Res. 100:425-449 (1984). W. Feller, An Introduction to Probability Theory and Its Applications, Vol. 2, Wiley, New York, 1966. D. Frankenberg, M. Frankenberg-Schwager, and R. Harbich, Interpretation of the shapes of survival curves in terms of induction and repair-misrepair of DNA double strand breaks, Br. J. Cancer 49cSuppl. 6):223-238 (1984). M. Frankenberg-Schwager, D. Frankenberg, D. Blocher, and C. Adamczyk, Repair of DNA double strand breaks in irradiated yeast cells under new growth conditions, Rudiat. Res. 82:498-510 (1980).
CELL SURVIVAL 10
11 12 13 14 15
16 17
18 19 20 21
22 23
24
PROBABILITY
UNDER
RADIATION
317
M. Frankenberg-Schwager, D. Frankenberg, and R. Harbich, Exponential or shouldered survival curves result from repair of DNA double strand breaks depending on postirradiation conditions, Radiat. Res. 11454-63 (1988). G. Iliakis and W. Pohlit, Quantitative aspects of repair of potentially lethal damage in mammalian cells, Int. J. Radiat. Biol. 36(6):649-658 (1979). A. M. Kellerer and H. H. Rossi, A generalized formulation of dual radiation action, Radiat. Res. 75:471-488 (1978). L. Le Cam, Stochastic models of lesion induction and repair in yeast, Math. Biosci., this volume. L. Le Cam and J. Neyman, Eds., Probability Models and Cancer, North-Holland, New York, 1982. N. F. Metting, L. A. Braby, W. C. Roesch, and J. M. Nelson, Dose rate evidence for two kinds of radiation damage in stationary-phase mammalian cells, Radiat. Res. 103:204-218 (1985). J. Neyman and P. S. Puri, A hypothetical stochastic mechanism of radiation effects in single cells, hoc. Roy. Sot. Land. B 213:139-160 (1981). F. Q. H. Ngo, E. A. Blakely, C. A. Tobias, P. Y. Chang, and L. Lommel, Sequential exposures of mammalian cells to low- and high-LET radiations. II. As a function of cell-cycle stages, Radiut. Res. 115:54-69 (1998). R. K. Sachs and L. R. Hlatky, Stochastic dose-rates in radiation cell survival models, Radiat. Emiron. Biophys. 29:169-184 (1990). R. K. Sachs, P.-L. Chen, P. Hahnfeldt, and L. R. Hlatky, DNA damage caused by ionizing radiation, Math. Biosci. (1992). C. A. Tobias, The repair-misrepair model in radiobiology: comparison to other models, Radiat. Res. 104:S-77-S-95 (1985). C. A. Tobias, E. A. Blakely, F. Q. H. Ngo, and T. C. H. Yang, The repair-misrepair model of cell survival radiation biology and cancer research, A. Meyn and R. Withers, Eds., Raven Press, New York, 1980, pp. 195230. R. C. Von Borstel and R. W. Rogers, Alpha-particle bombardment for the Habrobrazon egg. II. Response of the cytoplasm, Radiat. Res. 8:248-253 (1958). G. L. Yang and C. E. Swenberg, Stochastic models for cells exposed to ionizing radiation, in Modeling of Biomedical Systems, J. Eisenberg and M. Witten, Eds., Elsevier, Amsterdam, 1986, pp. 85-89. G. L. Yang and C. E. Swenberg, Stochastic modeling of dose-response for single cells in radiation experiments, Math. Sci. 16:16-65 (1991).