Abstract
The world without a disease is a dream of any human being. The disease spread if not controlled could cause an epidemic situation to spread and lead to pandemic. To control an epidemic we need to understand the nature of its spread and the epidemic spread model helps us in achieving this. Here we propose an epidemic spread model which considers not only the current infective population around the population but also the infective population which remain from the previous generations for computing the next generation infected individuals. A pushdown cellular automata model which is an enhanced version of cellular automata by adding a stack component is being used to model the epidemic spread and the model is validated by the real time data of H1N1 epidemic in Abu Dhabi.
1. Introduction
Computational models in epidemics provide insight into dynamics of the disease spread across a geographical region. Given a small amount of relevant real time data the model would give us after certain amount of time what would be the epidemic situation across the geographical region. Cellular automata (CA) models are considered to be very handy and efficient in handling the real time simulation problems in epidemiology. The reason behind the simplicity of this CA model is its ability to attain the global behavior from the local behavior by the interaction of its cells [1–3].
The power of CA has been utilized to solve wide variety of problems like solving Weyl, Dirac, and Maxwell equations [4], diffusion equation [5], and Poisson equations [6]. Epidemic spread modeling through cellular automata has been an importance to many of the researchers. An efficient epidemic spread model through CA is given by Hoya White et al. [7]; in this model the population is assumed to be constant and the rules for the CA are considered to be static. The movement of population and its effects during an epidemic spread are given by Sirakoulis et al. [8]. The population when grouped under patches and their movement during an epidemic spread are discussed in our earlier work [9]. The spatial pattern and its dynamics along with noise in epidemics are being provided by Sun et al. [10]. Stochastic model for epidemic spread with quarantine and vaccination strategies have been provided by Wang Jeffrey [11].
In all these works discussed above, the susceptible population is infected by the infected population available in the current time step. There are no considerations for the populations who are infected at the previous time steps. If we consider the infected population at current time as , then the population which was infected at the previous time steps is considered as , , and so forth. This work largely focusses on how the latent infected population, that is, the ones which are infected at , , has an effect in infecting the current set of susceptible population in a cell. This relates to a question: how could we segregate or handle the infected populations at or ? The solution is to append a stack to each cell of the cellular automaton and call the modified CA as pushdown cellular automaton or PCA [12]. The idea is to hold the infective population at different time steps at different levels of stack and use them for computation for the next time step accordingly.
In this research paper we begin introducing the cellular automata in Section 2.1, then in Section 2.2 we introduce the usage of cellular automata in epidemic spread simulation and the cellular automata with stack configuration would be handled in Section 2.3 followed by various modifications in stack configured cellular automata model in the next sections. The results and discussions are dealt with in Section 3 followed by conclusion and the future scope of the model.
2. Materials and Methods
2.1. Cellular Automata
Cellular automata is a universal Turing machine which is capable of solving many complex problems ranging from simulating forest fire spread, landmine detection, robotic movements, and object avoidance to understanding the spread of epidemics [15]. It consists of a grid of cells, which is the basic building component of a cellular automaton. Each cell could be in one of the finite amount of states. The next state of each cell depends on the state of the current cell and states of the neighboring cells. The neighborhood configuration can be many but mostly Von-Neumann neighborhood or Moore neighborhood is used depending upon which of the surrounding cells we consider for computing the next state of the current cell.
2.2. Cellular Automata for Epidemic Spread
A cellular automaton can be configured for simulating an epidemic spread. The following assumptions would be applied for configuring the cellular automata.
A cell could be in one of the following finite number of states: susceptible, infected, and recovered/removed (dead). We assume the epidemic model to be of SIR type with equal birth and death rate so the population as a whole remains constant.
The neighborhood configuration would be of Moore’s neighborhood which considers the top, down, left, right, and the diagonal neighbors as well.
The population count in each cell is assumed to be constant and no movement between the cells is permitted. This assumption is for the very basic model whereas in real time the movement of population between the cells is permitted.
The initial configuration of the cellular automata and the cellular status word for each cell is as described in our work [9].
The cell status word (CSW) (see Table 1) has five parameters, namely, Cell_id which is the current cell denoted by , CPV (critical population value) which is the fraction of infected population over total population in the cell, value which takes Boolean value 0/1: 0 means the cell is free of infective population; 1 means some percentage of cell is infected so we consider that the cell is infected, value which will also take either 0 or 1 depending upon whether the population inside the cell is recovered from the disease or not, and value which is the movement flag for a particular cell used to control the movement of the population from one cell to another. If the value is 0 it means no movement possible from the cell.
The critical population value is computed by where IPC represents the infected population count in the cell and TPC represents the total population count in the cell. The critical population value for a cell at the next time step would be the function of the CPVs of the current cell at time and the CPVs of the neighboring cells; that is,
Size of the grid or cellular array depends on the user and the type of geographical area under consideration. The bigger the cellular array is the longer it takes for the simulation to run.
2.3. Pushdown Cellular Automata for Time Variant Infective Population
A modified cellular automata with added stack to each cell is termed as pushdown cellular automata. The concept of time variant infective population comes to importance in the case of population movement across the cells. In each time step a certain percentage of population moves from one cell to another based on the configuration of the cellular array. During the movement a certain percentage of infective, susceptible, and recovered population moves from one cell to another and the count is considered to be random across all the three sets of populations. So some amount of infective population moves from current cell to another, and some remains within the current cell. So there can be two categories of infective population in the current cell, that is, the ones which remain from the previous time step and the ones which could come from another cell to the current one. The infective population which remain from the previous time step would be termed as and the infective population which come into the current cell would be termed as . The degree of infection caused by the newly infected infective population and the infection caused by the older ones would be different. So the cell update rule will be modified to accommodate the effect of these current infective population and the older ones. Here in this research work we consider three time steps , , and , respectively. So the stack size would be 3 (plus one to hold the bottom symbol $ which indicates the bottom of the stack).
2.4. Modification of CSW with Stack
The cell status word is modified to hold the stack contents for storing the critical population values of the previous time steps. The modified CSW will look like Table 2.
CPV part is replaced by SP (stack pointer) which points to the top of the stack attached to the current cell. Default value/initial value of the SP would be $ which indicates the bottom symbol of the stack and also the emptiness of the stack. The structure of the CSW during the three time steps is given by Tables 3, 4, and 5.
In every time step the current cell is updated based on the status of itself and the status of the neighboring cells.
2.5. Configuration of PDCA
The PDCA has the following components.
Finite set of states : neighborhood configuration would be Moore neighborhood (all cells surrounding the cell to be considered for update). The neighborhood status is divided into three categories based on the number of infected cells out of the eight neighborhood cells. The status is minimum denoted by “m” when the number of the infected cells surrounding the current cell is less than 4. The status is half denoted by “h” when the number of the infected cells surrounding the current cell is either 4 or 5. The status is full denoted by “f” when the number of the infected cells surrounding the current cell is more than 5. Given below is the various neighborhood status. White indicates susceptible state, red indicates infected state, and green indicates removed or recovered state (Figure 1).
Stack of size three is attached to each CSW of the cell with the following stack symbols:$:special symbol to indicate the bottom of the stack,:number of infective remaining from the previous time step,:number of infective population remaining from the previous to previous time step.
Transition function () is given by the following.
Initially all the cells are assumed to be in the susceptible state (S).
Each transition function will be of the following form: (current state of the cell, state of the neighborhood, stack content before transition) (state of the cell after transition, stack content after transition).(i) (pushing the $ on to the stack before starting the computation).
The biological meaning of $ can be considered as the residual infective population whose impact is not severe and remains in every cell:(ii).
When the cell is in susceptible state and the neighborhood is in the minimum status with the stack top being $ then the cell remains in susceptible state and the neighborhood also in minimum status and the number of infective population that remains in the cell would be pushed as on top of the stack. The biological meaning would be that current location is not infective (which means the number of infective population is less) and surrounding is also having less number of infective population; then in the next time step the cell is not infected and remains in susceptible state:(iii).
When the cell is in susceptible state and the neighborhood is in the half status with the stack top being $ then the cell remains in susceptible state and the neighborhood also in half status and the number of infective population that remains in the cell would be pushed as on top of the stack. The other choice of this state to the next time step would be changed to infective state, and the neighborhood status changes from half to full and the would be pushed on top of the stack. The choice would be selected at random as in real time cases this randomness helps us in introducing the noise factor into the system which is very appropriate since the epidemic situation is a nondeterministic one. As the number of infective population increases around the surroundings the current cell can change from susceptible state to infective state:(iv).
When the cell is in susceptible state and the neighborhood is in full status with stack top being $ then the cell changes its state to infective with the neighborhood remaining in full status and is pushed onto the stack. Likewise the remaining transition functions for various states of the PDCA are(i),(ii),(iii),(iv),(v),(vi).
2.5.1. Assumptions Regarding the Stack
The following assumptions were taken into account with regard to the stack.(i)Here we consider only two time steps before , so the stack size becomes three along with the bottom symbol $. If we increase the time steps to unlimited amount, this will lead to the proportionate amount of increase in the size of stack (memory) which will definitely affect the efficiency of the computation.(ii)There is no movement of the population which is inside the stack.
2.6. Setting Up the Model
The first step in designing the cellular automata model is defining the size of the grid. The size of the grid would be inversely proportional to the speed of computation [8, 9]. The grid would be taken with the background of the city of Abu Dhabi since we would be simulating our results with the data of H1N1 database of Abu Dhabi (Figure 2).
After mapping the grid with the map of the city the boundaries are drawn for the simulation area (Figure 3). We presume that the epidemic would not spread beyond these boundaries since it is surrounded by sea. Overlapping top and bottom rows and left and right rows of the grid are not a feasible assumption as the city is surrounded by the sea and three sides.
Initial conditions are as follows:(i)cellular array size defined by the user,(ii)number of persons in each cell 100,(iii) value 0 for all cells,(iv) value 0 for all cells,(v)infection period defined by the user as per the type of disease,(vi) value 0 no movement in initial configuration.
The critical population value (CPV) would be now calculated by the following way:
The CPV of a cell in the next time step would be calculated in the same way after transition into a new state based on the transition function.
2.7. Homogeneous PDCA without Movement of Population
In this configuration the number of persons in each cell is kept the same and the movement of population between the cells is also restricted by giving the value in CSW to 0. This configuration is not the right one to model the real time environment but necessary to understand the basic parameters during the simulation.
2.7.1. The Algorithm for Homogeneous PDCA (Equal Population in Each Cell with No Movement)
Step 0. Get the input from the user: the value of tr (time period for recovery), tg (number of generations), and ti (time period for immunity). Block the boundary cells so as to control the infection spread.
Step 1. Initialize the generation counter to 0 and push $ on to the stack for all the cells.
Step 2. Initialize all the cells with equal number of populations.
Step 3. Initialize the value of the CSW of all cells to 0 to control the movement of population.
Step 4. Infect the desired cell.
Step 5. Increment the generation counter by 1.
Step 6. Find out the region of cells which are infected and the region of cells which are not infected by checking the value and value.
Step 7. Based on the top symbol of the stack and the status of the neighborhood change the state of the current cell and modify the stack according to appropriate transition function.
Step 8. If value is 1 then go to 9.
Else calculate the CPV value based on the cellular automata update rule to find the fraction of infected population.
Step 9. If the number of generations is less than tg then go to step 4.
Step 10. Stop.
2.8. Heterogeneous Population with Movement
In this configuration, the population in each cell varies and moreover the movement of population between the cells is also possible. The movement of population from a particular cell is possible only if the cell’s CSW-M value is 1; otherwise movement cannot be done from that cell. If we consider a geographical area (city or town) to which the grid is mapped then all those cells which are in line with the streets or bridges of the city are enabled for movement. After deciding which cells allow the movement of population, we need to identify the amount of population that would be moving from the cell and the time of movement.
The amount of population which can move from one cell to another cell would be 5% to 7% and we can deploy two techniques to trigger the movement of the population. First one is that after some random amount of time or when the infected fraction of population reaches more than 15% of the total population the movement is done. The second choice would be more natural in the sense that the people tend to move from one place to another when there would be more infection around.
2.8.1. Algorithm for Heterogeneous Population with Movement
Step 0. Get the input from the user: the value of tr (time period for recovery), tg (number of generations), and ti (time period for immunity).
Step 1. Initialize the generation counter to 0 and set enable field of CSW to 1 for all those cells except the boundary cells.
Step 2. Initialize all the enabled cells with variable populations.
Step 3. Initialize the value of the CSW of all cells which are in the pathway to 1 in order to enable the movement.
Step 4. Push $ on the stack and infect the desired cells.
Step 5. Increment the generation counter by 1.
Step 6. Find out the region of enabled cells which are infected and the region of cells which are not infected by checking the value and value.
Step 6.1. Based on the top symbol of the stack and the status of the neighborhood change the state of the current cell and modify the stack according to appropriate transition function.
Step 7. Check if the infected percentage is greater than 15%. If yes then go to step 8; else go to step 9.
Step 8. Check whether the cell is movement enabled or not by checking value to be 1. If yes then go to step 8.1; else go to 9.
Step 8.1. If the cell is within the boundary limits transfer the 0 to 5% amount of population taken at random to the next cell in the neighborhood to the direction in which there is a pathway available.
Step 9. If value is 1 then go to 10.
Else calculate the CPV value by (3) based on the cellular automata update rule to find the fraction of infected population.
Step 10. If the number of generations is less than tg then go to step 4.
Step 11. Stop.
The PDCA could be slightly modified for considering the movement factor. The modification would be in terms of the neighborhood configuration; instead of using the traditional Moore neighborhood we can use the extended neighborhood concept so as to realistically model the situation. The extended neighborhood combines the traditional Moore neighborhood and randomly selects the cells across the neighborhood cells.
As we extend the neighborhood to broader radius say to 2 or 3 then the region of neighborhood increases for 8 to 16 or 24 cells, so the more the radius the less the probability of contact and infection. So the radius is confined to 1.
If the current cell is , refers to the Moore neighborhood and refers to the cells that are chosen randomly; then the extended neighborhood of and accordingly the cell update rule will also change for the extended neighborhood [16]:
This random choosing of cells in the neighborhood part and the random selection of individual for the movement give us the source for the noise addition in the model as the model we are considering is of nature which is not a deterministic one. The stack rules will remain the same and the transition functions based on stack will also be maintained similar to the previous one.
2.9. Important Parameters in the Epidemic Spread
Susceptible population denoted by is the population which is not infected but prone to infection when it is around the infected set of population. In this simulation environment we consider all the population initially to be susceptible. There may be many cases where certain population can be immune to the disease under consideration, so they may not be part of the susceptible set.
Infected population denoted by is the population which is actually infected by the disease under consideration. Certain amount of individuals from the susceptible set is infected by the infection rate denoted by “.” Here in our model we consider three different types of infective population; those are the ones present in the cell at current time , the ones which remain from the previous time step , and the ones from time step. The current ones are considered in the cell and in the next time step it is pushed onto the stack as . Since we have three different types of infective population, the infection rate also would be different for them; the infection rate will increase for the infected population if not treated or vaccinated.
We define three different infection rates for the three sets of infected population. They are , , and for , , and , respectively. We can define the relation between the three types of infection rates by
The above equation is valid in case there is no treatment or vaccination. So the infection rate increases as the time progresses and the infected individuals would be capable of infecting more numbers of susceptible population:
The above equation is applicable when there is a treatment or vaccination strategy is underway to the infected set of population. Initially when the treatment starts the effect will not be visible but after certain time period the infection rate starts reducing. That is the reason why is greater than both and .
Recovered populations denoted by are those which are recovered from the infected population and become immune to the disease. The rate of recovery from the infected population is given by . The recovery rate would be directly proportional to the level of vaccination strategy used for the infected population.
Reproduction number denoted by which is the important component in the epidemic environment denotes the amount of individuals infected by a single infected individual in a cell placed along with fully susceptible population: where is the infection rate, is the birth and death rate considered to be equal which means that the amount of birth is equal to the amount of death, and represents the recovery rate. When is less than 1 then the system is in infection-free stable state. When is greater than 1 the system is in endemic stable state. Always the objective is to reduce the to less than one by increasing the recovery rate. The above given equation is the general one and for estimating in the cellular automata paradigm with various radius of neighborhood the solution is given by where is the average number of connections per cell, is the probability that susceptible cells become infected, and is the neighborhood radius.
2.10. Controlling the Spread of Epidemics
There are various strategies for controlling the epidemics; quarantine and vaccination stand out amongst them.
2.10.1. Quarantine
There are two terms in the world of epidemic which needs to be clarified before dealing with quarantine. The terms are isolation and quarantine; what is the difference between these two? The clarification rightly given by Center for Disease Control (CDC-USA) [17] is as follows.
Isolation. It applies to person who is known to be ill with a contagious disease.
Quarantine. It applies to those who have been exposed to a contagious disease but who may or may not become ill.
Quarantine strategy is implemented by blocking 9 cells with a boundary and stating that the population outside the cell is susceptible to infection. The probability of the infection would be where is the escape rate, which means number of persons who can go through the quarantine defense [11]. This could be more realistic because we consider no quarantine system is full proof. The quarantine is analyzed with 10% and 25% escape rate.
The area to be quarantined can be identified using the full status of the stack given in the figure (Figure 4). The stack full status indicates that more numbers of infective population with different infection rates are accumulated when compared to the other region.
2.10.2. Vaccination
The most common and effective technique in dealing with control of epidemics is the vaccination technique. There can be many methods of vaccination: one is to continuously vaccinate for particular amount of time and the other one is to vaccinate for a small amount of time in intervals which is known as pulse vaccination. Here continuous vaccination after the 30th generation is done for another 10 generations. The epidemic spread simulation techniques help to reach cost effective technique in implementing the vaccination. The areas where the disease spread is more have to be identified and vaccinated accordingly. The important concept in vaccination is the efficacy of the vaccine; vaccine efficacy is defined as the reduction in effect of disease among people who have been vaccinated compared to the effect in unvaccinated people [18]. The influenza vaccine efficacy rate is considered to be 85% [19]. The probability of infection spread in the vaccinated area is where is the efficacy rate. Randomly selected cells which have large number of infective population starting from the 30th generation are being vaccinated and the results have been found.
Simulation environment: simulation is done using MATLAB Ver. 2013 (64 bit). The processor used is Intel Core-i5 and RAM capacity is 8 GB.
3. Results and Discussion
3.1. Homogeneous Model
Homogeneous model is the model in which there is a fixed amount of population in every cell. A 50 × 50 cellular grid (optimal grid size is the tradeoff between time, size, and accuracy) is used to simulate the epidemic spread as per the available surveillance data of H1N1 infection in Abu Dhabi [13, 14]. The black colored blocks in the grid are the boundary which is traced out from the Abu Dhabi city map. We have infected certain cells randomly in the left side of the grid and started the simulation and the snapshots of each of the 30 generations are provided here which clearly depicts how the infection spreads out in the beginning and is erased out in the end. The last two parts of the figure provide information about the cells in which the augmented stack is full and empty. The epidemic spread is simulated using the data of H1N1 infection in Abu Dhabi [13, 14] over the period of twelve months with the parameters of ; ; ; . Average comes out to be 1.52, , and .
The average infection period of a cell as a whole is around 15 to 20 generations (Table 6).
Figure 6 shows the curve of SIR model which matches with the ideal epidemic SIR curve and the number of infective population peak during the 4th and 5th month, respectively. Figure 7 gives us the picture of how the number of infective population increases during the 4-5th month along with the increase in the number of stacks getting full by , , .
Figure 8 shows us the number of confirmed cases with the effect of no quarantine, quarantine with 25% escape rate, and quarantine with 10% escape rate. 10% escape rate is considered to be optimum. Figure 9 gives the picture of how the number of infective population decreases with the random vaccination done every 30 generations across the cells in the grid.
3.2. Heterogeneous Population with Movement
Population movement is enabled in those CSWs where the stacks were full which is identified from Figure 5(k) (Figure 10). Simulation is done for another 300 generations and the infection spread reach is calculated from the start location to the various identified locations in the grid (Table 7). This analysis helps us in planning the vaccination, quarantine, and evacuation procedure which is a very important phase in epidemic situation.
(a) Initial status
(b) Generation 1
(c) Generation 30
(d) Generation 60
(e) Generation 90
(f) Generation 120
(g) Generation 150
(h) Generation 180
(i) Generation 210
(j) Generation 240
(k) Generation 270
(l) Generation 300
(m) Generation 330
(n) Areas where stacks were full
(o) Areas where stacks were empty
4. Conclusion and Future Scope
The epidemic spread model using pushdown cellular automata has been considered here in this work to examine the effects of the time variant infective population that remain in a cell for a certain time period. The effective use of the stack to hold these infective population of the previous time periods and effect of these infective population on the infection spread in the coming time steps are studied with two models: one with the homogeneous population across all the cells with no movements between the cells and the other one with different sizes of population in each cell with movement between the cells across the desired locations. The effect of quarantine and the vaccination process is also being investigated with the homogeneous model. The analysis helps us to understand the nature of the epidemic spread across the geographical region and to plan the quarantine and evacuation strategy in case of a pandemic. This model has been tested with the H1N1 data for Abu Dhabi over the period of twelve months from May 2009 to April 2010. Future work can be carried out in this area by considering 3D cellular automata which can take care of the population density like the people living in buildings on multiple floors. The other variation of considering two different types of influenza infection which has got two different recovery rates also can be explored.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.