BirthDeath Processes, Lattice Path Combinatorics,
Continued Fractions, and Orthogonal Polynomials
Fabrice Guillemin
France Telecom, CNET, Lannuon, Breizh
Algorithms Seminar
February 2, 1998
[summary by Philippe Flajolet and Fabrice Guillemin]
A properly typeset version of this document is available in
postscript and in
pdf.
If some fonts do not look right on your screen, this might be
fixed by configuring your browser (see the documentation here).
Abstract
Classic works of KarlinMcGregor and JonesMagnus
have established
a fully general correspondence between birthdeath processes
and continued fractions of the StieltjesJacobi type
together with their associated orthogonal polynomials.
This fundamental correspondence can be revisited
in the light of the otherwise known combinatorial correspondence
between weighted lattice paths and continued fractions.
For birthdeath processes,
this approach separates clearly the formal apparatus
from the analyticprobabilistic machinery and neatly delineates those
parameters that are amenable to
a treatment by means of continued fractions and orthogonal polynomials.
1 BirthDeath Processes
Consider a particle initially in state 0
that, at any given time, may change to another state 1 (where it stays),
with rate l.
This means that the probability of a state change in
an interval of time of length dt is l dt.
Then, the probability p_{0}(t) that the particle is still in state 0 at
time t satisfies
p_{0}(t+dt)p_{0}(t)=l p_{0}(t) dt
or p_{0}'(t)=l p_{0}(t), whose solution is an exponential distribution,
Similarly, a particle initially in state 0 that may change either to
state 1 with rate l or to state 1 with rate µ
will satisfy (p_{j}(t) is the probability of being in state j at time t)
p_{0}(t)=e 

,
p_{1}(t)= 

(1e 

),
p_{1}(t)= 

(1e 

).

The interpretation is obvious: the particle stays in state 0 for a random
amount of time with an exponential distribution of rate l+µ and
then changes to states 1,+1
with probabilities equal to to l/(l+µ)
and µ(l+µ).
In a general birthdeath process a particle can be in any state
in {0,1,2,...} and when in state j, it can only change to state
j+1 at rate l_{j} or to state j1 at rate µ_{j}.
By analogy with the model of
an evolving population (whose size is represented by the
state), the l_{j} are called
birth rates and the µ_{j} death rates.
The general problem is to understand the evolution of a process
given values (or properties) of its birth and death rates;
see [12, Ch. 4] for an excellent introduction.
Let p_{n}(t) be the probability of being in state n at time t.
An essential rôle is played by the coefficients
p_{n} = 
l_{0}l_{1}··· l_{n1} 

µ_{1}µ_{2}··· µ_{n} 

.

