Modelling Epidemic Spread using Cellular Automata

backporcupineAI and Robotics

Dec 1, 2013 (3 years and 8 months ago)

79 views

Modelling Epidemic
Spread using Cellular
Automata
Shih Ching Fu
This report is submitted as partial ful¯lment
of the requirements for the Honours Programme of the
Department of Computer Science and Software Engineering,
The University of Western Australia,
2002
Abstract
Most existing epidemic models make simplifying assumptions about the underly-
ing environmental conditions they are trying to emulate.Of particular interest is
the assumption of a homogeneous world where not only are all the hosts identical,
they are uniformly distributed over the landscape.Making such an assumption
limits the realism of these models.To improve realism we should incorporate
spatial heterogeneity.Cellular automata (CA) implicitly model space in their
structure making them an attractive prospect for incorporating spatial factors
into epidemic modelling.Using CA to model epidemics is not a new concept,
however previous models have usually focused on one or two epidemic spread
parameters in order to isolate and analyse their impact.As far as creating a pre-
dictive tool is concerned,we require a model that will take into account all the
major epidemic spread factors and produce a reasonable depiction of the future.
Not only can a CA model trivially implement many of the parameters found
in existing di®erential equation models,but they provide a useful tool for the
graphical visualisation of an epidemic's spread.The model proposed in this
dissertation is far from complete and requires extensive quantitative data for val-
idation purposes.Initial results suggest that through the use of simple localized
interaction rules,an accurate overall epidemic behaviour can be emulated.
Keywords:Cellular automata,epidemic,simulation.
CR Categories:F.1.1,I.6.8
ii
Acknowledgements
First and foremost I would like to thank my supervisor Professor George Milne.
His interest and enthusiasm in my project was motivation for me to work that
extra bit harder.I appreciate that he supervised a number of other Honours
students this year as well.I would also like to thank Dr Gareth Lee for proof
reading my thesis in Professor Milne's absence.
I would also like to thank our Head of Department,Professor Robyn Owens.
Not only did her Scienti¯c Communication unit teach me a lot about writing and
presenting,but the support she has provided the Honours students this year was
invaluable.
Thanks must also go out every member of the Honours Class of 2002.We
are a close knit group and their company made this year all the more enjoyable,
even at the expense of a little productivity!
Finally of course I would like to thank my family for their support and en-
couragement through out this year.
iii
Contents
Abstract ii
Acknowledgements iii
1 Introduction 1
2 Basic Viral and Epidemic Theory 3
2.1 What is a virus?...........................3
2.1.1 De¯nition...........................3
2.1.2 Infection life-cycle.......................3
2.1.3 Transmission.........................5
2.1.4 Mutation...........................5
2.2 What is an epidemic?.........................6
2.2.1 De¯nition...........................6
2.2.2 Transmission probability...................6
2.2.3 Basic reproductive number,R
0
...............7
2.2.4 Virulence...........................7
2.2.5 Periodicity...........................7
3 Cellular Automata 8
3.1 The Cell................................8
3.2 Update rules..............................9
3.3 Interaction neighbourhoods.....................9
4 Review of the Literature 10
4.1 Parameters that in°uence epidemic spread.............10
4.2 Di®erential equation models.....................11
4.3 MFT approximation.........................12
iv
4.4 Spatial models:CA..........................13
4.4.1 Variable population size...................14
4.4.2 Uneven population density..................14
4.4.3 Variable susceptibility....................15
4.4.4 Incubation and latency time.................15
4.4.5 Conclusion...........................15
5 My Epidemic Model 16
5.1 Cell De¯nition.............................16
5.2 World De¯nition...........................17
5.3 Adjustable Simulation Parameters..................17
5.3.1 Neighbourhood radius....................17
5.3.2 Motion probability......................18
5.3.3 Immigration rate.......................18
5.3.4 Birth rate...........................19
5.3.5 Natural death rate......................19
5.3.6 Virus morbidity........................19
5.3.7 Vectored infection rate....................20
5.3.8 Contact infection probability.................20
5.3.9 Spontaneous infection probability..............20
5.3.10 Recovery rate.........................21
5.3.11 Re{susceptible rate......................21
5.4 Cell update algorithm........................21
5.4.1 Movement phase.......................21
5.4.2 Infection and recover phase.................22
6 Control Scenario 23
6.1 World setup..............................23
6.2 Parameter settings..........................23
6.3 Results.................................24
6.4 Discussion...............................25
v
7 Viral Parameter Experiments 27
7.1 Residual immunity..........................27
7.1.1 World setup..........................27
7.1.2 Parameter settings......................28
7.1.3 Results.............................28
7.1.4 Discussion...........................30
7.2 High virulence.............................30
7.2.1 World setup..........................30
7.2.2 Parameter settings......................31
7.2.3 Results.............................31
7.2.4 Discussion...........................33
8 Spatial Experiments 35
8.1 Varied population density......................35
8.1.1 World setup..........................35
8.1.2 Parameter settings......................37
8.1.3 Results.............................38
8.1.4 Discussion...........................39
8.2 Host immigration...........................40
8.2.1 World setup..........................40
8.2.2 Parameter settings......................40
8.2.3 Results.............................40
8.2.4 Discussion...........................43
8.3 Corridors of spread..........................43
8.3.1 World setup..........................43
8.3.2 Parameter settings......................44
8.3.3 Results.............................44
8.3.4 Discussion...........................45
8.4 Barriers to spread...........................45
8.4.1 World setup..........................45
8.4.2 Parameter settings......................47
vi
8.4.3 Results.............................47
8.4.4 Discussion...........................49
9 Conclusion 50
10 Future work 52
10.1 Model calibration...........................52
10.2 Di®erent cell types..........................53
10.3 Other tools..............................53
A SimDemic 57
A.1 Habitable.java...........................57
A.2 World.java..............................58
A.3 Cell.java...............................59
A.4 The GUI................................59
B Original Honours Proposal 61
vii
List of Tables
4.1 Various epidemic models and the realismparameters they implement.13
6.1 Epidemic spread parameter values for the control scenario.....24
7.1 Parameter values for the residual immunity experiment......28
7.2 Parameter values for the virus morbidity experiment.......31
8.1 Parameter values for the population density experiment......37
8.2 Parameter values for the host immigration experiment......41
viii
List of Figures
2.1 Infection life-cycle...........................4
3.1 A general CA cell state transition..................8
3.2 Commonly used interaction neighbourhoods............9
6.1 Control scenario:ratio of infectives to total population over time.25
6.2 Control scenario:total population over time............26
7.1 Residual immunity:ratio of infectives to total population over time 29
7.2 Residual immunity:total population over time...........29
7.3 High virulence:infective ratio over time..............32
7.4 High virulence:total population over time.............32
7.5 High virulence:e®ects of immunity.................34
8.1 High to low population density gradient..............36
8.2 Low to high population density gradient..............36
8.3 Decreasing population density:ratio of infectives.........38
8.4 Increasing population density:ratio of infectives..........39
8.5 Population dynamics of the closed population experiment.....42
8.6 Population dynamics of the open population experiment.....42
8.7 Corridors of spread..........................44
8.8 Corridors of spread:lag map.....................46
8.9 Barriers to spread...........................47
8.10 Barriers to spread:lag map.....................48
A.1 A screen capture of the SimDemic application............60
ix
CHAPTER 1
Introduction
Public health issues are seeing greater visibility in the media;a particular concern
is viral spread through populated areas.Decreased worker productivity as a
result of viral illness costs industry millions of dollars every year [4].With recent
virus epidemics such as foot and mouth disease in the United Kingdom and
the threat of using viruses such as smallpox as biological warfare agents,the
monitoring of outbreaks is gaining importance for governments and public health
o±cials.It therefore becomes desirable to predict patterns of viral infection
under certain environments.It is hoped that the modelling of a phenomenon
such as virus spread will help us understand,predict,and ultimately control that
phenomenon's behaviour.
The majority of existing epidemic models utilize di®erential equations and do
not take into account spatial factors such as population density.These models
assume populations are closed and well mixed;that is,host numbers are constant
and individuals are free to move wherever they wish.When trying to devise more
realistic models it makes sense to incorporate spatial parameters that re°ect the
heterogeneous environment found in nature.An alternative to using deterministic
di®erential equations is to use a two-dimensional cellular automaton that relies
upon stochastic parameters.
The focus of this dissertation is to analyse how we can use a two-dimensional
CA to model the spread of a viral epidemic.Of particular interest is the manner
in which CA discretize space and I hope to show that spatial heterogeneity can be
easily incorporated into a CA epidemic model.The hypothesis is that the overall
behaviour of an epidemic can be produced from the summation of localized host
interactions.To demonstrate this hypothesis,it will be necessary to devise an
epidemic model and the show that this model can replicate the results of existing
non-CA models,as well as taking space into account.
Before beginning the modelling process,modellers must decide what level
of detail their models will examine.This is an important decision because too
much detail will produce a cluttered and unworkable model,whereas too little
1
detail provides no useful information.Chapter 4 outlines some of the existing
approaches that have been adopted for epidemic modelling and looks at which
epidemic spread factors,or level of detail,they have chosen to implement.After
a brief overview of basic epidemic theory and cellular automata in Chapters 2
and 3,I will describe the elements of my new composite epidemic model and the
epidemic spread parameters I have chosen to implement.Chapters 6,7 and 8
describe the experiments I have conducted to analyse the suitability of CA to
epidemic modelling,which also serve to partially validate the epidemic model I
devised as a part of this project.
2
CHAPTER 2
Basic Viral and Epidemic Theory
In contrast with medicine,epidemiology studies the health of entire populations
rather than just the individual.This chapter outlines some of the concepts as-
sociated with what is traditionally a statistical ¯eld and de¯nes some of the
terms you are likely to encounter in epidemiological articles and publications.
Although epidemics can stem from any infectious disease,I will focus solely on
viral epidemics.
2.1 What is a virus?
2.1.1 De¯nition
Andrewes [3] describes viruses as\on the borderline between life and death".
Under a microscope a virus resembles little more than a lifeless geometric crystal.
However,when a virus penetrates the cell of a living organismit is able to rapidly
replicate itself.
A virus comprises two parts:some genetic material and a protective protein
coating.After a virus has penetrated a living cell it rewrites the cell's DNA
and transforms it into an organism that produces hundreds or even thousands
of copies of itself.When these virus copies leave the infected cell they are once
again lifeless until they penetrate another cell.Disease arises from cell damage
caused by the genetic rewriting procedure.
2.1.2 Infection life-cycle
The aim of a virus entering a living cell is to replicate itself,however as the host
acquires antibodies to ¯ght the infection the virus will need to ¯nd another host.
Consequently,after being infected,an individual usually becomes infectious as
well as diseased.A diseased host is one that shows symptoms { these symptoms
3
Figure 2.1:The relationship between infectiousness and virus symptoms [21].
The labels above the line describe host infectiousness while the labels below the
line describe disease dynamics.Note that the infectious period can start before
or after the onset of symptoms.
are usually mechanisms to help the virus spread to new hosts,for example,cough-
ing.However,the disease cannot be too debilitating because if the host dies,so
does the contagion.The relationship between each of the phases of infection are
shown in Figure 2.1,with each period described below.
Latent period This is during the early stages of an infection where the virus is
yet to develop the ability to transfer to a new host.
Infectious period During this phase the virus is contagious and can be passed
onto other individuals through the virus'natural spread mechanisms.
Recovered or Removed As far as the virus is concerned,a host that has
gained a natural immunity or has died is no longer able to contribute to
the replication process.In either case,the virus cannot spread further.
Incubation period Early into an infection there may be no signs of infection
at all;this initial stage is known as the incubation period.It is during the
intersection of this period and the infectious period where viruses spread
the most [21].This is because hosts are unaware of their`carrier'status
and continue normal contact with other healthy hosts.
Symptomatic period This is the stage during the infection where there are
visible signs of infection,for example,an infected human would go and see
a doctor.For viral infections,treatment usually only comprises relieving
symptoms and isolation away from other healthy individuals.
4
2.1.3 Transmission
Within the body a virus can pass between cells via contact with adjacent cells.
However,outside of the body there are ¯ve main ways a virus can be transmitted
between hosts [3]:
² respiratory transmission,
² via food or faeces,
² mechanical transmission,
² via living carrier or agent,or
² vertical transmission.
Respiratory transmission includes the inhalation of droplets that contain the
virus.These droplets could be a result of the coughing and sneezing by infective
hosts typical of the in°uenza virus.The Hepatitis A virus can be contracted from
consumption of contaminated food or water.Mechanical transmission occurs
when virus particles enter through the skin such as through cuts.Transmission
due to a living carrier is also known as vectored infection.Vectored infection can
arise from tick bites or in the case of rabies,dog bites.Vertical transmission
refers to the transfer of contagion from parent to child during childbirth.
2.1.4 Mutation
What is often called viral mutation is really accelerated natural selection [3].
Viruses evolve like many other organisms but they have an advantage over more
complex lifeforms because they can multiply rapidly.By having millions of de-
scendants in a short space of time,the e®ects of natural selection are felt sooner
so that the resulting strains of virus are the ones that are the most successful at
duplicating and surviving.In general,the most successful parasites are the ones
that do not greatly harm their hosts,whereby reaching a state of`mutual toler-
ation'.A hidden virus will not spread and a highly virulent strain will kill the
host;the balance between the two extremes is attained through natural selection.
Successful mutations usually undermine the antibodies produced by recovered
hosts and make them resusceptible to further infection.For some viruses,such
as in°uenza and HIV,this rapid mutation makes vaccination di±cult.
5
2.2 What is an epidemic?
2.2.1 De¯nition
Aviral outbreak occurs when the number of cases of a particular virus or disease is
higher than the normally expected or endemic level of infection.Di®erent viruses
in di®erent regions of the world will have di®erent threshold values for what is
classed as an outbreak.In°uenza,with thousands of cases in the United States
every year,is considered endemic in many parts of the world [10].However,just
one case of smallpox in any country is considered an outbreak.An outbreak is
upgraded into an epidemic when it becomes prolonged and rapidly spreads to
neighbouring areas [21].An epidemic which spreads to cover continents is called
a pandemic.An example is the in°uenza pandemic of 1918 that killed roughly
40 million people across 4 continents [6].
2.2.2 Transmission probability
The transmission probability of infection refers to the chance that there is a
successful transfer of the virus from one host to another.Estimates of this
probability are useful to the epidemiologist in understanding the dynamics of
an epidemic [21].A good estimate of the transmission probability is found by
calculating the secondary attack rate.
Secondary attack rate
The secondary attack rate (SAR) is a measure of contagiousness and is de¯ned as
the ratio of individuals who develop an infection to the total number of susceptible
individuals.This is shown in Equation 2.1.It is deemed a secondary attack rate
because it refers to the infections that occur from a primary`source'.
SAR =
number of individuals that develop the disease
total number of susceptibles
(2.1)
It is calculated by identifying the infective hosts,tracking which healthy hosts
come in contact with them,and then noting which become infective as well.It is
important to note that the SAR is a value that is calculated in retrospect from
collected data rather than a metric that can be predicted.Consequently,it is a
static ¯gure that is averaged over the entire epidemic and does not provide an
instantaneous measure of epidemic spread.
6
2.2.3 Basic reproductive number,R
0
For viral outbreaks,the basic reproductive number,R
0
,is the mean number of
susceptibles that an infective host infects during its infectious lifetime.This num-
ber only includes secondary (direct) infections not tertiary ones.For example,if
R
0
= 6,we would expect an average of 6 secondary infections for each primary
infection.If R
0
= 1,then the number of infectives remains relatively constant
and the virus is deemed to be endemic.The basic reproductive number is a
function of three parameters,as shown in Equation 2.2.
R
0
= c £p £d
where c is the number of contacts per unit time,
p is the transmission probability,and
d is the duration of infectiousness
(2.2)
2.2.4 Virulence
Virulence,also known as virus morbidity,is a measure of how rapidly a virus
kills its host and is inversely proportional to R
0
.Highly virulent viruses will
have R
0
<< 1 and usually result in acute outbreaks where many die but the
virus does not spread far.The probability of dying from the infection before
recovering or dying from other causes is known as the case fatality ratio.
2.2.5 Periodicity
A common feature in epidemics is the damped oscillatory nature of the number of
infections over time [25].After an initial`boom'in the number of infections,the
number of infectives drops as the population acquires immunity.It might be ex-
pected that after reaching an endemic level there would be few future outbreaks,
but historic data has shown that the infective population oscillates.These oscil-
lations seem to continue inde¯nitely but become successively weaker.Epidemic
theory suggests that this periodicity is due to the turnover in host populations;
when the previous generation of immune hosts dies,a new injection of suscepti-
bles to perpetuate the outbreak is born.
7
CHAPTER 3
Cellular Automata
Cellular automata (CA) are dynamical systems characterized by their discretiza-
tion of time and space [9].Typically,a cellular automaton comprises an array
or lattice of automata that evolve over discrete time quanta.This lattice can be
n{dimensional,but is usually one or two dimensional.This chapter provides an
overview of CA and its inherently simple nature.
3.1 The Cell
The most basic component in a CA is the cell.Traditionally,each cell is a ¯nite
state automaton (FSA) that evolves according to a pre-de¯ned update rule.The
next state of a cell is a function of its present state and the current inputs as
shown in Figure 3.1.Classically,cells are square and placed side by side to form
a lattice,however,there are no formal restrictions on the size or shape of the
cells,their arrangement in the lattice,or whether all the cells must be identical.
Figure 3.1:A generalized state transition.A cell's next state depends on the
current states of its neighbours.
8
3.2 Update rules
The present state of a CA is de¯ned as the set comprising the current state of
all it cells.These states are in turn are governed by a global update rule [12].
Although called an update rule,it usually consists of a list of criteria rather than
just one.The next state of a cell is a function of its current state and the state
of its interaction neighbourhood.It is up to the implementer as to what the rules
contain and whether each cell must obey the same update function.The classic
example of a CA with a uniform update rule is John Conway's Game of Life [17].
3.3 Interaction neighbourhoods
As stated earlier,a cell's next state will depend on the current state of its neigh-
bours.Before deciding on its next state,a cell will interrogate its neighbours for
their present states and then evolve accordingly.Once again the size and shape
of the interaction neighbourhood is up to the implementer and will vary from
application to application.Figure 3.2 shows some of the more commonly used
neighbourhoods.
Figure 3.2:The black dot represents the target cell { the shaded cells are its
neighbours.The grid on the left shows the 4-connected,Von Neumann neigh-
bourhood.The centre grid shows the 8-connected,Moore neighbourhood.The
right-hand neighbourhood uses a hexagonal lattice to provide equal separations
between cells.
9
CHAPTER 4
Review of the Literature
Epidemic spread models are devised for a number of reasons.Epidemiologists
want to create simple models so they can test the impact of speci¯c parameters
on an epidemic's overall behaviour.Computer scientists might want to create
software packages that accurately depict the behaviour of epidemics as seen in
nature.The latter is the motivation for this dissertation.
There exist many models of epidemic spread,each with its own approach and
set of assumptions.However,these models all share one property:the virtual
world in which they run is an idealized one where noise and imperfections are
¯ltered out.This arises from the di±culty of incorporating all the variables we
see in nature into a simulation that has a reasonable execution time | hours
rather than days.When modelling a complex system there is a trade o® between
a model's degree of abstraction and its usefulness;that is,without devaluing the
results a model provides,which details can be left out?
In this chapter I examine some of the existing epidemic modelling techniques
and compare their levels of detail and realism.The majority of past approaches
have used ordinary and partial di®erential equations (ODE's and PDE's).I
examine those as well as mean ¯eld type (MFT) approximations [18] and cellular
automata (CA).
4.1 Parameters that in°uence epidemic spread
The main focus of my project is on the spatial behaviour of epidemics rather than
absolute numbers and densities of infected individuals.In order to objectively
compare the relative usefulness of modelling approaches I use a set of standard
criteria.These criteria are based on the following factors that Mollison [19]
describes as signi¯cant in determining epidemic spread:
² susceptible population size,
10
² homogeneity of population density,
² infection transmissibility,
² immunity levels of individuals,
² motion of individuals,and
² infection incubation time.
All of the above will in°uence whether an infection will rapidly propagate
through a population or head into extinction.Di®erent models have di®erent
assumptions regarding the above parameters and hence have di®erent success
rates in mimicking nature.The rest of this chapter examines three modelling
methodologies:ODE's and PDE's,MFT approximations,and CA to contrast
how each approach incorporates the above listed epidemic spread parameters.A
model's omission of a parameter does not necessarily imply that the model is
unrealistic,though it might mean that the model is designed to investigate the
impact of one particular parameter independently of the others.
4.2 Di®erential equation models
Deterministic approaches to modelling,such as those using ODE's and PDE's,
are poor at representing small populations compared to probabilistic models such
as CA [22].They are poor in the sense that the results from ODE's diverge to
unreasonable values as the population size is scaled down toward zero.This
divergence is due to the simplifying assumptions made in ODE models:
² Population sizes are constant:no births,deaths or immigration occurs in
the world.
² Populations are uniformly distributed over the world.
² Populations are well mixed,that is,there is homogeneous motion about the
world.
Most of the above assumptions arise fromregarding susceptible populations as
continuous entities rather than comprising discrete individuals.Boccara [7] states
that by recognizing that the spatial behaviour of an epidemic is\strongly linked
to the short range character of the infection process",continuous di®erential
11
equation models that neglect the individual are probably going to be misleading
on all scales,not just small.
However,Di Stefano [22] shows that the continuous nature of ODE approx-
imation is well suited to dealing with large populations because the e®ect of
the close contact between individuals becomes negligible compared with the epi-
demic's macroscopic behaviour.The local correlations are lost because individu-
als are able to roam all over the world.This spatial mixing e®ect is introduced
into PDE models through a di®usion term and become a function of the initial
and boundary value conditions.
The need for homogeneity in an ODE model means that the natural pro-
gression of an epidemic is not represented accurately.Variations in localized
population densities,variations in immunity and susceptibility,and variations in
incubation and sickness time are all attributes of natural epidemics but omitted
in ODE simulations.Given that many spatial properties of epidemics are not
realized by di®erential equation models,another paradigm should be chosen if
space issues are to be accurately addressed.
4.3 MFT approximation
Closely related to ODE and PDE models are those using mean ¯eld type (MFT)
approximation.Kleczkowski [18] proposed an epidemic model using MFT ap-
proximations examining the spread of childhood measles.
MFT approximation ignores localized correlations making it similar to a dif-
ferential equation model [8].MFT models assume that susceptible population
density is uniform over the world,which is also similar to ODE/PDE models.
Finally,MFT and ODE/PDE models share the assumption that hosts are capa-
ble of di®using around the world { that is,the population is well mixed.Despite
all these similarities,MFT approximations are signi¯cantly di®erent from dif-
ferential equation models in that their mixing parameter can be a probabilistic
variable.Unlike ODE's where either all individuals di®use or none at all,the
decision to move around the MFT world is independent among individuals.
Similar to CA,MFT approximations utilize a lattice structure to emulate
the spatial nature of epidemic spread.Each lattice site contains an individual
who can exist in one of several states.The set of possible existence states is
determined by the epidemic model being used.Many groups [7,8,13] have
used MFT approximations as a point of comparison with the CA models they
develop;particularly models investigating the e®ect of host motion.Noting that
MFT approximations neglect localized correlations,it appears meaningless to
12
compare them with CA models because CA models focus on the contact between
individuals.But as Kleczkowski [18] describes,MFT and CA models converge
when the MFT mixing parameter tends to in¯nity { that is,when the world
contains more disorder than correlation.This convergence is analogous to the
situation where di®erential equation and CA models converge when population
size tends to in¯nity.
In the context of exhibiting realistic spatial behaviour,MFT approximations
are better than modelling with di®erential equations because by determining
host movement probabilistically,they manage to partially portray the stochastic
°uctuations observed in nature.
4.4 Spatial models:CA
Epidemic spread in nature is a stochastic process,so it seems logical to use a
model that is probabilistic.According to Ahmed et al.[2],\[CAhave] a signi¯cant
role in epidemic modelling since it can be shown that [they are] more general than
ordinary and partial di®erential equations."This section examines existing CA
models and identi¯es which virus spread parameters they have chosen to include
and which they have chose to omit.
All of the models that I discuss in this section are based on the SIR model
superposed with CA.The letters in SIR correspond to the three states an individ-
ual can exist in:susceptible,infective and recovered [1].Susceptible individuals,
or susceptibles,are ones who can contract the pathogen from already infected
infective individuals.Infectives can later recover from this infection.There are
variants of this model that introduce other intermediate states.For example,in
SEIRthere is the`E'state,representing exposed but yet to be infected individuals.
Criteria
[13] [7] [22] [1] [8]
Wrap around world
£ £ £ £ £
Variable population size
£
Uneven population density
£ £ £
Movement of hosts
£ £ £
Immunity after recovery
£
Variable susceptibility
£
Includes incubation time
£
Includes latency time
£
Table 4.1:Various epidemic models and the realism parameters they implement.
13
The rest of this chapter looks at the models of Duryea et al.[13],Di Stefano et
al.[22],Ahmed and Agiza [1],and Boccara et al.[8] to see which epidemic spread
factors they have chosen to incorporate into their models.Such factors include
population density,host susceptibility and immunity,virus transmissibility,and
virus infection times.The inclusion of these parameters will directly a®ect the
realism of the simulations we generate { that is,they are parameters that add
heterogeneity to the otherwise idealistic world that designers build models in.
Table 4.1 shows the parameters each group has chosen to implement.
4.4.1 Variable population size
In nature,the population within a region is always changing.Internal events such
as births and deaths increase and decrease the population respectively.External
factors such as immigration and emigration have similar e®ects.Epidemic spread
is a®ected by this constant °ux in susceptible hosts but very few models include
this °ux.The reason perhaps lies in the di±culty to integrate such features into
a model,or more probably,these models have very speci¯c applications where
population variation is considered negligible.
Of the ¯ve groups mentioned here,only one,Boccara et al.[8],has chosen to
model population °ux.When the proportion of infectives in a population reaches
a constant steady-state value we say that the infection has become endemic.
Boccara et al.chose to include the death and birth rates to investigate how these
parameters a®ect the stability of such endemic states.Other models such as those
by Di Stefano [22] focus on the movement and the heterogeneity of susceptibility
in populations and choose to abstract out the population size parameter.
4.4.2 Uneven population density
The model of Boccara and Cheong [7] tries to introduce variations in population
density by allowing host movement.In e®ect,the cell update function is divided
into two phases:infection and motion.Once each cell has decided which of the
SIR states it evolves into,it chooses a destination cell and moves there.This has
the e®ect of generating a non-uniform landscape of susceptible hosts.
Ahmed and Elgazzar [2] also model variations in population densities by al-
lowing host movement;more speci¯cally,cyclic host movement.PDE models use
a di®usion term that creates random host motion whereas a cyclic cell to cell
mapping function in a CA model makes it possible to emulate regular periodic
host movement.This is analogous to someone going from home to work,and
then back home again.
14
4.4.3 Variable susceptibility
Susceptibility,immunity,and transmissibility relate to how easily a contagion can
pass between hosts.Probabilistic CA/SIRmodels can represent highly contagious
viruses by assigning a high probability of infection when susceptibles come in
contact with the contagion.Although these parameters are usually statically
de¯ned as in the model proposed by Ahmed and Agiza [1],the rule-based nature
of CA means that these parameters can be modi¯ed dynamically during run-
time.A population's innate immunity to infection has a signi¯cant impact on
whether a small outbreak grows into an epidemic or simply dies out;the model
by Ahmed and Agiza appears to solely address this susceptibility parameter and
neglects the others to perhaps isolate its e®ects.
4.4.4 Incubation and latency time
The last parameters in Table 4.1 to be discussed are incubation and latency time;
these two terms should not be confused.The distinction is outlined in Figure 2.1.
Incubation time is de¯ned as the length of time between an individual being
infected and an individual showing signs of disease.The time lapse between being
infected and becoming infective is known as latency.Both of these quantities
are modelled by Ahmed et al.[1],however they do not discuss the signi¯cance of
these times.Others fail to include these quantities in their models but Ahmed and
Agiza suggest that these parameters have an accelerating impact on an epidemic's
spatial spread.
4.4.5 Conclusion
Most of the models examined in this chapter have a speci¯c focus:they isolate
a particular epidemic spread parameter and see its e®ect over the epidemic as
a whole.As part of determining the suitability of CA to epidemic modelling,it
will be necessary to take these speci¯c models and try to compose them into a
single generalized composite model.By examining the merits of this composite
model we can endorse or reject the conjecture that CAis appropriate for epidemic
modelling.
15
CHAPTER 5
My Epidemic Model
The cell evolution in a cellular automaton follows an update function that takes
the state of a particular cell and its neighbours and determines the next state.
This chapter outlines the cell and update rule de¯nitions adopted in my epidemic
simulation as well as the epidemic parameters that are implemented.A more
detailed description of the software implementation of this model can be found
in Appendix A.
5.1 Cell De¯nition
The basic unit of computation in my simulation is the cell.Here,`cells'are
automata cells and not biological cells.However,rather than representing a par-
ticular area of space,the e®ective size of a cell will be determined by the epidemic
spread parameters set by the user.For example,by de¯ning a small value for the
mobility,the user is in e®ect simulating the increased di±culty of traversing a
cell that represents a large area compared to a cell that represents a small area.
The concept of a cell needs to be di®erentiated from the hosts that live inside
the cell,as illustrated by the following ¯ve cell attributes:
² carrying capacity,
² total population,
² susceptible subpopulation,
² infective subpopulation,and
² recovered subpopulation.
The carrying capacity of a cell is used as a mechanism to limit the movement
of hosts between cells.It is a mechanism used to prevent crowding within a
16
particular cell { comparable to a surface area.The number of newborns will also
be dependent on whether a cell has reached its carrying capacity.Although the
e®ect of the land's carrying capacity is not directly enforced in nature,for the
purposes of simulation,it is a straightforward way to encourage or attenuate the
motion of individuals between cells.
Traditionally,such as in the CA of Jon Von Neumann [9],a single automaton
occupied each of the cells that constituted the larger cellular automaton.Rather
than stipulate this,this model allows multiple individuals to dwell in one cell up
to the above-mentioned carrying capacity.Variable cell population has two main
purposes:¯rst it reduces the total number of cells and hence reduces computation
time;secondly it provides generality.If we set the carrying capacity to one we
can revert back to a traditional cellular automaton.
5.2 World De¯nition
A two{dimensional array of cells and the epidemic spread parameters that govern
their evolution constitute the world that the hosts`live'in.The cells are arranged
in a rectangular grid comprising square cells with external dimensions that may
or may not be square.The world boundaries serve as impenetrable barriers to
host movement,but conceptually can be thought of as political boundaries that
allow immigration into and out of this world into adjacent worlds.The adjustable
epidemic parameters that control cell evolution are described in the next section.
5.3 Adjustable Simulation Parameters
After researching existing epidemic models,particularly those examining virus
pathogens that can survive outside the bodies of hosts,I compiled the following
list of epidemic spread parameters.This is not an exhaustive list,but it contains
what I believe to be the most signi¯cant factors that account for the spatial
behaviour of an epidemic.Apart from the interaction radii,all the following
parameters are modelled using probabilities that directly impact the update rules
applied over the CA lattice.
5.3.1 Neighbourhood radius
This parameter determines the size of the interaction neighbourhood that a cell
interrogates for state information.My simulation uses a square interaction neigh-
17
bourhood whose area,n,in terms of cells is determined by an interaction radius,
r,as shown in Equation 5.1.
n = (2r +1)
2
(5.1)
There are two distinct interaction radii:motion and infection.The motion
radius de¯nes the greatest distance,measured in cells,a host can move during a
time step.The infection radius is slightly di®erent to the motion radius in that
it does not relate to hosts but to the virus pathogen.The infection radius de¯nes
the greatest distance the virus pathogen can travel outside the body of a host on
its own.This quantity is used to model the spread of a virus via natural vectors
such as airborne droplets in in°uenza or vermin as in bubonic plague.
5.3.2 Motion probability
The individuals in the world are permitted to move between cells.The motion
probability determines the frequency of this motion.This simulation assumes
homogeneous mixing within each cell,but the motion of hosts between cells is
limited by the motion probability parameter.For example,if the motion prob-
ability,p
move
,is set to 0:4,you would expect roughly two in ¯ve hosts to shift
from the cell they currently reside in to another cell within its motion neighbour-
hood.The success of a host's inter-cell movement is dependent on whether the
destination cell has reached its carrying capacity.In this model,the destination
cell is selected randomly from the surrounding interaction neighbourhood.
5.3.3 Immigration rate
This simulation can model an open or closed world,depending on the value of the
immigration parameter.For the purposes of this simulation,immigration refers
to susceptible individuals who are not already in the world coming in and ¯lling
any vacancies within the cell grid.This parameter determines the probability
that a cell will receive any immigrants.The actual number of immigrants will
be a random proportion of the number of spaces left in the cell.Immigration is
equally likely for all cells except for cells which cannot support any more hosts
(the cell is at its carrying capacity).
18
5.3.4 Birth rate
The birth rate parameter,as its name suggests,controls the addition of newborns
to the population pool at each time step.These increases are from births of
new individuals via reproduction between existing hosts.It does not include
population increases due to immigration.
The birth rate is given as a parameter,p
birth
,which is the probability that
a pair of hosts will produce an o®spring during a time step.The birth rate is
uniform across all the possible di®erent pair combinations such as SS,SI,SR,II,
IR,and RR,but all newborns are susceptibles.That is,infection or immunity is
not passed onto the o®spring.For many diseases this is a reasonable assumption
and makes the model simpler,but as discussed in Chapter 2 the contagion might
be passed onto a child via vertical transmission.Equation 5.2 shows the relation-
ship between the birth rate parameter and the total cell population increase per
time step.
% population increase per timestep ¼
p
birth
2
£100 (5.2)
Equation 5.2 is only an approximation because births are determined proba-
bilistically.
5.3.5 Natural death rate
Comparable to the birth rate parameter,deaths by natural causes are controlled
by a natural death parameter.Unlike the birth parameter,the death parameter
encompasses emigration from the world as well.This natural death parameter
does not discern between the host types.That is,susceptibles,infectives and
recovered are equally a®ected by this parameter.The average lifetime of a healthy
host is approximated by the reciprocal of the natural death probability.
5.3.6 Virus morbidity
The neighbourhood radius,natural birth,immigration,and natural death pa-
rameters mentioned above are all geographic or demographic factors that a®ect
epidemic spread.Virus morbidity and the parameters discussed beloware deemed
viral spread factors.
Virus morbidity is a measure of how rapidly a virus kills its host,if at all.
Unlike natural deaths,the virus morbidity parameter only a®ects infective hosts.
19
This parameter de¯nes the probability that an infective will die fromdisease dur-
ing a particular time step.For highly virulent viruses,the morbidity probability
will be close to unity,whereas very`mild'viruses will have this parameter value
very close to zero.
5.3.7 Vectored infection rate
Apart from entering an otherwise`clean'cell inside a mobile infective host,the
virus contagion can spread across cells using its natural spreading mechanisms.
Such spread is known as vectored infection.For this model,the velocity of inter-
cell infection is controlled through the p
input
parameter.The actual probability
of spread,p
vectored
,is a function of the vectored infection parameter provided by
the user,p
input
,and the density of susceptible hosts in the local interaction neigh-
bourhood.This is illustrated in Equation 5.3.A high value for this parameter
would represent a highly transmissible virus such as one that was airborne or
carried in bird droppings.
p
vectored
=
susceptible population
neighbourhood capacity
£p
input
(5.3)
5.3.8 Contact infection probability
Whilst the previous parameter handled the infection of hosts across cells,the con-
tact infection probability determines how susceptible hosts are infected through
contact with infectives within the same cell.Mixing within each cell is assumed
to be homogeneous { all hosts in a cell will come in contact with all the other
hosts during a time step.The contact infection parameter is a measure of virus
contagiousness { highly contagious viruses need only a small amount of exposure
between infective and susceptible hosts to pass on the contagion.
5.3.9 Spontaneous infection probability
To simulate the infection of individuals from factors external to the world,sus-
ceptibles may be spontaneously`struck down'with the virus.This parame-
ter implies that although a susceptible individual might have escaped infection
through vectored and direct contact means,it might still become infective from
outside means.Possibilities might include contracting the virus whilst on hol-
idays,breathing in vermin faeces,contact with virus particles adhered to the
tyres of vehicles,or other remote but still probable methods of infection.This
20
parameter is usually very small compared to the value of the other infection
probabilities.
5.3.10 Recovery rate
In this simulation,recovery corresponds to a state change from infective to re-
covered;it does not encapsulate death or emigration.Recovered individuals are
sometimes known as removed individuals because they no longer contribute to the
infection cycle;not only are they immune from contracting the contagion,they
cannot pass the contagion onto other susceptible hosts.The recovery rate param-
eter is the probability that at a particular time step an infective host becomes a
recovered host.The reciprocal of this parameter provides an approximation to
the mean infection duration.
5.3.11 Re{susceptible rate
The duration of a recovered host's immunity is determined by the re-susceptible
parameter.This probability controls the chance of a recovered individual as-
similating back into the susceptible population { in e®ect it is a recovered host
becoming a susceptible host again.For viruses that mutate rapidly,this parame-
ter will be close to unity,whilst viruses that o®er lifelong immunity after infection
will have a re{susceptibility close to zero.
5.4 Cell update algorithm
The CA cell update function is used to evolve each cell to its next state.The cell
update function takes as arguments all the parameters outlined in the previous
section along with the state information of the interaction neighbourhood of the
cell in question.The update of the world is done in two phases:¯rst the motion
of individuals between cells,then the evolution of individuals within cells.
5.4.1 Movement phase
1.Select a random cell from the world.
2.For each of the individuals in the cell,randomly select a neighbouring cell
and move the individual into it.This movement is dependent on the the
21
motion probability parameter described earlier and the destination cell not
being full.
3.Repeat from step one until all the cells in the world have been accounted
for.
5.4.2 Infection and recover phase
1.Select the ¯rst cell.
2.Deduct the`natural deaths'from the cell population.
3.Deduct the deaths from virus morbidity.
4.Add to the population any newborns.
5.Add any immigrant population.
6.Compute inter-cell (vectored) infections.
7.Compute intra-cell (contact) infections.
8.Compute spontaneous infections.
9.Compute host recoveries.
10.Compute re-susceptibles.
11.Repeat from step one for the next cell until all cells have been accounted
for.
22
CHAPTER 6
Control Scenario
Whenever a model of a natural phenomenon is devised,it needs to be calibrated
such that any numbers that we supply as input or receive as output can be readily
interpreted as`real-life'¯gures.However,creating a fully calibrated model is
beyond the scope of this project as I am only interested in examining whether
CA is a good starting point for a comprehensive model.To substitute calibration,
I have devised a control scenario which will be used as a comparison for all the
experiments discussed later.
6.1 World setup
The control world comprises a 100 £100 grid of cells,with each cell possessing
the characteristics described in Chapter 5.Each cell will start with a population
of 100 susceptibles and be capable of supporting a maximum of 200 individuals.
There is no initial source of infectives.The world starts`clean'and waits for
spontaneous outbreaks to spark an epidemic.A point source of infectives was
not used because later experiments rely on spontaneous infections to seed the
epidemics.This is to try and emulate how epidemics start in nature.
6.2 Parameter settings
The control scenario represents a homogeneously distributed population where
the total population is in a state of dynamic equilibrium { that is,the num-
ber of births roughly balances the number of deaths.However,the population
is expected to gradually increase due to immigration from outside worlds.A
non-zero motion probability means that at each time step,it is expected that
approximately 0.1% of individuals will attempt to shift out of their current cell
into another one in their interaction neighbourhood.The initial values of the
epidemic spread parameters described in Chapter 5 are summarized in Table 6.1.
23
Parameter
Value
Interpretation
Infection radius
1
adjacent spread only
Movement radius
1
adjacent movement only
Immigration rate
0.01
1% increase in population per time step
Birth rate
0.02
1% increase in population per time step
Natural death rate
0.01
1% decrease in population per time step
Virus morbidity
0.05
5% of infectives die per time step
Spontaneous infection rate
0.0001
1 in 10000 chance of outside infection
Vectored infection rate
0.2
20% chance of hostless inter{cell spread
Contact infection rate
0.4
2 in 5 chance of infection after close contact
Recovery rate
0.1
10% of infectives recover per time step
Resusceptible rate
0.001
residual immunity lasts ¼1000 time steps
Movement probability
0.001
1 in 1000 individuals attempt to shift cells
Table 6.1:Epidemic spread parameter values for the control scenario.
Whilst it is simple to de¯ne a typical host population,it is much more di±cult
to de¯ne a typical virus strain for the control scenario without losing generality.
To that end,although the values in Table 6.1 do not quantitatively represent
a particular virus,the ratios between the parameters need to be reasonable.
Essentially these numbers are just educated guesses at what a typical virus would
be.
6.3 Results
Figure 6.1 shows the °uctuation in the proportion of infectives in the world over
800 time steps.It is important to note that after approximately 100 epochs,
the proportion of infectives levels out.I will de¯ne this steady-state value as an
endemic level of infection in the population.
The plot in Figure 6.2 provides an indication of the dynamics of the total
population,not just the infective hosts.Starting from an initial population of
one million the population rises by approximately 12% and at about the ¯fti-
eth epoch the population starts to fall sharply to a much smaller value than it
started.A comparison of Figures 6.1 and 6.2 shows that the rapid decline in pop-
ulation corresponds to the sharp increase in infectives.The population reaches
a local minimum just before the hundredth epoch and then starts to regenerate.
This regeneration corresponds almost exactly with the point in Figure 6.1 when
the proportion of infectives bottoms out.As the population recovers it reaches
24
0
100
200
300
400
500
600
700
800
0
0.1
0.2
0.3
Control scenario
Time (epochs)
Proportion of infectives
Figure 6.1:The proportion of the total population that is infective during the
control scenario.Froman initial infective population of zero,the infection spreads
rapidly before peaking at ¼ 0.26.From there,the transient rapidly decays and
reaches an endemic proportion of ¼ 0.1.
a steady state value which is roughly one million { its initial size before the
epidemic.
6.4 Discussion
The shape of the graph in Figure 6.1 is one which will be seen often in later ex-
periments.The initial surge in infectives is a result of the virus spreading rapidly
through a purely susceptible population with no natural immunity.However,as
hosts recover and acquire immunity,the number of infectives quickly reduces to
an endemic level.For any subsequent epidemics to occur the population will need
to lose its natural immunity or be injected with a large number of infectives.
This control experiment provides some initial evidence for the suitability of
CA for epidemic modelling.The dynamics seen in Figure 6.1 are similar to those
presented by Boccara et al.[8] and Kleczkowski [18] in their models that use
mean ¯eld type approximations.However,they do not mention whether they
have used real life data to validate the results from their theoretical models.
The main purpose of this experiment is to provide a benchmark for the fol-
lowing experiments.As parameters and initial conditions are varied in later
25
0
100
200
300
400
500
600
700
800
0.85
0.9
0.95
1
1.05
1.1
1.15
x 10
6
Control scenario
Time (epochs)
Total population
Figure 6.2:The °uctuation of the total world population over 800 time steps.
After some rapid growth in the ¯rst 50 time steps the population declines just as
rapidly.This decline corresponds with the peak in Figure 6.1.From time = 100
onward,the population recovers to roughly its original size.
experiments,new results can be interpreted with respect to the results presented
here and their consequences discussed.
26
CHAPTER 7
Viral Parameter Experiments
This section investigates whether a CA epidemic model can encapsulate the ca-
pabilities of existing,non-spatially oriented models.Speci¯cally the viral param-
eters of immunity and virulence are examined.Each of these experiments has
the same starting state as the control scenario,except that some of the viral
parameters will be di®erent.
7.1 Residual immunity
A population's natural immunity to a viral contagion will limit the chances of a
small outbreak escalating into a large scale epidemic.There are many so-called
`childhood diseases'such as chicken pox and measles where exposure and recovery
provides the host with extended immunity.The duration of this immunity varies
from one virus strain to another because even though a host may produce anti-
bodies to ¯ght further infection,the virus may mutate into a more resilient strain.
This experiment does not examine the mutation rate of a virus (it assumes that
there is only one strain) rather it looks at how immunity can be implemented in
a CA model.
7.1.1 World setup
The world in this experiment is identical to the uniform world described in the
control scenario described in Chapter 6.Each cell has a carrying capacity of
200 and an initial population of 100 susceptibles.Being free of infectives,any
virus outbreaks in this world will be due to the non-zero spontaneous infection
parameter.
27
7.1.2 Parameter settings
The only parameter that is di®erent from its default value is the resusceptibility
probability.This probability is the chance that a recovered host will revert back
into a susceptible one.This experiment looks at resusceptibilities of 0,0.2,and
0.8 over 600 time steps.Table 7.1 gives the values of the other epidemic spread
parameters.
Parameter
Value
Infection radius
1
Movement radius
1
Immigration rate
0.01
Birth rate
0.02
Natural death rate
0.01
Virus morbidity
0.05
Spontaneous infection rate
0.0001
Vectored infection rate
0.2
Contact infection rate
0.4
Recovery rate
0.1
Resusceptible rate
varied
Movement probability
0.001
Table 7.1:Parameter values for the residual immunity experiment.Most of the
numbers here are the same as the defaults de¯ned in Chapter 6,except that the
resusceptibility rate will be varied.
7.1.3 Results
The population dynamics for each of the di®erent resusceptibility values are
shown in Figures 7.1 and 7.2.In Figure 7.1 notice that for any increase in re-
susceptibility the height of the peak,the proportion of infectives at steady state,
and the time taken to reach steady state all increase as well.
The graphs of Figure 7.2 show that although there is an initial increase for
all three values of resusceptibility,at steady state the plots settle at di®erent
values.For smaller a resusceptibility parameter,the ¯nal steady state population
is smaller also;the population does not appear to be capable of recovering to its
original size.
28
0
100
200
300
400
500
600
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
Residual immunity
Time (epochs)
Proportion of infectives
p = 0
p = 0.2
p = 0.8
Figure 7.1:The proportion of total population that is infective over 600 time
steps.The probability of resusceptibility is set to p = 0,p = 0:2,and p = 0:8.As
the resusceptibility probability is increased,the duration of the epidemic (signi-
¯ed by the width of the peak) and the level of endemic infection also increases.
0
100
200
300
400
500
600
0
2
4
6
8
10
12
x 10
5
Residual immunity
Time (epochs)
Total population
p = 0
p = 0.2
p = 0.8
Figure 7.2:The °uctuations in total population for the same world as in Fig-
ure 7.1.As the resusceptibility probability is increased from 0 to 0.2 to 0.8,there
is a reduction in the total steady-state population.Note that each graph reaches
approximately the same maximum,but have markedly di®erent populations when
the virus has decayed to its endemic level.
29
7.1.4 Discussion
A comparison of Figure 6.1 with Figure 7.1 shows that although their shapes
are similar,increases in the resusceptibility parameter result in a virus that is
able to maintain itself in greater numbers and for a longer period of time.This
is the expected result because resusceptibility controls how rapidly the pool of
susceptible hosts is replenished in the absence of external host in°uxes;higher
resusceptibility means more potential hosts for the virus.
This experiment,apart from exhibiting some of the epidemic behaviour ob-
served in nature,shows that we can use the resusceptibility parameter to directly
control residual immunity of a host population without needing to create a new
parameter that controls the duration of a host's immunity.
For the purposes of modelling,being able to use a single global parameter to
control the immunity of every host reduces the model's complexity.Additionally,
the resusceptibility parameter is stochastic so the bene¯t is twofold:each host
does not need to carry memory of how long it has been immune and there is no
need to actively introduce noise into any duration parameters.This bodes well
for a CA model because the cells in a CA model use localized information and
global rules to produce complex overall behaviour.
7.2 High virulence
Virulence,also called morbidity,is a measure of how rapidly a virus kills its host.
Viruses are unable to replicate on their own,so it is in their own survival interests
that they do not kill their host.However,many of the mechanisms for host to
host spread,such as coughing,rely upon symptomatic disease.Consequently,a
balance must be made to maximize the spread to new hosts through disease but
keeping the current host alive.This experiment investigates if a CA model can
accurately emulate highly virulent viruses.
7.2.1 World setup
Apart from saying the world setup is the same as the one in the previous exper-
iment,I will not discuss it further as it is described in Chapter 6.Note however,
that starting without any infective sources means that all infections must start
from spontaneous outbreaks.
30
7.2.2 Parameter settings
To isolate the e®ect of the virus morbidity parameter on epidemic spread,all of
the parameters described in Chapter 5 are set to their default values.For this
experiment virus morbidity is set to 0.5:it is expected that at each time step
half of the infective individuals will die.This represents a highly virulent virus
such as some strains of Ebola found in the developing countries of Africa [24].
The epidemic spread parameters are summarized in Table 7.2.
Parameter
Value
Infection radius
1
Movement radius
1
Immigration rate
0.01
Birth rate
0.02
Natural death rate
0.01
Virus morbidity
0.5
Spontaneous infection rate
0.0001
Vectored infection rate
0.2
Contact infection rate
0.4
Recovery rate
0.1
Resusceptible rate
varied
Movement probability
0.001
Table 7.2:Parameter values for the virus morbidity experiment.A virus mor-
bidity value of 0.5 represents a virus that on average kills within two epochs of
infection.
7.2.3 Results
The graph in Figure 7.3 shows that after a rapid increase in the proportion of
infectives (numbers reach approximately 8% of the total population),there is a
similarly steep decay back down to approximately 1.5%.From about t = 100
onwards,the graph becomes very jagged where values °uctuate between 2.5%and
4.5%.There are distinct troughs and peaks in the graph where a local maximum
is followed by a local minimum less than 10 time steps later.This jagged pattern
is superposed with a larger pattern with a period of approximately 50 time steps,
seen more clearly after t = 500.However,there does not appear to be any clear
periodicity elsewhere.
31
0
100
200
300
400
500
600
700
800
0
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
Highly virulent virus
Time (epochs)
Proportion of infectives
Figure 7.3:A virus that has high virulence (or morbidity) will kill a high propor-
tion of infective hosts before they can recover.This plot shows a virus infecting
the control population described in Chapter 6 with a resusceptibility of 0.001 and
morbidity of 0.5 executed over 800 time steps.
0
100
200
300
400
500
600
0
1
2
3
4
5
6
7
8
9
10
11
x 10
5
Highly virulent virus
Time (epochs)
Total Population
Figure 7.4:A similar plot to Figure 7.3 except that it shows total population not
just infectives.Once again there appears to be an underlying oscillatory nature
to the curve.
32
Figure 7.4 does not show the same jagged characteristics of Figure 7.3.How-
ever,after the initial population decay the graph seems to settle into an oscilla-
tory pattern.To show that this oscillatory pattern is not a result of coincidence,
I have repeated this same experiment but with resusceptibility probabilities of
0,0.2,and 0.8.The results are presented in Figure 7.5.Notice that at least
qualitatively,each set of three graphs looks very similar.
7.2.4 Discussion
A highly virulent virus that kills its host quickly is unlikely to successfully spread
to a new host [3].This e®ect is shown by the much reduced numbers of Figure 7.3
in relation to the control scenario of Chapter 6.The jagged nature of this graph,
not present in the control scenario,is due to the sudden death of many infec-
tives during a particular time step.This periodicity has been observed in closed
populations such as those on found on the Faroe Islands [20].The quasi-periodic
nature of the curve results from the regular injection of new susceptibles (new-
born or immigrant) that provide the epidemic with new hosts only to later have
these new hosts die from the infection a few time steps later.
The similarities in the graphs of Figure 7.5 suggest that changing the residual
immunity of hosts has no signi¯cant e®ect when dealing with a highly virulent
virus.The e®ect of the immunity parameter appears to be negligible compared
with the in°uence of the virulence parameter.This makes sense because the virus
kills so swiftly that recoveries are few.In nature it is probably di±cult to ¯nd
a virus that kills rapidly and yet provide the survivors with lifelong immunity.
This is because both of these characteristics work against virus spread and would
not be deemed advantageous mutations during natural selection.
This experiment shows how some emergent behaviours of epidemics such as
periodicity have been captured in a simple CA model.Speci¯cally,this period-
icity is encapsulated into a single parameter:the virus morbidity.Although it is
the ability of CA to incorporate space into epidemic models that is under ques-
tion,showing that a CA model can reproduce the statistical results as observed
in nature provides us with some con¯dence in CA as a modelling approach in
general.
33
0
100
200
300
400
500
600
0
0.05
0.1
Highly virulent virus
p = 0
0
100
200
300
400
500
600
0
0.05
0.1
Proportion of infectives
p = 0.2
0
100
200
300
400
500
600
0
0.05
0.1
Time (epochs)
p = 0.8
0
100
200
300
400
500
600
0
5
10
x 10
5
Highly virulent virus
p = 0
0
100
200
300
400
500
600
0
5
10
x 10
5
Total population
p = 0.2
0
100
200
300
400
500
600
0
5
10
x 10
5
Time (epochs)
p = 0.8
Figure 7.5:Plots of the same experiment as Figures 7.3 and 7.4,but with resus-
ceptibility values of p = 0,p = 0;2,and p = 0:8.It is interesting to observe the
similarity in shape to each of the curves regardless of the variation in resuscepti-
bility.
34
CHAPTER 8
Spatial Experiments
The focus of this dissertation is on the spatial aspects of epidemics,particularly
how spatial parameters a®ect the way they spread.This chapter contains a
series of experiments that show how a CA model can emulate several emergent
behaviours of an epidemic yet use only one rule set.The results of the experiments
in this chapter provide a basis for deciding the suitability of CA to epidemic
spread modelling.
8.1 Varied population density
A virus relies upon a supply of susceptible hosts to survive and replicate.But
before a virus can infect a new host it must come in contact with it.An indi-
cator to whether an outbreak will spread is the basic reproductive number,R
0
,
discussed in Chapter 2.According to Equation 2.2,R
0
is directly proportional to
the number of contacts per unit time which itself is proportional to the local pop-
ulation density.Therefore,in regions of high local population density we expect
an epidemic to spread more rapidly than in low density regions.This experiment
examines the e®ect of varied local population density and how CA rules can be
used to emulate those e®ects.
8.1.1 World setup
This experiment is broken down into two sub-experiments.To juxtapose the
e®ect of local population density on virus infection velocity I have used two
synthetic landscapes { one with a negative density gradient from top to bottom,
and another with a positive density gradient from top to bottom.These are
shown in Figures 8.1 and 8.2.Those ¯gures also show the initial`line source'of
infectives at the top of the grid.It is assumed that during the early stages of
the experiment the epidemic spreads through the top of the grid and during the
later stages spreads through the lower part of the grid.
35
Figure 8.1:The population distribution for the high to low density scenario.
On the left is the population distribution of susceptibles (S),infectives (I),and
recovered (R).The shading indicates a negative population density gradient from
top to bottom.The ¯gure on the right shows the line source of infectives denoted
by the dark row of cells at the top of the grid.These infectives start o® the
epidemic spread.
Figure 8.2:The population distribution for the low to high density scenario.Its
features are similar to Figure 8.1,except that the gradient direction is reversed.
36
8.1.2 Parameter settings
To isolate the e®ect of local susceptible population density on the spread of an
epidemic,most of the other parameters in this experiment are set to zero.This
is to remove any skew that may arise from including several factors at once.
The vectored and contact infection parameters are set to unity to accelerate the
spread of the virus;doing this should not adversely a®ect the results because it
is not the absolute value of the infection velocity that is of concern,but how it
changes with respect to population density.The recovery probability is set to
zero so that the epidemic spreads in one direction only { that is,velocity cannot
be negative.
Additional to using a ramped population density,this experiment will look
at the e®ects of increasing the interaction neighbourhood size of the local pop-
ulation density,as de¯ned in Equation 8.1.Consequently,changing the size of
the interaction neighbourhood will also change the local population density.The
values of the other parameters are summarized in Table 8.1.
local population density =
neighbourhood population
neighbourhood capacity
(8.1)
Parameter
Value
Infection radius
varied
Movement radius
0.0
Immigration rate
0.0
Birth rate
0.0
Natural death rate
0.0
Virus morbidity
0.0
Spontaneous infection rate
0.0
Vectored infection rate
1.0
Contact infection rate
1.0
Recovery rate
0.0
Resusceptible rate
0.0
Movement probability
0.0
Table 8.1:The parameter values for the population density experiment.Most
of the parameters are set to zero except for the ones that relate directly to
population density.For this experiment,the infection radius is varied and its
e®ect examined.Note that both the high to low density scenario and the low to
high density scenario use the same parameter values.
37
8.1.3 Results
Each of the plots in Figures 8.3 and 8.4 show the proportion of the total pop-
ulation that is infective over the course of 100 epochs or iterations.There are
two features of note for each plot of each ¯gure:the plot's gradient and the time
before the population becomes saturated with infectives.
For the high to lowdensity scenario in Figure 8.3,we can see that the gradients
of the plots start o® much steeper than they ¯nish.That is,as the population
density decreases so does the gradient.A similar e®ect can be seen for the
low to high density scenario in Figure 8.4.Each of the curves starts with a
shallow gradient,or small infection velocity,and as the epidemic reaches the
higher density regions,the gradients steepen.
When the infection radius is increased the population becomes saturated with
infectives more rapidly.That is,the`wave'of infection reaches the bottom of the
grid in less time.
0
10
20
30
40
50
60
70
80
90
100
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
High density to low density
Time (epochs)
Proportion of infectives
r = 1
r = 2
r = 5
r = 10
Figure 8.3:E®ect of population density on epidemic spread velocity.For the
earlier time steps,we see the gradient of the graph is much steeper than for the
latter time steps.Each of the graphs corresponds to a di®erent infection radius,
r.Notice that for increasing infection radius the gradients become steeper.
38
0
10
20
30
40
50
60
70
80
90
100
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Low density to high density
Time (epochs)
Proportion of infectives
r = 1
r = 2
r = 5
r = 10
Figure 8.4:Similar scenario to Figure 8.3,but for a di®erent CA starting state.
This plot appears to be opposite to the previous scenario:here the transient
¯nishes much steeper than it starts.Each of the graphs corresponds to a di®erent
infection radius,r.
8.1.4 Discussion
The infection velocity of the epidemic is proportional to the gradient of the plots
in Figures 8.3 and 8.4.As expected,as the epidemic spreads through a high
density region,the infection velocity is greater than for low density regions.In
nature,this e®ect is attributed to the increased number of contacts between an
infective host and other susceptible hosts,as shown in Equation 2.2.
As the size of the infection interaction neighbourhood is increased there is a
narrowing e®ect on the graph:the virus propagates through the population more
rapidly,but still exhibits the velocity changes noted earlier.This acceleration in
infection is mostly due to the larger`reach'a contagion has to infect hosts in
other cells.For example,viruses that are airborne or spread by birds would have
a high`reach'or transmissibility.
From a modelling perspective,the two curves of Figures 8.3 and 8.4 are ex-
actly as predicted by epidemic theory { that is,high densities promote rapid
spread.This experiment has shown that a CA model can take into account the
local population pro¯le and produce the appropriate epidemic spread velocity.
This partially justi¯es CA as a good epidemic modelling paradigm in that we
have shown that the ruleset I devised for my model produces the desired results.
39
However,as yet we have no assurances that these rules will work for more complex
landscapes.
8.2 Host immigration
Whilst the experiment that examined variations in population density had the
density pro¯le remain static for the duration of the experiment,this experiment
looks at the e®ect of dynamic variations in population density and composi-
tion.Of particular interest are the variations due to host immigration.Although
not purely a spatial factor,immigration a®ects the population density and dis-
tribution over the landscape,indirectly in°uencing the spatial behaviour of an
epidemic.
8.2.1 World setup
This experiment uses the exact same initial setup as the control scenario described
in Chapter 6 where the population is uniformly distributed over the landscape.
Each cell starts with a population of 100 susceptibles and has room to ¯t 100
more SIR hosts.
8.2.2 Parameter settings
Once again this experiment is broken down into two sub-experiments.However,
instead of using di®erent starting landscapes,this experiment uses two di®erent
sets of spread parameters.Most of the parameters are initialized to the default
values de¯ned in Chapter 6 except of course for the immigration rate.To simulate
a closed population,the immigration rate is set to zero { there are no incoming
hosts from outside the world.For an open population,this parameter is set to
0.05,meaning that approximately 5% of cells will receive susceptible hosts from
an external source to ¯ll up some of its empty space.The values of the other
parameters are shown in Table 8.2.
8.2.3 Results
The results for the closed population are displayed in Figure 8.5.It shows the
familiar shape of an initial virus outbreak that decays back down to a smaller
steady state value.The simulation was run for 800 epochs and for four di®erent
40
Parameter
Value
Infection radius
1
Movement radius
1
Immigration rate (closed)
0.0
Immigration rate (open)
0.05
Birth rate
0.02
Natural death rate
0.01
Virus morbidity
0.05
Spontaneous infection rate
0.0001
Vectored infection rate
0.2
Contact infection rate
0.4
Recovery rate
0.1
Resusceptible rate
varied
Movement probability
0.001
Table 8.2:The epidemic spread parameters used in the host immigration ex-
periment.Most of the values here are the same as in Table 6.1 except for the
immigration probability which is 0 for the closed scenario and 0.05 for the open
scenario.Part of this experiment will look at the e®ects of residual immunity
controlled by the resusceptibility parameter.
values of resusceptibility:0,0.1,0.5,and 0.9.As the chance of resusceptibility
is increased there is a corresponding increase in the maximum ratio between
infectives and total population.
In each case,the proportion of infectives decays down to a value very close to
zero.However,it is interesting to note that although for values of resusceptibility
equaling 0.1,0.5,and 0.9 the infective population decayed to zero by t = 300,
when the resusceptibility was set to zero,the infective population did not reach
zero until much later.
The results for the open population shown in Figure 8.6 are much more con-
sistent.Although the simulation was run for 800 epochs,only the ¯rst 400 are
shown here as the infections rapidly reached steady state.The main feature to
note is that all of the plots reach a non-zero steady state value.Increases in
the resusceptibility parameter have resulted in a corresponding increase in the
proportion of infectives at steady state,as well has removing the distinct peak
that represents the epidemic.
41
0
100
200
300
400
500
600
700
800
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
Closed population
Time (epochs)
Proportion of infectives
p = 0
p = 0.1
p = 0.5
p = 0.9
Figure 8.5:A closed population does not allow outside individuals to enter the
world.From the graph above we see that after the initial outbreak peaks and
decays to zero,the infection does not become endemic.Each of the plots corre-
sponds to a resusceptibility probability of p = 0,p = 0:1,p = 0:5,and p = 0:9.
Notice that for increasing resusceptibility,the number of hosts that become in-
fected during the epidemic (the hump) increases dramatically.
0
50
100
150
200
250
300
350
400
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Open population
Time (epochs)
Proportion of infectives
p = 0
p = 0.1
p = 0.5
p = 0.9
Figure 8.6:An open population allows the free movement of individuals in and
out of the world.These graphs show a plateau structure as the resusceptibility
probability is increased from 0 to 0.9.In contrast with Figure 8.5,all of the plots
reach an endemic non-zero value.
42
8.2.4 Discussion
Historical evidence has shown that for a population in which individuals are
always coming and going there is a greater chance that an infection will become
endemic after an outbreak [21].Conversely,for isolated communities,it is usually
the case that after an initial outbreak,the immunity acquired by the survivors
limits any chance of a virus continuing to be active [20].The reason for the e®ect
mentioned by Scott and Duncan [21] is that in an open population,the virus can
rely upon a constant in°ux of new susceptible hosts.On the other hand,remote
populations such as those on islands or in isolated villages provide the virus with
no new hosts and usually results in the virus'demise.
The results presented in this section partially portray the e®ect of immigra-
tion on endemic infection,but other anomalies make it hard to draw any other
conclusions.Quantitatively,the results are very extreme:an endemic proportion
of 70% is highly unlikely { such large ¯gures clearly show the limitations an un-
calibrated model.Qualitatively,we are still able to see the expected increase in
endemic infection with an increase in resusceptibility,so I believe there is still
enough evidence to justify the development of a composite CA epidemic model.
8.3 Corridors of spread
This experimental scenario is the ¯rst of two that will use my CA model to try
and visualize how an epidemic would spread through a more complex virtual
landscape than that seen earlier.
8.3.1 World setup
This experiment tries to reproduce what might be a real life landscape with
imaginary town centres and transport links between them.I use`towns'and
`roads'because they are human in°uenced and attract high population densities.
The scenario is not restricted to cultural features;other natural features such as
rivers might lead to development along their shores.Conversely,other geograph-
ical features such as mountain ranges or swamps will limit development.What
is important is their directed and linear shape rather than the wide sweeping
population densities used in previous experimental scenarios.Figure 8.7 shows
the virtual landscape used in this experiment.
Each cell has a carrying capacity of 1000,with three`towns'already at this
maximum.Two of these towns,the north-west one and the south-east one,
43
Figure 8.7:The synthetic landscape that I devised to show how an epidemic is
likely to spread along features of high density before spreading across open space.
There are three points of high population density { two of which are connected by
a transport link which has its own settlement developed on either side of it.The
north-west and south-east`towns'both have 100 infectives and 900 susceptibles.
start with a 1:9 infective to susceptible ratio.The`transport links'comprise a
three cell wide bar of susceptibles.Cells on the central axis of this bar have
an initial susceptible population of 100,whilst the cells on either side have an
initial susceptible population of 75.The rest of the landscape comprises cells
with 10 susceptibles in them.The town in the south-west corner contains 1000
susceptibles and no infectives.
8.3.2 Parameter settings
The parameter settings for this scenario,because I am looking at the directed
spread of the epidemic with respect to population density,are the same as for
the high to low,and low to high scenarios discussed earlier.Table 8.1 lists
the parameter values:essentially all of the parameters are zero except for the
infection radius,contact infection probability,and vectored infection probability
which are all equal to one.
8.3.3 Results
The results of this scenario are presented in the lag map shown in Figure 8.8.The
lag map is basically a series of snapshots taken at t = 0;20;40;60;80;100;200;300.
44
Each cell is represented by a coloured square:red squares contain at least one
infective host and white squares contain only susceptible and recovered hosts.
Figure 8.8 shows the tendency of the epidemic to follow the lines of population
density to produce the`fuzzy cross'pattern.
After 300 epochs,the top left outbreak has reached all four edges of the map
but the bottom right outbreak is yet to reach any.Notice that the epidemic
spreads outward along the arms or`roads'before ¯lling up the space between the
roads.
8.3.4 Discussion
This particular scenario provides some evidence of the power of CA models for
visualizing an epidemic;di®erential equation models provide no such capability.
In Figure 8.8 it becomes obvious to see where the infective towns were,where the
densest populations are situated,as well as identify where the infection velocity
is greatest.Features such as leap frogging or in¯lling where the epidemic appears
to skip regions of land before coming back to ¯ll in the space can be seen in the
lag map.This behaviour was documented in the foot and mouth disease epidemic
in Great Britain in 2002 [14].
Knowing beforehand the probable direction that an epidemic will take might
help public health o±cials e±ciently direct containment measures and medical
services to deal with infections.
8.4 Barriers to spread
The last experimental scenario in this chapter tries to show how a CA model
can simulate the e®ects of erecting barriers to slow or stop virus spread.As seen
in the foot and mouth disease epidemic in Great Britain during 2001,a key to
slowing down disease spread is restricting movement.Other measures included
the culling of livestock or the inoculation of livestock.I try to reproduce these
measures in this scenario by varying the initial state of the CA lattice.
8.4.1 World setup
Figure 8.9 shows the starting distribution with two squares that represent`cattle'
farms.The top left farm has a four square wide barrier surrounding it,whereas
the bottom left farm has a one square wide barrier surrounding it.Barriers are
45
Figure 8.8:A lag map showing the state of the epidemic at t =
0;20;40;60;80;100;200;300.Notice that the outbreak to the north-west is able
to cover a greater distance than the outbreak in the south-east because it has
access to the road link and the population associated with that link.Notice that
the spread from t = 20 in the top left of the map appears asymmetric.This is
probably an artifact of the stochastic nature of this model.
46
implemented as cells with zero carrying capacity.In Figure 8.9,barriers are
represented by black squares and host occupation is represented by blue shading.
Figure 8.9:In this scenario there are two sources of infection,each with a di®erent
sized barrier around it restricting host movement and providing no hosts for any
viruses that try to cross.
8.4.2 Parameter settings
Apart from increasing the infection radius to two,there are no additional param-
eter changes for this scenario { the default values used in the control scenario are
used again.Having an interaction radius of two means that the barrier around
the infectives in the bottom right corner will still be able to disperse.
8.4.3 Results
The resultant epidemic spread is shown in the lag map of Figure 8.10.Notice
that the virus is able to elude containment in the bottom right`farm'.Despite
the farm itself becoming free of infectives,the neighbouring country side has
become infective.The upper left farm is no longer contaminated and neither is
its surrounding hinterland.
47
Figure 8.10:This lag map shows that bu®er zones that are too narrow provide
no resistance to the spread of a virus.
48
8.4.4 Discussion
Much like the previous experiment,this scenario provides an opportunity to
demonstrate the visualization capabilities of a graphical CA model.There are no
new epidemiological conclusions to be drawn from either of the last two experi-
ments { they show exactly what our intuition would suggest.What is important
is that other statistical models do not appear to illustrate such behaviours very
well.As far as sample data is concerned,very little is available about the spatial
behaviour of epidemics;much of epidemiology looks at statistical data and try-