On the stochastic theory of compartments: Solution for n-compartment systems with irreversible, time-dependent transition probabilities

B U L L E T I N OF I~ATHEMATICAL BIOLOGY

VOLUME 36, 1974

ON THE STOCHASTIC THEORY OF COMPARTMENTS: SOLUTION FOR n-COMPARTMENT SYSTEMS WITH IRREVERSIBLE, TIME-DEPENDENT TRANSITION PROBABILITIES

[ ] M. C~mDENAS and J. tI. MATIS Institute of Statistics, Texas A & M University, College Station, Texas 77840

The bivariate distribution of a two-compartment stochastic system with irreversible, time-dependent transition probabilities is obtained for any point in time. The mean and variance of the number of particles in any compartment and the covariance between the number of particles in each of the two compartments are exhibited and compared to existing results. The two-compartment system is then generalized to an n-compartment catenary and to an n-compartment mammillary system. The multivariate distributions of these two systems are obtained under two sets of initial conditions: (1) the initial dintribution is known; and (2) the number of particles in each compartment of the system at time t = 0 is determined. The moments of these distributions are also produced and compared with existing results. I. Introduction. Several recent p a p e r s h a v e i n t r o d u c e d a stochastic t h e o r y for c o m p a r t m e n t a l models b a s e d on a finite p o p u l a t i o n t h e o r y . Matis a n d H a r t l e y (1971) solved t h e c o m p l e t e m u l t i v a r i a t e d i s t r i b u t i o n for a n n - c o m p a r t m e n t stochastic m o d e l w i t h c o n s t a n t flow i n t e n s i t y coefficients a n d T h a k u r et al. (1973) were able to v e r i f y t h e m a r g i n a l distributions w i t h o u t resorting to t h e c o m p l e t e c u m u l a n t g e n e r a t i n g function. I n o t h e r r e c e n t work, T h a k u r et al. (1972) e x a m i n e d a o n e - c o m p a r t m e n t m o d e l w i t h a t i m e - d e p e n d e n t i n p u t process b y using g e n e r a t i n g functions, a n d P u r d u e (1974) solved t h e o n e - c o m p a r t m e n t m o d e l w i t h t i m e - d e p e n d e n t i n p u t a n d o u t p u t processes w i t h o u t a p p e a l to t h e g e n e r a t i n g f u n c t i o n techniques. 489

490

M. CARDENAS AND J. H. MATIS

Other directly related work has been reported as stochastic migration processes. Chiang (1968) derives the solution to a general n-dimensional process with constant migration and death rates but time-dependent immigration rates to the compartments. The single dimensional time-dependent birth-deathimmigration process has also been widely studied (see e.g. Adke, 1969; Bailey, 1964). I n addition, Purl (I968) investigated a special time-dependent two compartment birth-death process without immigration and gives an infinite series representation of the solution. The multiple compartment model which in addition to time-dependent immigration incorporates time-dependent migration and death rates, however, remains yet to be solved. Such a model is the natural generalization of the previous work above. With the abundant current research establishing the widespread existence of biological rhythms in nature from the cellular to the population level (see e.g. Sollberger, 1965) such time-dependent models become highly relevant. The present paper first obtains the complete bivariate distribution for a two compartment model where the transition coefficients are time-dependent but irreversible and where the probability distribution of particles at t = 0 is known. The solution is then generalized to other models with irreversible transitions, and the complete explicit solution is given for the general n-compartment, time-dependent eatenary and mammillary systems.

II. Two-Compartment Model. 1. Definitions and Hypothesis. Let the stochastic variables Xl(t) and X2(t) denote the number of particles at time t in compartments 1 and 2 respectively as X,o(t}

Xzo[t}

~,ol(t)

),oz(t)

Figure I. Two-compartment mode] with irreversible, t i m e - d e p e n d e n t transition probabilities

illustrated in Figure 1. Also let AXi(t ) = X~(t + At) - Xi(t), for i = 1, 2. The transition probabilities are defined as follows: Prob. {a single unit moves from compartment 1 to compartment 2 in the interval (t, t + At)} = Xl(t)221(t ) At + o(At);

