P1: JzG
WY047-Holman 0470022159pre September 18, 2004 2:45
ii
A First Course
in Stochastic Models
Henk C. Tijms
Vrije Universiteit, Amsterdam, The Netherlands
Copyright
c
2003 John Wiley & Sons Ltd, The Atrium, Southern Gate, Chichester,
West Sussex PO19 8SQ, England
Te lephone (+44) 1243 779777
Email (for orders and customer service enquiries):
Visit our Home Page on www.wileyeurope.com or www.wiley.com
All Rights Reserved. No part of this publication may be reproduced, stored in a retrieval system or
transmitted in any form or by any means, electronic, mechanical, photocopying, recording, scanning
or otherwise, except under the terms of the Copyright, Designs and Patents Act 1988 or under the
terms of a licence issued by the Copyright Licensing Agency Ltd, 90 Tottenham Court Road, London
W1T 4LP, UK, without the permission in writing of the Publisher. Requests to the Publisher should
be addressed to the Permissions Department, John Wiley & Sons Ltd, The Atrium, Southern Gate,
Chichester, West Sussex PO19 8SQ, England, or emailed to , or faxed to (+44)
1243 770620.
This publication is designed to provide accurate and authoritative information in regard to the subject
matter covered. It is sold on the understanding that the Publisher is not engaged in rendering
professional services. If professional advice or other expert assistance is required, the services of a
competent professional should be sought.
Other Wiley Editorial Offices
John Wiley & Sons Inc., 111 River Street, Hoboken, NJ 07030, USA
Jossey-Bass, 989 Market Street, San Francisco, CA 94103-1741, USA
Wiley-VCH Verlag GmbH, Boschstr. 12, D-69469 Weinheim, Germany
Contents
Preface ix
1 The Poisson Process and Related Processes 1
1.0 Introduction 1
1.1 The Poisson Process 1
1.1.1 The Memoryless Property 2
1.1.2 Merging and Splitting of Poisson Processes 6
1.1.3 The M/G/∞ Queue 9
1.1.4 The Poisson Process and the Uniform Distribution 15
1.2 Compound Poisson Processes 18
1.3 Non-Stationary Poisson Processes 22
1.4 Markov Modulated Batch Poisson Processes 24
Exercises 28
Bibliographic Notes 32
References 32
2 Renewal-Reward Processes 33
2.0 Introduction 33
2.1 Renewal Theory 34
2.1.1 The Renewal Function 35
2.1.2 The Excess Variable 37
2.2 Renewal-Reward Processes 39
2.3 The Formula of Little 50
2.4 Poisson Arrivals See Time Averages 53
2.5 The Pollaczek–Khintchine Formula 58
2.6 A Controlled Queue with Removable Server 66
2.7 An Up- And Downcrossing Technique 69
Exercises 71
Bibliographic Notes 78
References 78
3 Discrete-Time Markov Chains 81
4.5.3 First Passage Time Probabilities 170
4.6 Transient Distribution of Cumulative Rewards 172
4.6.1 Transient Distribution of Cumulative Sojourn Times 173
4.6.2 Transient Reward Distribution for the General Case 176
Exercises 179
Bibliographic Notes 185
References 185
5 Markov Chains and Queues 187
5.0 Introduction 187
5.1 The Erlang Delay Model 187
5.1.1 The M/M/1 Queue 188
5.1.2 The M/M/c Queue 190
5.1.3 The Output Process and Time Reversibility 192
5.2 Loss Models 194
5.2.1 The Erlang Loss Model 194
5.2.2 The Engset Model 196
5.3 Service-System Design 198
5.4 Insensitivity 202
5.4.1 A Closed Two-node Network with Blocking 203
5.4.2 The M/G/1 Queue with Processor Sharing 208
5.5 A Phase Method 209
CONTENTS vii
5.6 Queueing Networks 214
5.6.1 Open Network Model 215
5.6.2 Closed Network Model 219
Exercises 224
Bibliographic Notes 230
References 231
6 Discrete-Time Markov Decision Processes 233
6.0 Introduction 233
References 338
9 Algorithmic Analysis of Queueing Models 339
9.0 Introduction 339
9.1 Basic Concepts 341
viii CONTENTS
9.2 The M/G/1 Queue 345
9.2.1 The State Probabilities 346
9.2.2 The Waiting-Time Probabilities 349
9.2.3 Busy Period Analysis 353
9.2.4 Work in System 358
9.3 The M
X
/G/1 Queue 360
9.3.1 The State Probabilities 361
9.3.2 The Waiting-Time Probabilities 363
9.4 M/G/1 Queues with Bounded Waiting Times 366
9.4.1 The Finite-Buffer M/G/1 Queue 366
9.4.2 An M/G/1 Queue with Impatient Customers 369
9.5 The GI /G/1 Queue 371
9.5.1 Generalized Erlangian Services 371
9.5.2 Coxian-2 Services 372
9.5.3 The GI /P h/1 Queue 373
9.5.4 The Ph/G/1 Queue 374
9.5.5 Two-moment Approximations 375
9.6 Multi-Server Queues with Poisson Input 377
9.6.1 The M/D/c Queue 378
9.6.2 The M/G
/
c Queue 384
9.6.3 The M
conviction that theory is better understood when the algorithms that solve the
problems the theory addresses are presented at the same time. This textbook tries
to recognize what the computer can do without letting the theory be dominated
by the computational tools. In some ways, the book is a successor of my earlier
book Stochastic Modeling and Analysis. However, the set-up of the present text is
completely different. The theory has a more central place and provides a framework
in which the applications fit. Without a solid basis in theory, no applications can be
solved. The book is intended as a first introduction to stochastic models for senior
undergraduate students in computer science, engineering, statistics and operations
research, among others. Readers of this book are assumed to be familiar with the
elementary theory of probability.
I am grateful to my academic colleagues Richard Boucherie, Avi M andelbaum,
Rein Nobel and Rien van Veldhuizen for their helpful comments, and to my stu-
dents Gaya Branderhorst, Ton Dieker, Borus Jungbacker and Sanne Zwart for their
detailed checking of substantial sections of the manuscript. Julian Rampelmann
and Gloria Wirz-Wagenaar were helpful in transcribing my handwritten notes into
a nice Latex manuscript.
Finally, users of the book can find supporting educational software for Markov
chains and queues on my website />CHAPTER 1
The Poisson Process and
Related Processes
1.0 INTRODUCTION
The Poisson process is a counting process that counts the number of occurrences
of some specific event through time. Examples include the arrivals of customers
at a counter, the occurrences of earthquakes in a certain region, the occurrences
of breakdowns in an electricity generator, etc. The Poisson process is a natural
modelling tool in numerous applied probability problems. It not only models many
real-world phenomena, but the process allows for tractable mathematical analysis
as well.
The Poisson process is discussed in detail in Section 1.1. Basic properties are
k
,n= 1, 2, .
A First Course in Stochastic Models H.C. Tijms
c
2003 John Wiley & Sons, Ltd. ISBNs: 0-471-49880-7 (HB); 0-471-49881-5 (PB)
2 THE POISSON PROCESS AND RELATED PROCESSES
Then S
n
is the epoch at which the nth event occurs. For each t ≥ 0, define the
random variable N(t) by
N(t) = the largest integer n ≥ 0 for which S
n
≤ t.
The random variable N(t) represents the number of events up to time t.
Definition 1.1.1 The counting process {N(t), t ≥ 0} is called a Poisson process
with rate λ if the interoccurrence times X
1
,X
2
, have a common exponential
distribution function
P {X
n
≤ x}=1 − e
−λx
,x≥ 0.
The assumption of exponentially distributed interoccurrence times seems to be
restrictive, but it appears that the Poisson process is an excellent model for many
real-world phenomena. The explanation lies in the following deep result that is
only roughly stated; see Khintchine (1969) for the precise rationale for the Poisson
j=0
e
−λt
(λt)
j
j!
,t≥ 0. (1.1.1)
The Erlang (k, λ) distribution has the probability density λ
k
t
k−1
e
−λt
/(k − 1)!.
Theorem 1.1.1 For any t>0,
P {N(t) = k}=e
−λt
(λt)
k
k!
,k= 0, 1, . (1.1.2)
That is, N(t) is Poisson distributed with mean λt.
Proof The proof is based on the simple but useful observation that the number
of arrivals up to time t is k or more if and only if the kth arrival occurs before or
at time t. Hence
P {N(t) ≥ k}=P {S
k
≤ t}
= 1 −
k−1
The memoryless property of the Poisson process
Next we discuss the memoryless property that is characteristic for the Poisson
process. For any t ≥ 0, define the random variable γ
t
as
γ
t
= the waiting time from epoch t until the next arrival.
The following theorem is of utmost importance.
Theorem 1.1.2 For any t ≥ 0, the random variable γ
t
has the same exponential
distribution with mean 1/λ. That is,
P {γ
t
≤ x}=1 − e
−λx
,x≥ 0, (1.1.3)
independently of t.
4 THE POISSON PROCESS AND RELATED PROCESSES
Proof Fix t ≥ 0. The event {γ
t
>x} occurs only if one of the mutually exclusive
events {X
1
>t+x}, {X
1
≤ t, X
1
+X
n
≤ t, S
n+1
>t+x}=
t
0
P {S
n+1
>t+x | S
n
= y}λ
n
y
n−1
(n − 1)!
e
−λy
dy
=
t
0
P {X
n+1
>t+ x − y}λ
n
y
n−1
(n − 1)!
0
e
−λ(t+x−y)
λdy
= e
−λ(t+x)
+ e
−λ(t+x)
(e
λt
− 1) = e
−λx
,
proving the desired result. The interchange of the sum and the integral in the second
equality is justified by the non-negativity of the terms involved.
The theorem states that at each point in time the waiting time until the next arrival
has the same exponential distribution as the original interarrival time, regardless
of how long ago the last arrival occurred. The Poisson process is the only renewal
process having this memoryless property. How much time is elapsed since the last
arrival gives no information about how long to wait until the next arrival. This
remarkable property does not hold for general arrival processes (e.g. consider the
case of constant interarrival times). The lack of memory of the Poisson process
explains the mathematical tractability of the process. In specific applications the
analysis does not require a state variable keeping track of the time elapsed since the
last arrival. The memoryless property of the Poisson process is of course closely
related to the lack of memory of the exponential distribution.
Theorem 1.1.1 states that the number of arrivals in the time interval (0,s) is
Poisson distributed with mean λs. More generally, the number of arrivals in any
time interval of length s has a Poisson distribution with mean λs.Thatis,
P {N(u + s) − N(u) = k}=e
−10/3
(10/3)
k
k!
= 0.3528.
The answer to question (b) follows from the memoryless property stated in Theo-
rem 1.1.2 and is given by
P {γ
5
> 5}=e
−5/3
= 0.1889.
In view of the lack of memory of the Poisson process, it will be intuitively clear
that the Poisson process has the following properties:
(A) Independent increments: the numbers of arrivals occurring in disjoint intervals
of time are independent.
(B) Stationary increments: the number of arrivals occurring in a given time interval
depends only on the length of the interval.
A formal proof of these properties will not be given here; see Exercise 1.8. To
give the infinitesimal-transition rate representation of the Poisson process, we use
1 − e
−h
= h −
h
2
2!
+
h
3
3!
2
(t), t ≥ 0} are indepen-
dent Poisson processes with respective rates λ
1
and λ
2
, where the process {N
i
(t)}
corresponds to type i arrivals. Let N(t) = N
1
(t) + N
2
(t), t ≥ 0. Then the merged
process {N(t), t ≥ 0} is a Poisson process with rate λ = λ
1
+ λ
2
. Denoting by Z
k
the interarrival time between the (k − 1)th and kth arrival in the merged process
and letting I
k
= i if the kth arrival in the merged process is a type i arrival, then
for any k = 1, 2, ,
P {I
k
= i | Z
k
= t}=
(a) It will be obvious that the process {N(t)} satisfies the properties (A) and (B).
To verify property (C) note that
P {one arrival in (t, t +t]}
=
2
i=1
P
one arrival of type i andnoarrival
of the other type in (t, t + t]
= [λ
1
t + o(t)][1 − λ
2
t + o(t)]
+ [λ
2
t + o(t)][1 − λ
1
t + o(t)]
= (λ
1
+ λ
2
)t + o(t) as t → 0.
Property (D) follows by noting that
P {no arrival in (t, t +t]}=[1 − λ
1
∞
t
P {Y
2
>Y
1
>t| Y
1
= x}λ
1
e
−λ
1
x
dx
=
∞
t
e
−λ
2
x
λ
1
e
−λ
1
x
dx =
+ λ
2
)t]. Hence
P {I
k
= 1,Z
k
>t}=P {I
k
= 1}P {Z
k
>t},
showing that P {I
k
= 1 | Z
k
= t}=λ
1
/(λ
1
+ λ
2
) independently of t.
(b) Obviously, the process {N
i
(t)} satisfies the properties (A), (B) and (D). To
verify property (C), note that
P {one arrival of type i in (t, t + t]}=(λt)p
i
+ o(t)
k + m
k
p
k
1
p
m
2
e
−λt
(λt)
k+m
(k + m)!
= e
−λp
1
t
(λp
1
t)
k
k!
e
−λp
2
t
(λp
2
1
+ λ
2
and the probability that a newly arriving parker is a long-term
parker equals λ
1
/(λ
1
+ λ
2
).
Example 1.1.2 A stock problem with substitutable products
A store has a leftover stock of Q
1
units of product 1 and Q
2
units of product 2.
Both products are taken out of production. Customers asking for product 1 arrive
according to a Poisson process with rate λ
1
. Independently of this process, cus-
tomers asking for product 2 arrive according to a Poisson process with rate λ
2
.
Each customer asks for one unit of the concerning product. The two products serve
as substitute for each other, that is, a customer asking for a product that is sold
out is satisfied with the other product when still in stock. What is the probability
distribution of the time until both products are sold out? What is the probability
that product 1 is sold out before product 2?
To answer the first question, observe that both products are sold out as soon as
−1
k=0
Q
1
+ Q
2
− 1
k
λ
2
λ
1
+ λ
2
k
λ
1
λ
1
+ λ
2
Q
1
+Q
inventory system with back ordering. In this model customers asking
for a certain product arrive according to a Poisson process with rate λ. Each cus-
tomer asks for one unit of the product. The initial on-hand inventory is S. Each
time a customer demand occurs, a replenishment order is placed for exactly one
unit of the product. A customer demand that occurs when the on-hand inventory
is zero also triggers a replenishment order and the demand is back ordered until
a unit becomes available to satisfy the demand. The lead times of the replenish-
ment orders are independent random variables each having the same probability
distribution with mean τ . Some reflections show that this
(
S − 1,S
)
inventory sys-
tem can be translated into the M/G/∞ queueing model: identify the outstanding
replenishment orders with customers in service and identify the lead times of the
replenishment orders with the service times. Thus the limiting distribution of the
number of outstanding replenishment orders is a Poisson distribution with mean
λτ . In particular,
the long-run average on-hand inventory =
S
k=0
(
S − k
)
e
−λτ
(λτ )
k
k!
k
D
k
. By the above argument, the number of type k
customers present at time t is Poisson distributed with mean (λp
k
)D
k
. Thus, by the
independence property of the split Poisson process, the total number of customers
present at time t has a Poisson distribution with mean
s
k=1
λp
k
D
k
= λµ.
This proves (1.1.6) for the case that the service time has a discrete distribution
with finite support. Any service-time distribution can be arbitrarily closely approx-
imated by a discrete distribution with finite support. This makes plausible that the
insensitivity result (1.1.6) holds for any service-time distribution.
Rigorous derivation
The differential equation approach can be used to give a rigorous proof of (1.1.6).
Assuming that there are no customers present at epoch 0, define for any t>0
p
j
(t) = P {there are j busy servers at time t},j= 0, 1, .
Consider now p
(t + t), dividing by t and letting t → 0, we find
p
0
(t) =−λ(1 −B(t))p
0
(t)
p
j
(t) =−λ(1 −B(t))p
j
(t) + λ(1 − B(t))p
j−1
(t), j = 1, 2, .
Next, by induction on j, it is readily verified that
p
j
(t) = e
−λ
t
0
(1−B(x)) d x
λ
t
0
(1 − B(x))d x
the vehicles to the regions in such a way that all regions provide, as nearly as
possible, a uniform level of service to the customers. The service level in a region
is measured as the long-run fraction of time that all vehicles are occupied (it will
be seen in Section 2.4 that the long-run fraction of delayed customer orders is also
given by this service measure).
Let us assume that the parameters are such that each region gets a large number
of vehicles and most of the time is able to directly provide a vehicle for an arriving
customer order. Then the M/G/∞ model can be used as an approximate model
to obtain a satisfactory solution. Let the dimensionless quantity R
i
denote
R
i
= λ
i
E(S
i
), i = 1, ,F,
12 THE POISSON PROCESS AND RELATED PROCESSES
that is, R
i
is the average amount of work that is offered per time unit in region i.
Denoting by c
i
the number of vehicles to be assigned to region i,wetakec
i
of
the form
c
i
k=c
i
e
−R
i
R
k
i
k!
≈ 1 −
c
i
− R
i
√
R
i
,i= 1, ,F,
where (x) is the standard normal distribution function. By requiring that
c
1
− R
1
√
R
i
F
i=1
R
i
.
This value of k is the guideline for determining the allocation (c
1
, ,c
F
) so that
each region, as nearly as possible, provides a uniform service level. To illustrate
this, consider the numerical data:
c = 250,F = 5,λ
1
= 5,λ
2
= 10,λ
3
= 10,λ
4
= 50,λ
5
= 37.5,
E(S
1
) = 2,E(S
2
≈ 63.05 and
c
5
≈ 90.98. This suggests the allocation
(c
∗
1
,c
∗
2
,c
∗
3
,c
∗
4
,c
∗
5
) = (16, 34, 46, 63, 91).
Note that in determining this allocation we have used the distributions of the
processing times only through their first moments. The actual value of the long-run
fraction of time during which all vehicles are occupied in region i depends (to
a slight degree) on the probability distribution of the processing time S
i
.Using
simulation, we find the values 0.056, 0.058, 0.050, 0.051 and 0.050 for the service
level in the respective regions 1, 2, 3, 4 and 5.
The M/G/∞ queue also has applications in the analysis of inventory systems.
Example 1.1.4 A two-echelon inventory system with repairable items
An approximate analysis of this inventory system can be given by using the
M/G/∞ queueing model. Let (S
0
,S
1
, ,S
N
) be a given design for which S
0
spare parts have been assigned to the depot and S
j
spare parts to base j for
j = 1, ,N such that S
0
+ S
1
+···+S
N
= J . At the depot, failed items arrive
according to a Poisson process with rate
λ
0
=
N
j=1
λ
j
(1 − r
j
. On average λ
0
failed items arrive at
the depot per time unit and on average a failed item at the depot waits W
0
time
units before a replacement is shipped. Thus the average number of failed items at
the depot waiting for the shipment of a replacement equals λ
0
W
0
. This heuristic
argument shows that
L
0
= λ
0
W
0
.
This relation is a special case of Little’s formula to be discussed in Section 2.3.
The relation W
0
= L
0
/λ
0
leads to an explicit formula for W
0
,sinceL
defined as the repair time in case of repair at the base and otherwise as the time
until receipt of a replacement from the depot. Thus the average service time of a
customer at base j is given by
β
j
= r
j
µ
j
+ (1 − r
j
)(τ
j
+ W
0
), j = 1, ,N.
The situation at base j can only be modelled approximately as an M/G/∞ queue.
The reason is that the arrival process of failed items interferes with the replacement
times at the depot so that there is some dependency between the service times at
base j . Assuming that this dependency is not substantial, we nevertheless use the
M/G/∞ queue as an approximating model and approximate the limiting distri-
bution of the number of items in service at base j by a Poisson distribution with
THE POISSON PROCESS 15
mean λ
j
β
j
for j = 1, ,N. In particular,
the long-run average number of back orders outstanding at base j
≈
0
,S
1
, ,S
N
can be
calculated.
1.1.4 The Poisson Process and the Uniform Distribution
In any small time interval of the same length the occurrence of a Poisson arrival is
equally likely. In other words, Poisson arrivals occur completely randomly in time.
To make this statement more precise, we relate the Poisson process to the uniform
distribution.
Lemma 1.1.4 For any t>0 and n = 1, 2, ,
P {S
k
≤ x | N(t) = n}=
n
j=k
n
j
x
t
j
1 −
=
1
P {N(t) = n}
n
j=k
P {N(x) = j, N (t ) −N(x) = n − j}
=
1
e
−λt
(λt)
n
/n!
n
j=k
e
−λx
(λx)
j
j!
e
−λ(t−x)
[λ(t − x)]
n−j
(n − j)!
=
n
p
(1 − y)
q
dy = 1,p,q= 0, 1, .
The right-hand side of (1.1.7) can be given the following interpretation. Let
U
1
, ,U
n
be n independent random variables that are uniformly distributed on
the interval (0,t). Then the right-hand side of (1.1.7) also represents the probability
that the smallest kth among U
1
, ,U
n
is less than or equal to x. This is expressed
more generally in Theorem 1.1.5.
Theorem 1.1.5 For any t>0 and n = 1, 2, ,
P {S
1
≤ x
1
, ,S
n
≤ x
n
| N(t) = n}=P {U
(1)
≤ x
1
denotes the number of passengers joining this crossing and the random variable S
k
represents the arrival epoch of the kth passenger. By conditioning, we find
E(total waiting time)
=
∞
n=0
E(total waiting time | N(T) = n)P {N(T) = n}