Indeed, a classical result asserts that the process is ergodic
(the expected time to return from each state to itself is finite)
if and only if
(The first condition ensures the existence of an invariant measure for the
embedded discretetime Markov chain; the second one
guarantees that, in the continuoustime process,
the particle is not absorbed at infinity in
finite time.)
In that case, one has
where these quantities represent the long run probability of being in
state n.
More puzzling is the nonstationary behaviour of the process
that is described by the infinitedimensional
differential system
p'
_{j}(
t)=
l
_{j1}p_{j1}(
t)(
l_{j}+µ
_{j})
p_{j}(
t)+µ
_{j+1}p_{j+1}(
t),
p_{j}(0)=
d_{j,0}.
(1)
Although finitedimensional versions are ``easy'' and reduce to combinations of
exponentials, it is precisely the infinitedimensional character of
the system that renders its analysis interesting.
In a series of important papers, Karlin and McGregor [10, 11]
have
developed a general connection between
the fundamental system (1)
and an associated
family of orthogonal polynomials.
Later, Jones and Magnus
constructed a direct
continued fraction representation; see [8, 9].
This summary is an account of Guillemin's lecture (see [5, 6]),
as well as of later developments.
The point of view that is adopted here
consists in relating the combinatorial theory of lattice paths to
birthdeath processes in the following way:
(i) trajectories of birthdeath processes are precisely
lattice paths; (ii) lattice paths have generating functions
expressed as continued fractions; (iii) the Laplace transform
expresses the main parameters of birthdeath processes as weighted
lattice paths to which the combinatorial theory applies.
2 Lattice Paths and Continued Fractions
It is known that the formal theory of continued fraction expansions
for power series is identical to the combinatorial theory of
weighted lattice paths; see [1, 2, 4].
Define a path u=(U_{0},U_{1},...,U_{n}) to be a sequence of points in
the lattice N×N such that if U_{j}=(x_{j},y_{j}), then
x_{j}=j and y_{j+1}y_{j}=1. If successive points are
connected by edges, then an edge can only be
an ascent (a: y_{j+1}y_{j}=+1),
a descent (b: y_{j+1}y_{j}=1),
or a level step (c: y_{j+1}y_{j}=0).
Thus a path is always nonnegative
and by a horizontal translation, one may always assume that
x_{0}=0. A path can be encoded by a word
with a,b,c representing the three types of steps.
What we call the standard encoding is
such a word in which each step a,b,c is subscripted by the
value of the ycoordinate of its associated point.
For instance,
w=a_{0}a_{1}a_{2}b_{3}c_{2}c_{2}a_{2}b_{3}b_{2}b_{1}a_{0}c_{1}
encodes a path that connects the source U_{0}=(0,0) to the destination
U_{12}=(12,1). We freely identify a path u defined
as a sequence of points,
its word encoding w, and the corresponding monomial.
We consider various geometric conditions that may be imposed on paths:
H_{k,l} is the collection of all paths that connect
a source at altitude k to a destination at altitude l,
H^{[£ h]} denotes paths of height (maximal altitude)
at most h, etc.
Theorem 1
The collection H_{0,0} of all paths has generating function
Proof.
It suffices to observe that (1f)^{1}=1+f+f^{2}+···
generates symbolically all the sequences with components f.
For instance, in H_{0,0}, the expressions
generate successively
paths composed from c_{0} level steps only,
paths of height at most 1 without c_{1} steps, all paths
of height at most 1. The complete continued fraction representation is
easily built by stages in a similar fashion.
In particular, the collection of all paths from
level 0 to level 0 with height at most h is
a rational
fraction, whose numerators and denominators,
P_{h},Q_{h}, each satisfy the recurrence
y_{h+1}=(1c_{h})y_{h}a_{h1}b_{h}y_{h1},
with Q_{1}=P_{0}=0, Q_{0}=P_{1}=1.
(Linear fractional transformations are 2×2 matrices in disguise!)
Wellknown path decompositions, like those based
on first or last time
at which levels are reached, can then be used provided they are
combinatorially ``unambiguous''. This and simple manipulations
on linear fractional transformations give access
to many geometric constraints in addition to (2) and (3).
We cite here some representative identities from [1, 2],
H_{0,h1}^{[<h]}= 
a_{0}a_{1}··· a_{h1} 

Q_{h} 

,
H_{0,k}= 

( 
Q_{k}H_{0,0}P_{k} 
) 
,
(4) 
H_{k,l}= 
Q_{k} 

a_{0}··· a_{k1}b_{1}··· b_{l} 


( 
Q_{l}H_{0,0}P_{l} 
) 
,
(5) 
where the latter holds provided k£ l.
The forms (2), (3) (4), (5) can be
converted into bona fide counting generating functions of
paths weighted multiplicatively by means of the
combinatorial morphism,
c(a_{k})=a_{k}z,
c(b_{k})=b_{k} z,
c(c_{k})=g_{k}z.
In that case,
the continued fraction (2) becomes
the general fraction of the Jtype (for Jacobi);
see [7, 9, 13].
3 The Connection
We illustrate here
in its simplest form the manyfaceted connection between birthdeath
processes and continued fractions.
It was apparently first stated explicitly by Jones and Magnus but it is
implicit in earlier works of Karlin and McGregor.
The connection goes through
the probabilities p_{i,j}(t) of being in state j
at time t starting from state i and
the Laplace transforms,
P_{i,j}(s) = 
ó õ 

p_{i,j}(t) e^{st} dt.

Theorem 2
The Laplace transform of the probability of
return to the origin satisfies
We offer here two proofs. A third proof that is based on
``uniformization of time'' can also be given but is omitted in this note.
Proof.[Proof 1]Take the Laplace transform of the fundamental system (1)
(so that p_{j}(t)=p_{0,j}(t))
and use the induced relations on the ratios P_{0,r}/P_{0,r+1}.
This proof is the most direct but the least
illuminating from a structural standpoint.
In particular, this proof does not provide an immediate grasp
on the question of deciding which parameters are amenable
to continued fraction representations.
Proof.[Proof 2]
Examine the times at which the (continuous time) birthdeath
process {L_{t}} changes states.
This defines an embedded (discrete time) Markov chain {Y_{n}}.
Then the set of
trajectories of the chain {Y_{n}} is exactly the family of lattice
paths of Section 2. The method
consists in splitting the probabilities by conditioning according to
all legal trajectories.

