The Effect of Patchy Host Distribution on the Dynamics and Persistence of Directly-Transmitted Pathogens : A Cellular Automata Study

The effect of fragmented host distributions on the transmission dynamics of directly-transmitted pathogens was explored via stochastic automata simulation. Sixteen diverse population distributions varying in shape and density were used as a substrate for simulated outbreaks. Extended neighborhoods (80 cells), with probability of infection weighted by proximity to an infective source were used to define the overall probability of transitions from susceptible to infected. A static probability defined transitions from infected to recovered. The duration of active transmission as well as the proportion of each population infected per outbreak was averaged over a series of 30 simulations per parameter set. The level of aggregation for each population, measured in terms of the Moran Coefficient (MC) of spatial autocorrelation, was found to affect both the intensity of an outbreak and its length of persistence. Denser populations produced the most cases and lasted longer than those that were sparser. Elongated distributions, measured as the ratio between perimeter and area (PA) reversed some of the trends of increasing density. Long, narrow distributions produced fewer cases and were less persistent than populations composed of more compact clusters but with similar MC. Thus, both the shape and density of host distribution patterns affected the incidence rate, duration of epidemics and the percent of the population infected. Certain patterns of habitat fragmentation, thus, may put more hosts at risk of becoming infected than others.


Introduction
Fragmented populations of endangered flora and fauna may be exposed to greater stresses than populations that are larger and more continuous [1].Increasing the amount of edge or transition area between habitats may increase exposure to predation, parasitism and environmental stressors.From the perspective of certain plant hosts, some of these effects may be positive, such as decreased herbivory when patches become isolated due to the destabilizing effects of habitat fragmentation on local populations of insect herbivores [2,3].Associations between host spatial structure and pathogen dynamics have been observed with rusts and smuts of various plants [4][5][6][7], where patchiness appears to prevent the perpetuation of outbreaks.
For pathogens that rely on close contact for transmis-sion, the configuration of populations in space may greatly affect the rate and stability of epidemic outbreaks.Such effects have been repeatedly noted in the natural documentation and simulation of measles outbreaks [8,9].Distributions of organisms that are organized into discontinuous clusters display especially complex and unpredictable dynamics.The manner in which such ensembles of populations or "metapopulations" are arranged can affect the persistence of certain populations both free-living or parasitic.In metapopulation models of free-living organisms, dispersal can either stabilize [10,11] or destabilize [12] populations depending on how far, and with what frequency organisms disperse.Parasitoid populations seem to persist less stably in the face of strong dispersal in structured metapopulations than do free-living organisms [13].
The mathematical modeling of infectious disease transmission still generally relies on non-spatial methods based on ordinary differential equations (ODEs) pioneered early in the 20th century [14].These techniques treat time and state as continuous entities and assume that populations are homogeneous and completely mixing.Various extensions of the basic method have been developed to recreate some level of natural complexity by simulating spatial or social discontinuity.Stochastic models also allow for the effects of chance on small populations or subsets of larger ones.
These extensions generally rely on subdividing populations into smaller groups with shared characteristics.Such methods become very quickly limited in how much spatial structure they are able to represent before becoming computationally intractable.
Cellular-automata offer an opportunity to study the effect of host spatial distribution on the dynamics of infectious disease transmission in a manner that captures the full complexity of the spatial structure of natural population distributions [15].
Cellular automata are discrete mathematical systems that allow "bottom-up" simulations of complex systems by applying simple sets of local rules to a matrix of integers over time [16].These rules take into account both the state of individual cells and those that neighbor them [17].Rule sets can be probabilistic or deterministic.Neighborhoods may include only adjacent cells, or any number of concentric groupings or layers outward and upward.Thus, automata are ideally suited for simulation of complex systems that depend on local interactions between discrete entities, and useful for understanding dynamics of biological systems [18].
Automata have been widely applied to ecological systems including forest gap dynamics, forest fires [19], and to explore fundamental properties of spatial distributions in eco-systems [20].It the life sciences it has been used, for example, to study dynamics of neural populations [21] and to model proliferation of tissue growth [22].Epidemics modeling with cellular automata has been previously implemented to study the spread of disease among moving individuals [23,24].Applications include simple models of measles outbreaks [25], effect of vaccination [26], to explore drug therapy for HIV [27] and to study transmission of vector borne diseases [28].
Here we apply methods derived from stochastic automata to explore the effect of host organization in specific spatial pattern on the persistence and intensity of outbreaks of directly-transmitted pathogens.