ON T H E STOCHASTIC T H E O R Y O F C O M P A R T M E N T S

Prob. {a unit (t, t + Prob. {a unit (t, t +

491

in compartment i leaves the system in the interval At)} -- X,(t)Aoi(t)At + o(At), for i = 1, 2; enters the system through compartment i in the interval At)} = am(t ) At + o(at), for i -- 1, 2.

Clearly the probability of two or more units moving from one compartment to another (here the exterior of the system is considered as a third compartment) in At is o(At). Therefore there are only five events with probability of first order magnitude of At. The probabilities of these five events are:

Prob. {Axe(t) = kl, Axa(t) = ka I x~(t), xdt)} { al0(t)At + o(At) aol(t)X~(t)At + o(At) = a2,(t)Xl(t)At + o(At) a2o(t)At + o(At) }to2(t)X2(t)At + o(At)

fork1 fork1 for kz for/cl for/c 1

= = = = =

1, k 2 = - 1 , k~ - 1, kg, O,/ca = O,/ca =

O; = O; = 1; 1; --1.

There are several methods in the literature which one could employ to derive a differential equation for a generating function. I n particular, one can use the "random variable" technique (see Bailey, 1964) to obtain the differential equation ~M(01, Oa, t) _ [)tot(t)(e_e, _ 1) + }~2~(t)(e-e,+e= - 1)] ~M(8~, 02, t) St ~O~

+ [ao2(t)(e-0~ - 1)] a M ( 0 . 02, t)

a0a

+ [al0(t)(e el - ~) +

;~2o(t)(e<

-

1)]M(01, 02, t)

(1)

for the joint moment generating function, M(01, 02, t), of X l ( t ) and Xa(t ). The differential equation for the cumulant generating function may now be obtained by dividing (1) by M(O 1, Oa, t). This results in ~K(01, 02, t) = [Ao~(t)(e_Ol _ 1) + }~a~(t)(e-e~ +°~ - 1)] ~K(01, 02, t) St 801 + [Ao2(t)(e-°~

-

1)]

~K ( OI, 02, t) 802

+ [}~lo(t)(e°* -- 1) + A~o(t)(e°= -- 1)].

(2)

2. The Joint Cumulant Generating Function. The linear partial differential equation (2) m a y be solved with characteristic theory (see Ford, 1955). A sketch of the steps of the solution for this simple model is given in order that it

492

M. C A R D E N A S A N D J. I t . M A T I S

m a y serve as a guide in t h e later n - c o m p a r t m e n t models. the t h e o r y one m u s t first solve the subsidiary equations,

I n aeeordance w i t h

dt dO1 dO2 -i- = ,~o~(t)(1 - e-°~) + ~21(t)(1 - e-°~+°2) = Ao2(t)(1 - e-°:) dK

(3)

= 2lo(t)(e°l -- 1) + i2o(t)(e°2 - 1)" L e t t i n g V~ = e°l - 1 a n d V2 = e°2 - l, (3) can be rewritten as dV1 = [2o~(t)V1 + 22~(t)(V1 -

V2)] dt

d V~ = )~o2(t)V2 dt d K = [~lo(t)Vx + )12o(t)V2] dt. The system of differential equations (4) can be solved sequentially. solution is

(4) The

V2 = C2 exp [ f ~ 2o2(Z) dz]

Now substituting back the value of V1 a n d V2 into (5) and solving sequentially for the arbitrary constants C 1, C2 and C 3 one obtains

C l = u l ( O l , 02, t ) = ( e ° l -

1) e x p { - ; [ ~ o l ( Z ) + ) ~ 2 1 ( z ) ] d z }

+ (eO2-1) f: z2dt ) exp {- ft o2(z) dz - ; ~ [;~ol(Z) + ~2~(z)] dz} dt~

(6)

ON THE STOCHASTIC THEORY

OF COMPARTMENTS

493

C 8 = u3(el, e2, t) = K(01, 02, t) - (e°l - 1) f : 21o(t2)

× exp(-ft:

21(z)]dz}dt2 2