The first observation is that,
given a lattice path w=w_{1}w_{2}··· w_{n},
the probability p_{0,0}(t w) of being back to 0 at time t
having followed the path w is
Pr{L_{t}=0 w}
=Pr 
ì í î 
S 

+S 

+... +S 

£ t,
S 

+S 

+... +S 

+S 

> t 
ü ý þ 
,

where S_{q}_{j} is the random variable that represents the sojourn time at
the state q_{j} determined by w_{1}··· w_{j},
while the righthand side involves q_{n+1} that ranges over all
legal ``continuations'' of w (in the case of
H_{0,0}, one has w_{n+1}=a_{0} and q_{n+1}=0).
As seen already, the sojourn time at some state e
is exponential with parameter (l_{e}+µ_{e}) so that its
Laplace transform is (l_{e}+µ_{e})/(s+µ_{e}+l_{e}).

The second observation is that the probability of a path in the
embedded chain
is the product of the individual transition probabilities, namely
l_{j}/(l_{j}+µ_{j}) and µ_{j}/(l_{j}+µ_{j}).
The different sojourn times are independent by the nature of the process
(the strong Markov property satisfied by {L_{t}}).
Also, sums of independent random variables
correspond to products of Laplace transforms.
Thus, the Laplace transform of the probability in the continuous model of
following a path w has a product form;
for instance, to w=a_{0}a_{1}b_{2}a_{1}, there corresponds the transform

æ ç ç è 









ö ÷ ÷ ø 
·

æ ç ç è 

l_{0}+µ_{0} 

s+l_{0}+µ_{0} 


l_{1}+µ_{1} 

s+l_{1}+µ_{1} 


l_{2}+µ_{2} 

s+l_{2}+µ_{2} 


l_{1}+µ_{1} 

s+l_{1}+µ_{1} 


ö ÷ ÷ ø 
.

Thus, the Laplace transform P_{0,0}(s) is, apart
from a fudge factor of 1/(s+l_{0}), a
sum over all paths lattice from zero to zero weighted multiplicatively
by the probabilistic morphism,
c'(a_{j})= 

,
c'(b_{j})= 

,
(6) 
with c'(c_{j})=0.
In other words, one has
P_{0,0}(s)=c'(H_{0,0})1/s+l_{0},
and the statement follows.
The same method applies to the computation of transition probabilities,
the analysis of maximum height, and so on. For instance,
the probability of reaching state k has
P_{0,k}(s)= 


( 
A_{k}(s)P_{0,0}(s)B_{k}(s)

) 
,

where A_{k}/B_{k} is the kth convergent of the continued fraction that
represents P_{0,0}, so that A_{k},B_{k} are simple variants of
c'(P_{k}),c'(Q_{k}).
Orthogonality.
In the case of paths,
the reciprocals of the Q_{h} polynomials,
Q_{h}(z)=z^{h}c(Q)(z^{1}) are
formally orthogonal with
respect to a measure defined its moments,
L 
[z^{n}] º 
ó õ z^{n} dµ (z) = H_{0,0,n}.
(7) 
Formal aspects of paths and orthogonality are detailed in
Godsil's book [3].
A similar orthogonality property then holds for the probabilistic
counterparts A_{h},B_{h} of the P_{k},Q_{k} polynomials.
This provides alternative expressions
of various probabilistic quantities
in terms of scalar products involving the measure µ of (7).
One can rederive in this way, via the combinatorial theory,
a number of formulæ
originally discovered by Karlin and McGregor. For instance,
one has
p_{m,n}(t)=p_{n} 
ó õ 

e^{tx}
q_{m}(x)q_{n}(x) dµ(x),