Host Population Distribution
Patchy population distributions were created by first seeding 100 cells randomly in a blank matrix then applying a cellular automata-based growth algorithm that used various rules governing the probability of an empty cell becoming occupied based on the number of neighboring cells that were already occupied.In some cases, probabilities were weighted more heavily along certain axes in order to create long, narrow aggregations.These population growth models were suspended once the population reached 10,000.
To produce more variation in patterns, differing numbers of initial seeds were used to initiate growth and were allowed to expand for different durations until the population size limit was reached.To produce more linear patterns, probabilities of growth were weighted to favor expansion along particular axes.The intensity of the weighting affected the length and width of the resulting pattern which is shown in Figure 1.For each distribution Moran Coefficients (MC) and perimeter/area ratios (PA) were calculated to characterize aggregation as presented in Table 1.The Moran Coefficient (MC) is a powerful method to assess spatial dependence or auto correlations [29,30].MC = 0 indicates total randomness, while MC = 1 indicates maximum spatial clustering; values in between indicate various degrees of aggregation.See Figure 1.

Initial Conditions
All simulations were encoded using the interpreted "J" programming language (Version 4.62, Iverson Software), which is optimized for array manipulation.Simulations were structured as follows: Assumptions and rules: 1) The environment is represented by a square matrix (224 × 224) of 50,176 integers.
2) Population distributions of hosts are represented by 10,000 variously distributed cells.
Patchy population distributions were created by first seeding 100 cells randomly in a blank matrix then applying a cellular automata-based growth algorithm that used various rules governing the probability of an empty cell becoming occupied based on the number of neighboring cells that were already occupied.In some cases, probabilities were weighted more heavily along certain axes in order to create long, narrow aggregations.These population growth models were suspended once the population reached 10,000.
3) Codes representing the current infective state of individual cells are as follows: 0 = unoccupied 1 = susceptible and uninfected 2 = infected and infectious 3 = immune 4) Time (t) advances in discrete units of one week per time step.5) During each time step (t), an infection can be transferred from a source (i) to adjacent or nearby contacts up to 4 cells away from the source, as described in Figure 2.
The basic structure of an automaton is as follows: The future composition of a matrix is thus a function of the current state of that matrix   t X and external factors   t Y such as treatment or vaccination that might be imposed.In the present simulation, only the internal state of the matrix is considered.
6) The probability of a set of contacts giving rise to an infection is calculated using a state transition probability function as in Figure 1.where: inf P = the probability that an uninfected host becomes infected.
f = an index of the relative transmissibility of a particular pathogen., , , to a set of randomly generated numbers between 0 and 1 using the Mersenne Twister algorithm.
8) Births and mortality from other causes are not considered due to the short term (up to several months) that was assumed for the duration of an outbreak.
Each simulation was terminated until extinction of the pathogen occurred and all infected individuals either became immune or dead, however long that required.For each of the eight distributions, simulations were replicated 10 times and averaged.

