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-

## Comments 0

Log in to post a comment