where the q_{k} polynomials
(closely related to the B_{k} and Q_{k})
satisfy the recurrence l_{n}q_{n+1}
+(xl_{n}µ_{n})q_{n} +µ_{n} q_{n1}=0.
4 So What?
The original motivation for the talk comes from
the need to elucidate the behaviour of certain queueing systems
in the context of telecommunication applications.
For instance, the single server queue (M/M/1)
is modelled by l_{j}=r,
µ_{j}=1, while the infinite server queue (M/M/¥)
corresponds to l_{j}=r, µ_{j}=j.
(Models of population growth lead to considering different types
of weights, like l_{j}=(j+1)r, µ_{j}=j.)
More specifically, the problem is to quantify parameters
of some simple statistical multiplexing scheme that
describe
the quality of service on an ATM link. The
relevant model is that of the
M/M/¥ queue and parameters are to be analysed,
like the duration q of an excursion
above some level c,
the volume V of lost information, or the number of bursts C in
a busy period.
Each parameter leads to a specific continued fraction
representation.
By Theorem 2, the basic continued fraction
of the M/M/¥ process is
This is recognizable as an instance of Gauß's
continued fraction associated to a quotient
of contiguous hypergeometric functions.
The numerator and denominator polynomials
are the PoissonCharlier polynomials that are
orthogonal with respect to the Poisson measure.
The quantity V (area) leads to
challenging asymptotics questions
both for the M/M/¥ queue
and for the M/M/1 queue.
A simple modification of the basic techniques
of this note shows that the bivariate
Laplace transform with
(s,u) ``marking'' (t,V) is obtained by the modified
morphism,
c''(a_{j})= 

,
c''(b_{j})= 

.

In the case of area under the M/M/1 queue,
quotients of continuous Bessel functions make an
appearance.
Stripped of its probabilistic context,
the
corresponding
problem of tail estimation then admits a purely analytic formulation:
Problem. Let A(x) be a function whose Laplace transform is
with J_{n} a Bessel function, and r>0 a parameter.
Show that, for some constants c_{1},c_{2}, one has
ó õ 

A(y) dy 
~ c_{1} x^{1/4}e 

,
(x®+¥).

Under plausible analytic
or probabilistic conjectures, precise (and useful!)
quantitative conclusions can be drawn.
See the papers by Guillemin and Pinchon [5, 6] for
full developments.
References
 [1]

Flajolet (P.). 
Combinatorial aspects of continued fractions. Discrete
Mathematics, vol. 32, 1980, pp. 125161.
 [2]

Flajolet (P.), Françon (J.), and Vuillemin (J.). 
Sequence of operations analysis for dynamic data structures. Journal of Algorithms, vol. 1, 1980, pp. 111141.
 [3]

Godsil (C. D.). 
Algebraic Combinatorics. 
Chapman and Hall, 1993.
 [4]

Goulden (Ian P.) and Jackson (David M.). 
Combinatorial Enumeration. 
John Wiley, New York, 1983.
 [5]

Guillemin (Fabrice) and Pinchon (Didier). 
Continued fraction analysis of the duration of an excursion in an
M/M/¥ system. Journal of Applied Probability, vol. 35,
n°1, 1998, pp. 165183.
 [6]

Guillemin (Fabrice) and Pinchon (Didier). 
Excursions of birth and death processes, orthogonal polynomials, and
continued fractions. 
Preprint, 1998. To appear in Journal of Applied
Probability. 21 pages.
 [7]

Henrici (Peter). 
Applied and Computational Complex Analysis. 
John Wiley, New York, 1977. 3 volumes.
 [8]

Jones (William B.) and Magnus (Arne). 
Application of Stieltjes fractions to birthdeath processes. In
Saff (E. B.) and Varga (Richard S.) (editors), Padé and rational
approximation. pp. 173179. 
New York, 1977. Proceedings of an International
Symposium held at the University of South Florida, Tampa, Fla., December
1517, 1976.
 [9]

Jones (William B.) and Thron (W. J.). 
Continued Fractions: Analytic Theory and Applications. 
AddisonWesley, 1990, Encyclopedia of
Mathematics and its Applications, vol. 11.
 [10]

Karlin (S.) and McGregor (J. L.). 
The differential equations of birthanddeath processes, and the
Stieltjes moment problem. Transactions of the American Mathematical
Society, vol. 85, 1957, pp. 489546.
 [11]

Karlin (Samuel) and McGregor (James). 
The classification of birth and death processes. Transactions of
the American Mathematical Society, vol. 86, 1957,
pp. 366400.
 [12]

Karlin (Samuel) and Taylor (Howard). 
A First Course in Stochastic Processes. 
Academic Press, 1975, second edition.
 [13]

Perron (Oskar). 
Die Lehre von der Kettenbrüchen. 
Teubner, 1954, vol. 2.
This document was translated from L^{A}T_{E}X by
H^{E}V^{E}A.