Results
The spatial characteristics (shape, density) of host population distributions were found to greatly affect the dynamics of infectious disease outbreaks.Epidemic curves generated for each distribution (Figures 2(a)-(p)), were diverse in amplitude and period (Figure 3).The least aggregated population (Figure 2(a) produced a peak incidence of less than 100 daily cases during the simulated outbreak, while the most aggregated conditions (Figures 2(o)-(p)) produced a peak incidence of 5000 daily cases.
Each of the other distribution patterns produced incidence peaks that were intermediate between these two extremes.These incidence peaks reflect the cumulative morbidity for each of the simulated outbreaks.
The mean number of infections occurring in an outbreak tended to be positively associated with the intensity of spatial autocorrelation (Figure 4) as represented by the first order Moran Coefficient (R2 = 56.7%).The average magnitude of an outbreak was also strongly, but negatively associated (R2 = 41.3%) with the perimeter/ area ratio (Figure 5).
Variation in the total number of infections occurring in an outbreak between simulations tended to be negatively associated (R2 = 61.4%)with the degree of spatial autocorrelation (Figure 8).Perimeter/area ratio (Figure 9) showed a milder, positive association (R2 = 40.2%)with the average magnitude of outbreaks.The duration of outbreaks was less strongly associated with either the Moran Coefficient (R2 = 30.5%,)or perimeter/area ratio (R2 = 27.9%).
Spatial statistics that account for degree of aggregation and exposure did not always suffice to explain the variation in outcomes between simulated outbreaks.Other aspects of shape, such as degree of elongation and fenestration may have impacted outbreak dynamics.For ex- ample, cumulative number of infections was lower in elongated patches (MC = 0.84) with a higher perimeterarea ratio (Figure 2(m)) relative to a more clustered distribution (MC = 0.858) with more rounded patches and a lower proportion of individuals along boundaries (Figure 2(i)).The persistence of outbreaks of infectious disease was also affected by the shape of host distribution patterns.Generally, the time required for epidemic outbreaks to resolve and for transmission to cease increases as the degree of aggregation increases.Spatial discontinuities between individuals and groups of indi-viduals increase the probability that transmission will be interrupted and that epidemics will abort.This general trend of extended transmission with increasing density can again be modified by altering the shape of particular distributions, particularly when the force of infection is relatively high.A population of elongated patches (PA = 0.581) persists longer than a population of similar density but lower perimeter-area ratio (PA = 0.449).Thus, both population density and the shape of population distributions can determine the intensity and persistence of outbreaks of infectious diseases.

Discussions
The degree of fragmentation and the connectivity of population distributions contribute to the dynamics of the spread of infectious pathogens.The shape of population distributions can affect both the duration of persistence of an epidemic outbreak.While the intensity of contacts reflected by density was also a critical factor, differences in shape had the ability to negate trends that were otherwise apparent across a gradient of density values.
Shape modulates the intensity of infectious contacts in a population by altering the number of individuals exposed on the edges of a distribution.When the perimeter to area ratio is small, the population behaves more like a homogeneous distribution.For a given density, when the perimeter to area ratio is large, simulations tended to terminate more quickly, and left a greater number of individuals uninfected.
Habitat fragmentation may either protect populations or render them more vulnerable to disease outbreaks, depending on the nature of the pathogen and the arrangement and size of habitat patches.Protection may be provided in certain cases by fragmentation by interrupting transmission between patches.In natural populations of a Swiss wetland plant (Primula farinosa), outbreaks of smut fungus (Urocystis pimulicola) were more likely to occur in larger, more contiguous habitats [7].
Conditions might be considered as well in which fragmentation allows a pathogen to persist by imposing asynchrony on its infection dynamics.While an active outbreak is occurring in one patch, other patches are recovering and increasing their immunologic vulnerability through the recruitment of young individuals.Certain pathogens, especially those which either kill or render hosts permanently immune, might be expected to cause their own extinction in host populations that are distributed uniformly and contiguously.
These results may be important in terms of evaluating trends in local and global biodiversity.Models similar to this may define more precisely the vulnerability of spatially discontinuous populations to pathogens, particularly populations of threatened flora and fauna that have become isolated through habitat fragmentation.Population viability analyses (PVAs) have not been reliable in defining these risks, partly because of the complexity that homogeneous models fail to capture but which can be highly relevant to outcomes.The manner in which hosts are distributed can fundamentally affect the dynamics of disease spread through a population.

Figure 1 .
Figure 1.Host population distribution patterns used as the basis for simulated epidemics.The first level Moran Coefficient of spatial autocorrelation is indicated below each matrix.

Figure 2 .
Figure 2. Diagram of stochastic automaton model structure.
d = the number of infective neighbors at each distance a-d. a d W  = the relative effect of distance on relative transmissibility portraying a geometric decline.7) For each time step, transitions from an uninfected state to infected are determined by comparing the probability of transmission to each uninfected cell   inf P

Figure 7 .Figure 8 .
Figure 7.The relationship of perimeter to area ratio to the total duration of outbreaks for each population distribution pattern (n = 30).

Figure 9 .
Figure 9.The relationship of perimeter to area ratio to the coefficient of variation of the total duration of outbreaks for each population distribution pattern (n = 30).