[~01(~) "~-

- (e°2 - 1)(fl ;~2o(t2) exp [- ftt2 ;~o2(z) dz] dtl

+ flft'; 1o(t2); 21(tl)exp{-f'l; o2(z)dz - f t t : [2~01(z)+ ;~21(z)] dz} dt 1 dt2)According to the theory, the three independent integrals can now be equated by the functional relationship u3(01, 02, t) = 9~[u~(01, 02, t), u2(01, 02, t)]

(7)

and the function ~o is identified from the initial conditions. If the joint cumulant generating function of XI(0) and X~,(0) is k n o w n to be k(01, 02), t h e n evaluating (7) at t = 0 yields

~(o~, 02) = ~[(e°l - 1), (e°~ -- 1)].

(s)

Letting Yl = (e°l - 1) and y~. = (e°~ - 1), (8) becomes /ella (y~ + 1), In (Y2 + 1)] = ~(yl, Y2).

(9)

I n view of (9) equation (7) can be written as

%(0~, 02, t) = k[ln {u~(O~, 02, t) + 1}, In {u2(0~, 02, t) + 1}].

(10)

Finally substituting (6) into (10) and rearranging terms one obtains K(0~, 02, t) = k[ln {1 + (e°~ - l)p1~(t) + (e°2 - 1)p12(t)}, In {1 + (e°~ - 1)p~,2(t)}] + (e°l - 1)a~(t) + (e°~ - 1)32(t) where pll(t) = exp

-

[;~ol(z) + 221(z)] dz}

(11)

494

M. CARDENAS AND J. H. MATIS

lal2(t)----; ;~21(tl) exp {- ;t:~o2(z)dz - ;i [~01(z)+ ~21(z)]dz} dtl

x exp

(; -

Ao2(Z)dz 1

f)

[Io1(Z) + 121(z)] dz

2

;

dtl dt2.

This solution m a y be verified b y s u b s t i t u t i o n into (2). 3. M o m e n t s of the Distribution. The m o m e n t s of the distribution can be f o u n d b y differentiating ( 1 1 ) w i t h respect to the a p p r o p r i a t e parameter(s) a n d setting 01 = 02 = 0. If/~1,/t2, a n , a22, ~12 are the means, variances, a n d covariance of XI(0 ) a n d X2(0 ) respectively, t h e n the m o m e n t s a t a n y time t are given b y

E[XI(t)]= th1911(t)+ 31(t) E[X2(t)] = t~2p:2(t) + t~p~2(t) + ~2(t) Vat [x~(t)]

= (~

Var

=

[X2(t)]

-

((Y22 -

t~)[p~(t)] ~ + t~p~l(t) + .~(t) 1~2)[Ps2(t)] 2 + 2c~12P22(t)P~2(t)

+ ((rii -- /~i)[Pi2(t)] 2 + /~2/O2z(t) -{- /tiP12(t ) + ~2(t) Coy [X~(t), X2(t)] = (r~2pil(t)p22(t) + (aii - th)P~2(t)Pii(t) •

(12)

I n the special case where X i ( 0 ) = X1 a n d X2(0) = X2, i.e. the initial values are known, the m o m e n t s at time t are given b y E[Xi(t)] = X l P i i ( t ) + ~i(t) E[X2(t)] = X2P22(t) + Xlp~2(t) + 32(t) Var [Xl(t)] = X ~ p n ( t ) [ 1 - p n ( t ) ]

+ 3~(t)

Var [X2(t)] = Xlp~2(t)[1 - P12(t)] + X2p22(t)[1 - p22(t)] + 39.(t) Coy [Xl(t), X2(t)] = - X l p 1 2 ( t ) p n ( t ) .

(13)

This agrees with t h e previous results in Matis (1970) in the special ease of constant t r a n s i t i o n rates. T h e m o m e n t s for t h e case of one c o m p a r t m e n t can be

ON T H E S T O C H A S T I C T H E O R Y O F C O M P A R T M E N T S

495

derived from (12) b y setting all variables with a subscript of 2 equal to zero. In this case E[XI(t)] = /~1 exp [ - ;

~01(z) dz ]

+ j: 1o(t2)exp Var [Xl(t)] = (c~11 - /xl) exp

Aol(Z)dz] at2

-2

~ol(z) dz

(14)

0

This is the previous result reported b y Purdue (1974) for the one-compartment system.

I I I . n-Compartment Catenary System. 1. Definitions and Hypothesis. The two compartment model will now be generalized to a complete, irreversible, time-dependent catenary system illustrated in Figure 2. The number of particles in the ith compartment of the ),io (t)

kzo(t)

),no[t)

X~z(t)

koi(t)

.Xoz(t }

,kn, n_i(t )

kon(t)

Figure 2. n-Compartment catenary model with irreversible, time-dependent transition probabilities system is the stochastic variable X~(t), and let AXi(t) = X~(t + At) - X~(t) for i = 1. . . . , n. Then the transition probabilities are defined by: Prob. {A single unit moves from compartment i to compartment i + 1 in the interval (t, t + At)} = 2~+l.~(t)X~(t) At + o(At), for i -- 1. . . . . n-l; Prob. (A unit enters the system through compartment i in the interval (t, t + At)} = ,~o(t) At + o(At), for i = 1 , . . . , n;

496

IV[. CARDENAS AND J. H. MATIS

P r o b . (A unit in c o m p a r t m e n t i leaves t h e s y s t e m in the i n t e r v a l (t, t + At)} = ~o~(t)Xi(t) At + o(At), for i = 1 , . . . , n. T h e p r o b a b i l i t y of t w o or m o r e t r a n s i t i o n s in At is clearly o(At). T h e r e f o r e t h e r e are 3n - 1 e v e n t s w i t h p r o b a b i l i t y of first order m a g n i t u d e of At. T h e s e e v e n t s h a v e probabilities as follows: P r o b . {AXe(t) = k~ . . . . .

:{

AXe(t) = k n ] X~(t) . . . . . X~(t)}

)~0(t)At + o(At) for k, = 1 a n d all o t h e r k's are zero, i = 1 . . . . , n; Ao~(t)X,(t)At + o(At) for k~ = - 1 a n d all o t h e r k's are zero, i = 1. . . . . n; ~, ~_ l(t)X,_ l(t)At + o(At) for k, = 1 a n d k~_ ~ = - 1 a n d all o t h e r k's are zero, i = 2 . . . . . n.

( l i v e n t h a t t h e particles in t h e s y s t e m act i n d e p e n d e n t l y , t h e differential e q u a t i o n for t h e joint e u m u l a n t g e n e r a t i n g f u n c t i o n is hence f o u n d to be

~K(O~,..., 0~, t) ~t = ~.~ {Ao,(t)(e-°~- 1) + A,+l,~(t)(e-°~+e~+~ - 1)} i=1

~K(O~ . . . . .

0~, t)

~0~

+ ~

~=1

h~o(t)(e°~- 1),

(15)

w h e r e An+ l.n(t) -- O.

2. The Joint Cumulant Generating Function and its Moments.

The procedure T h e necessary s u b s i d i a r y

u s e d in solving (2) will a g a i n be used in solving (15). e q u a t i o n s for solving (15) are dt 1

d0i Ao,(t)(1 - e°0 + hf+l.,(t)(1 - e-°,+°,+l)

dK(O 1. . . . . On, t) ~ = ~ ~jo(t)(e °, S u b s t i t u t i n g V, = (e°* -

?to~(t)Vi+ ~+1 i(t)(V~

~K ~t

i=1

~

for i

1,..

n.

(16)

1) for i = 1 . . . . . n in (16) one o b t a i n s

~V, ~t

1)

Ajo(t) V i

V,+I),

for i = 1 , . .

n

(17)

T h e s y s t e m of differential e q u a t i o n s (17) can be solved s e q u e n t i a l l y (with t h e use o f m a t h e m a t i c a l induction). T h e solution is a s y s t e m of e q u a t i o n s in t h e

ON T H E STOCHASTIC T H E O R Y O F COMPARTMENTS

497

n + 1 constants of integration. One c a n t h e n solve this s y s t e m of e q u a t i o n s s e q u e n t i a l l y (again w i t h t h e use o f m a t h e m a t i c a l i n d u c t i o n ) f o r t h e a r b i t r a r y constants, obtaining

C i = u~(01,..., On, t) = (e °, - 1) e x p {-gt(t, 0)}

+ (e°,+~ - 1) f~ k~+l,~(tn-~) exp {-g,(tn_,,

O)

g~+l(t, tn-i)} dtn-~

-

+

(e o, - 1)

..

l=/+2 x exp

f

n-~" - g ~ ( t n _ i, O) -

F , ( t n - l + l . . . . , t~_~) n-t+2

l-1 ~_~ gin(in_m, t n _ m + l ) m=~+l

gz(t, tn_z+l)}dtn_~+l . . . . ,dtn_ ~ f o r i = 1 . . . . , n -

-

Cn-1 = Un-l(O1 . . . . .

2

On, t) = (e°~-~ -- 1) e x p { - - g n _ l ( t , 0)}

+ (e°~ -- 1) f~ kn.n_~(t~) e x p {--g~_~(t~, O) -- g~(t, t~)} dt~ Cn = u ~ ( O ~ , . . . , 0,, t) = (e°- -

1) e x p {-gn(t, 0)}

C,+1 = un+l(01 . . . . , On, t) = K(OI . . . . , On, t) -

~

1) f~ ~jo(tn_j+l) e x p { - g j ( t ,

(e°J-

tn_]+l)}dtn_j+l

J=l

-

(e°J+~ -

1)

,~jo(tn_j+l),~j+l,i(tn_j)

J=l

n-I+1

x e x p { - g j ( t n _ j, tn_j+l) - gj+l(t, tn_~.)} dtn_ J d Q _ j + 1 +

(e °, -- 1) ]=I

l=j+2

... n-]+l

× Flj(tn-l+ 1. . . . .

n-]

,~jO(tn_j+l) n-t+2

tn-¢) l-1

× eXlO{--gAtn_j, tn_j+~) -

~

gm(tn-m,t~-m+~)

m=]+l -

gz(t, tn-z+~)} dtn_~+~. • .dtn_j+~

where

g~(x, y) =

K

[2o~(z) + k~+l,~(z)] dz

for i = 1,..., n

(18)

498

M. C A R D E N A S A N D J . H. M A T I S

and I~j(z~ . . . .

, z~_ j _ ~) = ;~, ,_ ~ ( z ~ ) . . . ~. + ~, ~.(z~_ ~._ ~).

T h e n the general integral of (15) is ¢[u1(01 . . . . . 0n, t) . . . . , un(O1. . . . , On, t)] = u~+1(01 . . . . , 0~, t) and the function ¢ is identified from the initial conditions. dition is K(O~ ....

, On, o) = ]c(O~ . . . .

(19)

If the initial con-

(20)

, On),

t h e n using (19) and (18) one can verify t h a t the joint c u m u l a n t generating function of X~(t) . . . . , X=(t) at any time t is given by

I=1

. . . ,

In [1 + (e°o - 1)p=(t)]} + ~ (eo,_ 1)~(t) (21) i=1

where

p~(t) = exp {-g~(t, 0)}, for i = 1 . . . . , n, p~.,+l(t) =

)~+l,~(tn-~) exp {-g~(tn_~, O) - g~+l(t, t~_i) } dtn_~, 0

Fj~(tn-j+l . . . . . tn-~)

loi. j(t) . . . . n-~ × exp

f

n-t+2

--g~(tn-t, O) --

i--1

~

gm(tn_m, t n _ m + l )

m=i+l -

gj(t, t~_j+l)} dtn-s+l-- .dtn_~, for/ = 1,...,n-

2andj_>

5~(t) = f : 21o(tn) exp {-gl(t, tn)} dtn, ~2(t) = f j A2o(tn_l)exp {-g2(t, tn_l)dtn_l + f j ~tt. ?ho(tn)~2~(tn_l) × exp {-g~(tn_l, t~) - g2(t, t~_~)} dt~_l dtn,

i + 2

ON T H E STOCHASTIC T H E O R Y OF COMPARTMENTS

499

3i(t) = f{ Z~o(t~-~+l) exp {-g~(t, t,-~+l)} dt,-~+l

+

;f

~ - 1, o(tn-~+2)~.~- l(t~-i+ 1)

× exp {-g~_l(tn_~+l, tn_~+e) - g~(t, t~_~+l)} dt~_~+~ dtn_~+ 2 +

""" j=l

x exp

~jo(t~- j + 1)Pij'(t~- ~+ 1. . . . . t~ _ j) n--I+2

n-J+l

(

-gj.(t~_j, tn_.~+l ) --

~

~m(tn_m,

t n _ m + l ) -- ~'i(t, t n _ i + l )

m=]+l

x dt~_~+l...dt~_~.+l,

}

for i = 3 . . . . , n.

The moments of the distribution are found b y differentiating (21) and are given b y

E[X~(t)] = ~ t~lPu(t) + ~,(t) for i = 1 , . . . , n, l=l

Var [X,(t)] = k=x ~ { ~z=~a~Pu(t)P~(t) + i~P~(t)[1- p~(t)]} + ~(t)

f o r i = 1. . . . , n ,

J

rain U, ]]

f o r i , j = 1. . . . . n

(22)

where t~, ~, are the m e a n a n d variance of X,(0) for i = 1. . . . , n and a,~.is the covariance of X,(0) andX~(0) for i , j = 1 , . . . , n a n d i # j. Note t h a t previous equations (12) are a special case of this.

3. Special Case: Initial Number of Particles in Each Compartment is Known. I f the initial n u m b e r of particles in each c o m p a r t m e n t is known, say X~(0) = X~ for i - 1. . . . . n, t h e n (20) becomes

K(O,,. . ., ~n, O) = ~ X~Oi

(23)

i=1

and the joint c u m u l a n t generating flmction at a n y time t is found to be

K(O 1. . . . ,O,,t) = ~ X ~ l n [ 1 + ~ (e° , - 1)pkl(t)] + ~ (e°~--1)3~(t) k=l

l=I¢

~=1

(24)

500

M. C A R D E N A S A N D J. H. MATIS

and its m o m e n t s are i

E[X~(t)] = ~ Xzpli(t ) + 3~(t) f o r i = 1. . . . . n, l=1 i

Xzpli(t)[1 - pz,(t)] + 3~(t) for i = 1. . . . . n,

Var [X~(t)] = ~ 1=1

rain [i, j]

Cov[Xi(t), Xj(t)] = -

Xzpli(t)plj(t ) for i , j = 1. . . . . n a n d i ¢ j.

~

(25)

1=1

I t has been p o i n t e d out in Marls (1970) t h a t equations of the form (24) i m p l y t h a t t h e r a n d o m v e c t o r IX1 . . . . . Xk] at a n y time t is d i s t r i b u t e d as a sum of multinomials and Poissons. N o t e t h a t t h e same distribution t h e o r y holds in this general t i m e - d e p e n d e n t case.

I V . n-Compartment Mammillary System. 1. Definition and Hypothesis. T h e two c o m p a r t m e n t model will now be generalized to a complete, irreversible, t i m e - d e p e n d e n t m a m m i l l a r y system, illustrated in Figure 3. The n u m b e r of particles in each of the i c o m p a r t m e n t s

;~z°(t)---~I

f I

X I(t)

kol(t}

IX2,Ct' IX3,{t' X20(t)

X:~o(t)

Iknl(t)

•

• Xon(t}~ kno(t]

Figure 3. n-Compartment mammillary model with irreversible, time-dependent transition probabilities of the s y s t e m at time t is the stochastic variable X~(t) for i = 1 . . . . . n, with T h e transition probabilities are:

Xl(t ) representing the central c o m p a r t m e n t .

Prob. {a single u n i t moves from c o m p a r t m e n t 1 to c o m p a r t m e n t i in the interval (t, t + At)} = Xl(t)~l(t)At + o(At), for i = 2 , . . . , n; Prob. {a unit enters t h e system t h r o u g h c o m p a r t m e n t i in the interval (t, t + At)} = ~o(t)At + o(At), for i = 1 , . . . , n; Prob. (a unit in c o m p a r t m e n t i leaves t h e s y s t e m in t h e interval (t, t + 'At)} = X~(t)2o~(t)At + o(At), for i = 1 . . . . . n.

ON THE

STOCHASTIC

THEOI~Y

OF

COMPARTMENTS

There are 3n - 1 events w i t h probability of first order m a g n i t u d e of At. probabilities of these events are

501

The

Prob. {AXI(t ) = k~ . . . . , AXn(t ) = k n ] X l ( t ) . . . . , Xn(t)} A~0(t)At + o(At) for k~ = 1 and all other k's are zero, i = 1 , . . . , n; Ao,(t)X,(t)At + o(At) for k, = - 1 and all other k's are zero, i = 1. . . . . n; A,l(t)Xl(t)At + o(At) for k 1 = - 1 a n d k, = 1 and all other k's are zero, i = 2 . . . . , n. Given t h a t the particles in the system act independently, the differential equation for the joint e u m u l a n t generating function is found to be

~K(O~ . . . . . ~t

On, t)

= [2ol(t)(e-°~ - 1) + ~ 2n(t)(e-°~+5 - 1)] j=2

eK(O~,. . ., 0,, t) 901

X

+ ~ toj(t)(e-5 - 1) ~K(01 . . . . . 0~, t) s=2 ~Oj. n

+

a

o(tl(e o, -

(26)

1).

i=1

2. The Joint C u m u l a n t Generating Punction and Its Moments. is solved b y the same procedure employed in solving (2). I f the generating function of XI(0), X 2 ( O ) , . . . , Xn(O) is k n o w n and 02. . . . . 0~), t h e n the joint cumulant generating function at

E q u a t i o n (26) joint cumulant given b y k(01, a n y time t is

i=1

In [1 + (e°2 - 1)p22(t)] . . . . . In [1 + (e°,~ - 1)pn~(t)] ~ n

+

(e 0, -

1) j(t)

(27)

i=1

where

Pll(t) = e x p f - f ~ p~(t) =

[~ol(Z)+ j~2 ~j~(z)] dz}

exvI-f£

2oi(z) dz}

fori=

2. . . . , n

502

M. C A R D E N A S A N D J. H . M A T I S

plj=

~ kjl(tl) e x p { - ; 1

[)~01(z) + ~

k~l(z)] dz - ~ttl 2oi(z) d z } d t l forj = 2,...,n

forj = 2,...,n. The moments of the distribution at any time t are obtained b y differentiating the eumulant generating function and are given b y

E[X~(t)] = I~zpzz(t) + 3~(t) E[Xdt)] =/~zlol~(t) + /hp~(t) + ~(t)

for i = 2 , . . . ,

n

Var [X~(t)] = ((r~ - i~)p~(t) + tz~p~(t) + S~(t) Var [X~(t)] = (au - I~,)P~(t) + Iz,P**(t) + as(t)

+ 2%,pu(t)pl~(t ) + ((r~ - izl)p~(t) +/~pl,(t), f o r i = 2 . . . . ,n, Coy [X~(t), X,(t)] = (q~p~,(t)pz~(t) - i~p~t(t)p~(t) + ~,p,(t)p~(t) for i = 2 , . . . , n Coy [X~(t), Xj(t)] = a~jp~(t)p~(t) + ~r~p~(t)pl~(t ) + aljpz(t)pl~(t )

+ a~p~(t)p~(t) - t~p~(t)p~j(t) for i , j = 2 , . . . , n a n d i

#j,

where t~,au, a~j fori, j = 1 , . . . , n ,

i # j,

are the moments of the distribution at time t = 0. I t can easily be verified t h a t the m o m e n t s derived here are the same as the m o m e n t s for the t w o - c o m p a r t m e n t system (12) if one considers only the central c o m p a r t m e n t (as the first c o m p a r t m e n t in the t w o - c o m p a r t m e n t system) and any one of the n - 1 peripheral compartments. s

3. Special Case: Initial Number of Particles in Each Compartment is Known. Suppose t h a t the n u m b e r of particles in each c o m p a r t m e n t at time t = 0 is

ON T H E STOCHASTIC T H E O R Y OF COMPARTMENTS

known to be X~(0) = X~, for i = 1. . . . . n. generating flmction is given b y

2(01,...,

503

E q u i v a l e n t l y the initial cumulant

=

(2s)

x,e,. t=l

I n this case the e u m u l a n t generating function at a n y time t is K(01 . . . . , O~,t) = X l l n [ 1

+ ~

(e ° ~ - 1)iolj.(t)]

j=l

+

X, ln [1 + (e°, - 1)p,~(t)] + ~=2

(e°, - 1)~j(t)

(29)

"=

a n d the moments are E[Xl(t)] = X l p l l ( t ) + 81(t) E[X~(t)] = X~p~(t) + X~p~(t) + 3~(t) for i = 2 . . . . . n Var [Xl(t)] = Xzp~l(t)[1 - pll(t)] + 31(t ) Var [X,(t)] = X,p,,(t)[1 - pu(t)] + X~p~,(t)[1 - p~,(t)] + 3,(t) f o r i = 2. . . . . n Coy [Xdt), Xj(t)] = - X l p l , ( t ) p l j ( t )

for i, j = 1 , . . . , n a n d i ¢ j .

This is a generalization of a time-independent result given in Matis et al. (1974). V. Discussion. Explicit solutions have been given for the n - c o m p a r t m e n t catenary and m a m m i l l a r y systems with irreversible rates. One m i g h t readily conceive several extensions to the foregoing c o m p a r t m e n t a l systems. One extension would be to consider the most general ease of an n-compartm e n t stochastic system w i t h irreversible, time-dependent transition probabilitie_s. This describes a system in which a particle in the first c o m p a r t m e n t can move into a n y one of the remaining n - 1 c o m p a r t m e n t s and a particle in the second c o m p a r t m e n t can move into a n y one of the remaining n - 2 compartments, and so on. This system was not considered because the equations are tedious and because the authors could not conceive a n y practical applications for it. A realistic and useful extension would be to consider reversible transition probabilities. The general solution is a rather ambitious undertaking, b u t particular cases seem fruitful.

504

M. CARDENAS AND J. H. MATIS

T h i s r e s e a r c h w a s c a r r i e d o u t i n p a r t i a l f u l f i l l m e n t o f t h e r e q u i r e m e n t s for t h e degree o f D o c t o r o f P h i l o s o p h y , for M. C a r d e n a s . T h e a u t h o r s are i n d e b t e d t o t h e r e v i e w e r s for t h e i r c o m m e n t s a n d s o m e h e l p f u l s u g g e s t i o n s o n n o t a t i o n .

LITERATURE Adke, S. R. 1969. "A Birth, Death and Migration Process." J. Appl. Prob., 6, 687-691. Bailey, N. T. J. 1964. The Elements of Stochastic Processes with Applications to the Natural Sciences. New York: Wiley. Chiang, C . L . 1968. Introduction to Stochastic Processes in Biostatistics. New York: Wiley. Ford, L. I~. 1955. Differential Equations. New York: McGraw-Hilh Matis, J. I-I. 1970. "Stochastic Compartmental Analysis: Model and Least Squares Estimation from Time Series Data." Dissertation, Texas A & M University. - a n d H. O. g a r t l e y . 1971. "Stochastic Compartmental Analysis: Model and Least Squares Estimation from Time Series Data." Biometrics, 27, 77-102. , M. Cardenas and 1~. L. Kodell. 1974. " O n the Probability of Reaching a Threshold in a Stochastic MammiUary System." Bull. Math. Biology, 36, 445-454. Purdue, 1) . 1974. "Stochastic Theory of Compartments." Bull. Math. Biology, 36, 305-309. Puri, P. S. 1968. "Interconnected Birth and Death Processes." J . Appl. Prob., 5, 334-349. Sollberger, A. 1965. Biological Rhythm Research. Amsterdam: Elsevier. Thakur, A. K., A. l~escigno and D. E. Schafer. 1972. "On the Stochastic Theory of Compartments: I. A Single Compartment System." Bull. Math. Biophysics, 84, 53-65. , and - - - 1973. "On the Stochastic Theory of Compartments: II. Multi-Compartment Systems." Bull. Math. Biology, 35, 263-271. I~¢C~rVED 2-7-74