The role of spatial heterogeneity in the evolution of local and global infections of viruses

Viruses have two modes spread in a host body, one is to release infectious particles from infected cells (global infection) and the other is to infect directly from an infected cell to an adjacent cell (local infection). Since the mode of spread affects the evolution of life history traits, such as virulence, it is important to reveal what level of global and local infection is selected. Previous studies of the evolution of global and local infection have paid little attention to its dependency on the measures of spatial configuration. Here we show the evolutionarily stable proportion of global and local infection, and how it depends on the distribution of target cells. Using an epidemic model on a regular lattice, we consider the infection dynamics by pair approximation and check the evolutionarily stable strategy. We also conduct the Monte-Carlo simulation to observe evolutionary dynamics. We show that a higher local infection is selected as target cells become clustered. Surprisingly, the selected strategy depends not only on the degree of clustering but also the abundance of target cells per se.


Author summary
Viruses such as human immunodeficiency virus and measles virus can spread through physical contact between infected and susceptible cells (cell-to-cell infection), as well as normal cell-free infection through virions. Some experimental evidences support the possibility that high ability of cell-to-cell infection is selected in the host. Since the mode of spread affects the evolution of life history traits, it is important to reveal what condition favors high ability of cell-to-cell infection. Here we address what level of cell-to-cell infection is selected in different target cell distributions. Analysis of ordinary differential equations that keep track of dynamics for spatial configuration of infected cells and the Monte-Carlo simulations show that higher proportion of local infection is selected as target cells become clustered. The selected strategy depends not only on the degree of clustering but also the abundance of target cells per se. Our results suggest viruses have more chances to evolve the ability of local infection in a host body than previously thought. In particular,

Introduction
Viruses have evolved various mechanisms to spread within a host body and between hosts. There are two modes for viral spread in a host body, one is to release infectious particles (virions) from the infected cells into the extracellular medium, and the other is to infect directly from an infected cell to an adjacent cell. The mode of viral spreading depends on the type of virus, their target cells and tissues. For example, viruses that lyse the host cell rely on the release of virions as the only way of spreading. In contrast, viruses that exit host cells by budding or some forms of exocytosis have a potential to spread directly from cell to cell. Conceptually the simplest mechanism of cell-to-cell spread is the fusion of infected and uninfected cells. To enter into a host cell, some viruses have proteins that cause membrane fusion, and these fusion proteins are expressed on the cell surface after viral replication is initiated. Thus, fusion proteins on the infected cell may cause membrane fusion to the adjacent uninfected cell, resulting in a single giant cell (syncytium). For instance, vaccinia virus forms two different forms specific to each mode: mature virus released after lysis of infected cells, and double membraneenveloped extracellular virus that remains associated with the producer cell surface and spreads by cell-to-cell [1,2]. Influenza virus have the potential to spread in a cell-to-cell manner but inherently release virions [3,4]. Other more sophisticated mechanisms of cell-to-cell virus spread also exist (for more examples, see [5]). Both cell-free and cell-to-cell modes of viral spread have their own advantages and disadvantages [5]. Since virions are much smaller and more resistant to environmental change than infected cells, they can disperse farther from the infected cell and even outside the infected host. However, cell-free infection takes a longer time for the virus to encounter a target cell and to engage attachment and entry receptors because the new infection event depends on diffusion and kinetic processes. This is particularly disadvantageous for viruses that bind to receptors that have low expression on host cells and/or those that must engage multiple receptors in order to enter the cell. Immunological barriers to free virions such as antibodies, complement, defensins and macrophages are also factors that discourage cell-free mode. In contrast, viruses that use cell-to-cell infection can avoid many of such obstacles. Another advantage of cell-to-cell infection is the efficiency of new infections: once an infection has occurred, the cell-to-cell mode of viral spread eliminates the rate-limiting step of diffusion. The disadvantage of cell-to-cell infection is a locality of new infection: as infections progress, "self-shading" occurs whereby host cells near infected cells decrease. Cell-mediated immunity, mainly caused by killer T cells, also has a large impact to cell-to-cell infection because infected cells are clustered.
The present study focuses on the intra-host evolution of cell-free and cell-to-cell infection, which we refer to as global and local infection respectively, with spatial structure of target cells. Since the mode of spread of a pathogen affects the evolution of its life history traits, such as virulence [6][7][8][9][10][11][12], it is important to reveal what level of global and local infection is selected. Previous studies assumed that there is a trade-off between global and local infection [6][7][8][9][10][11][12]. Our present study also uses this assumption that denotes the retention of virions on the infected cell surface and polarization to the side of cell-cell contact promote the efficient local infection but interfere with global infection [13]. Because virions could be transmitted to different susceptible cells during global infection, local infection can "waste" some virus particles by putting them all into one cells. As an example, previous study that addressed the evolution of global versus local spread, assuming spatial host dynamics whereby reproduction of host individuals is always local and, as a result, the host population is spatially structured [14]. They showed that the spatial structure is important for the evolution of local infection but this conclusion is limited to particular spatial structures generated by infection dynamics and host reproduction. Such self-organized structures are only a part of possible spatial structures; spatial distribution of target cells in a body, for example, have far more diversity than those self-organized by a particular epidemiological model. Another example of previous theoretical studies include local infection through virological synapses in retroviruses [15,16]. These studies considered viruses that spread via cell-free and cell-to-cell infection, and viral strategy is defined as the number of viruses passed per virological synapse. Since total number of virions produced before an infected cell dies is limited, putting virions all into one cell can "waste" those virus particles. The strategy selected in evolutionary competition can be an intermediate number of viruses passed per synapse (i.e. evolved viruses make use of local infection) depending on the viral kinetics. However the effect of spatial structure was only implicitly discussed. To date, the relationship between a measure of spatial structure and an evolutionary outcome remains unclear in the literature.
In the ecological context, our focus is interpreted as the evolution of short-vs. long-range dispersal. Harada [17] modeled population dynamics in a lattice-structured habitat and assumed a linear trade-off between global and local dispersal. In this model, the assumption that vacant sites are always available also causes the similar limitation as in Kamo and Boots [14]. Hiebeler [18,19] assumed the mixture of suitable and unsuitable habitats in the lattice space in which individuals reproduce globally and locally. In these works, pairwise invasibility was examined by Monte-Carlo simulation but the evolution of the proportion of global reproduction along an adaptive dynamics framework was not discussed.
The purpose of this paper is to analyze how the evolutionarily stable proportion of global (or local) infection is related to spatial heterogeneity. We model the evolution of the proportion of global and local infection in a spatially structured SIS model in which some sites are occupied by target cells and other sites are occupied by non-target cells on a lattice space. For simplicity, we assume that whether a site is occupied by a target cell or a non-target cell will never change. Although it may seem unrealistic to assume that an infected cell directly return to susceptible state (SIS model), this would be justified if a susceptible cell fills a blank immediately after an infected cell dies. The manner of target and non-target cell distribution over the lattice space was parameterized by the frequency of target cells and the pair frequency of target cells. We also assume viruses have an ability to establish persistent infection like human immunodeficiency viruses (HIV) that have enough time for within-host evolution and thus adaptive dynamics can be applied. The dependence of evolutionary outcomes on parameters that denote spatial structure is the main focus of our research, which has not been clearly shown in previous studies. At first, the infection dynamics are modeled by pair approximation and the endemic condition for a virus with a certain proportion of global infection is calculated. Next, the evolutionarily stable strategy is obtained by an adaptive dynamics framework. Finally, we conduct the Monte-Carlo simulation and compare results with analytical results.

Model
Spatial structure of cells. We consider a two-dimensional lattice space in which each site is occupied by either a target cell or a non-target cell for a particular virus. The state of target cells dynamically change between susceptible (S) and infected (I) whereas non-target cells (O) never change.
Let x σ be the global density of the sites with state σ in the population (σ is either S, I, or O), and p σσ , be the probability that a randomly chosen site has state σ and one of its randomly chosen nearest neighbors has state σ 0 . Throughout this paper, we call p σσ , pair density. The conditional probability that a randomly chosen σ site has a σ 0 at its nearest neighbor is denoted by q σ 0 /σ = p σσ 0 /x σ . By definition, we have a following relationship about p σσ 0 Note that x σ = ∑ σ 0 p σσ 0 . We define two parameters to characterize the spatial configuration: x C , the global density of target cells, and p CC , the pair density of target cells. One of the main interests of this study is the dependence of the evolutionary outcome on x C and p CC . In a pair of two neighboring target cells, three possible states should be considered (SS, SI, and II), where the density of IS state is exactly the same as SI state due to Eq (1). The pair density of two target cells in the population is fixed to p CC , and hence On the other hand, there are two states (SO and IO) in a pair of target and non-target cells.
The sum of these densities is constant and because x C = p CC + p CO , Although there are six possible states in a pair, only three of them are independent because we have three constants, the total fraction of sites (x S + x I + x O = 1), the global density of target cells x C , and the pair density of target cells p CC . The degree of cell clustering is denoted by p CC =x 2 C , where p CC =x 2 C ¼ 1 is equivalent to zero spatial correlation and called complete spatial randomness (CSR). Target cells are aggregated more than CSR when p CC =x 2 C > 1, whereas cells are more uniform than CSR when p CC =x 2 C < 1 (See Fig 1A-1C). Infection dynamics. The spatially structured infection dynamics is modeled by pair approximation and the main scheme is the same as Boots and Sasaki [8]. Let global and local infection rate be β G and β L respectively, and we define the proportion of global infection G (0 G 1) as a trait that viruses can change through mutation. Then, for a susceptible cell, the probability of getting infected by global infection per unit time is β G Gx I . In case of local infection, there are two possibilities: infection between cells that we are focusing on or infection from an infected cell standing at the other 3 neighboring sites. The probabilities of getting infected by local infection per unit time are β L (1−G)θ for the former case, and β L (1−G)(1−θ)q I/S for the latter case, where θ is the inverse of the number of nearest neighbors for each site (θ = 1/4 for our twodimensional lattice situation). For convenience, we use following expressions, where g = β G G, and l = β L (1−G)(1−θ). Using Eq (3), the probability that a susceptible cell in a SI pair gets infected per unit time interval is expressed as ϕ + ψ. Infected cells either die and are immediately filled by new susceptible cells, or recover to susceptible state with rate α. To express infection dynamics by pair approximation, we choose x I ( p SI + p II + p IO ), p SI and p IO as three independent variables, For more explanation, see equations (A1) and (A2) in S1 Text.
To consider the evolution of virus trait value G, we use invasion analysis which examines whether a mutant strain can increase or not when it arises at a resident equilibrium. We can obtain analytical expressions of the resident's equilibrium densities for a given value of G from implicit relationship (A6)-(A8) in Appendix A. A resident equilibrium exists as long as the parameters satisfy the endemic condition that is obtained in the next section. At the resident equilibrium, mutants (denoted by subscript J) with proportion of global infection is G 0 , are introduced with a small initial frequency. In the full dynamics of residents (with G) and mutants (with G 0 ), there are ten possible pair densities but we only have to consider seven pairs due to the three constraints, Eq (2) and x S + x I + x O = 1. The dynamics are obtained by adding four more equations to the Eq (4): where g 0 , l 0 , and ψ 0 are as defined before but the definition of G is replaced by G 0 , and ϕ 0 = g 0 x J + l 0 p SJ /x S ,. From the dynamics of resident and mutant, the linearized dynamics of a mutant strain around the resident equilibrium are obtained. If the mutant can increase its density, which is judged by inequality (B9) in S1 Text, we infer that the mutant can invade the resident population. To find an ESS, we assume 101 different possible strains with equally divided G values from G min = 0 to G max = 1, and pairwise invasibilities are checked between two neighboring strains. If a strain prevents the invasion by adjacent strains, we regard it as a local ESS.

Monte-carlo simulation.
To confirm these analytical predictions, we also conducted Monte-Carlo simulations in a 100 ×100 lattice space with periodic boundary conditions. At first, two parameters x C and p CC are determined. Next, according to the parameters, Metropolis-Hastings algorithm generates the spatial structure as follows: 1) each site is made to be occupied by a target cell with probability x C , 2) exchange the states of two randomly chosen sites if current p CC approaches to the goal by exchanging, 3) repeat step 2 until current p CC becomes sufficiently close to the goal. For infectious dynamics of a single strain, each target cell changes its state in a short time interval with some probability that is expected by ordinary differential equations Eq (4). In the simulation of evolutionary dynamics, we also assume 101 different possible strains and each strain can change by mutation to a strain of adjacent G value at a certain rate. Other infectious dynamics, occurrence of new infection and transition from infected to susceptible, is similar to the case of a single strain.

Condition for endemic equilibrium
For the first step, we checked whether a strain with a certain G value can be endemic or not by using the stability analysis of the disease free equilibrium. A similar analysis is done by Hiebeler [18] but the endemic condition was calculated for extreme cases (G = 0 or G = 1) in that study. Here we showed that the endemic condition is also obtained for a virus strain with intermediate G value. In addition to the stability analysis, we also obtained the next generation matrix [20] from the linearized dynamics and calculated basic reproductive number (R 0 ), as, The derivation is shown in Appendix A (S1 Text). When G = 1, R 0 is β G x C /α ρ 1 , which is consistent with the result from the SIS model without spatial structure. In this case, the infection becomes endemic when ρ 1 > 1 that means an infected cell infects more than one susceptible cell. When G = 0, R 0 is β L (1−θ)q C/C / α ρ 0 . ρ 0 > 1 is consistent with the "dyad heuristic" of Levin and Durrett [21], that is, a pair of infected cells will reproduce more than one pair of infected cells. For general values of G, the endemic condition is obtained by stability analysis of the disease free equilibrium, By solving (7) with respect to α, the result becomes consistent with R 0 > 1. For the intermediate value of G, the condition (7) is rewritten by using ρ 1 and ρ 0 , Fig 1 shows the region of G that satisfies endemic condition (7) with changing the recovery rate α. The difference among Fig 1D-1F is the degree of spatial correlation, p CC =x 2 C (but it is not exactly same as the spatial correlation ðp CC À x 2 C Þ=ðx C À x 2 C Þ). Of course, viruses cannot be endemic when the recovery rate α is too high. With increasing α, highly locally infecting strains drop out first when host cells distribute like CSR or when cells distribute more uniform than CSR. On the other hand, strains with intermediate G value are more resistant to the increase of α than other strains when host cells distribute with a positive correlation. Results of the simulation (Fig 1G-1I) indicate that the probability of survival in several trials show similar dependence as predicted by condition (7).

Invasibility analysis and evolutionarily stable strategy
Using the invasibility analysis, we drew a pairwise invasibility plot (PIP) ; Fig 2A is an example of a single parameter set in which intermediate value of G is evolutionarily stable strategy (ESS). In all parameter region examined, there is a unique ESS and ESS strategy is any of completely global (G = 1), a mixture of local and global infections (an intermediate G), or completely local (G = 0). In the present model, there are no other patterns like the evolutionary branching, or more than two evolutionary singular points in the present model. The dependence of ESS on the parameters is shown in Fig 2B-2D. When target cells are relatively less clustered like CSR, ESS G is independent of α (red line in Fig 2B). In contrast, when target cells are relatively clustered, increase of α makes the ESS proportion of global infection higher (green and blue lines in Fig 2B). This is because high recovery rate increases the density of disease-free clusters, which makes global infection beneficial in accessing the isolated clusters. This dependence is also observed with β G < β L (S1 Fig). In Fig 2C and 2D, the dependence of ESS on the degree of clustering for target cells, p CC /x C 2 , is shown. In general, the higher the degree of target cell clustering becomes, the more local infection is optimal. When the rates of global and local infection are equal (β G = β L ), the threshold below which the completely global infection becomes ESS is CSR, p CC =x 2 C ¼ 1 (Fig 2C). In addition to the dependence on p CC /x C 2 , the ESS proportion of global infection also depends on the fraction of cells x C alone ( Fig 2C). If x C becomes higher with fixing p CC =x 2 C , the ESS level of global infection becomes lower. The reason is that in this alteration, the conditional probability that a randomly chosen target cell has a target cell at its nearest neighbor, q C/C = p CC /x C , becomes higher. Thus, local infection becomes more efficient in finding susceptible cells than global infection. The threshold point at which the completely global strain cannot be ESS does not change by altering x C , which is analytically shown in the next section. When the two infection rates differ, the ESS G value tends to prefer the infection mode of better efficiency (Fig 2D). It should be noted that the ESS proportion of global infection is not always a R 0 maximizing strategy (Fig 2E). When p CC =x 2 C is very high, the ESS proportion is much lower than the G value that maximizes R 0 . In most epidemiological or infection dynamics without structure, the ESS trait is to maximize R 0 [22,23]. When cells are spatially clustered, however, the ESS is not always maximizing R 0 , which may be due to a "self-shading" problem as discussed in a previous study [8].

Invasivility to the completely global strain
Here we consider the special case in which a mutant strain with a certain G 0 (< 1) invades a completely global resident strain (G = 1). The endemic equilibrium of the completely global strain is obtained from Eq (4) (for the calculation, see equations (A9) in S1 Text), wherex I ;p SI andp IO are pair densities at the equilibrium. In this case, we can analytically obtain the condition for a mutant strain to increase its density around the resident's endemic The right-hand side of (8) denotes the degree of spatial correlation, and this is the reason why we choose p CC =x 2 C as the horizontal axis of Fig 2C and 2D. The left-hand side of (8) represents an increasing function of the mutant's proportion of global infection G'. Therefore if b G =b L > p CC =x 2 C holds, the completely global strain prevents any kind of mutant from invading and it is an ESS. Especially when β G = β L , the threshold of the spatial configuration at which the completely global strain can be the ESS corresponds to complete spatial randomness (CSR). It means that when cells distribute with a negative correlation, the completely global strain is the ESS, but when cells distribute with a positive correlation, some degree of contact infection can be beneficial. As predicted in this section, the threshold point moves to β G /β L when β G 6 ¼ β L (Fig 2D).

Monte-Carlo simulation of evolutionary dynamics
The mean G in the population quickly converges to a certain level, and fluctuates around that level ( Fig 3A). As long as we use the same parameters, the evolutionary outcomes are similar to each other regardless of the initial condition and the number of trials. We found that evolutionary branching never occurs and the distribution of strains is always unimodal (see S2 Fig). Fig 3B and  3C shows evolutionary outcomes with changing the degree of spatial correlation p CC =x 2 C like Fig  2C and 2D, where 20 trials are conducted for each parameter set. In general, the results are similar to the pair approximation, notably, 1) high p CC =x 2 C prefers local infection, 2) evolutionary outcome depends on both p CC =x 2 C and x C , and 3) increasing x C promotes the local infection. The difference between the results from pair approximation (Fig 2C and 2D) and those from simulations (Fig 3B and 3C) is that the mean value of G in the simulations does not converge to extreme values (the completely global or the completely local). This may be because the population is always polymorphic as a result of mutations. In addition, there are two other reasons for the region in which the completely local is the ESS according to the pair approximation. The first reason is the effect of finite population size. Especially, since the completely local strain spreads only in a contiguous cluster, the number of available hosts is smaller than for other strains. Therefore, diffusing to other clusters becomes adaptive. The second one is the limitation of pair approximation. In this approximation, we approximate q σ/σ 0 σ@ the conditional probability that a randomly chosen nearest of a σ'σ'' pair has a σ site by q σ/σ 0 the conditional probability that a randomly chosen nearest neighbor of σ' site is a σ site. When viruses are too biased toward local infection, infected cells tend to form a large cluster and the approximation does not work well. For these reasons, there is a discrepancy between pair approximation and simulation.
To check whether these results are specific to our method of generating a spatial structure, we also conducted the same evolutionary simulations on the several different deterministically generated structures each having the same global and pair densities (x C and p CC ) as those in randomly generated structure. the evolution of local and global infections based on the randomly generated spatial configurations for given singlet and doublet densities, x C and p CC .

Nonlinear trade-off
In the above sections, we assumed a linear trade-off, that is, the proportion of local infection decreases at the same amount as the proportion of global infection G is increased. However, this should not be always the case in reality. When the local infection decreases in proportion  to G 0.5 , the result is quite different (Fig 5). In this case, the PIP shows bistability in which the evolutionary outcome depends on an initial state (Fig 5B). This prediction by pair approximation is also confirmed by simulation (Fig 5C).

Discussion
In this paper, we investigated the effect of population structure on an evolutionarily stable proportion of global infection. Before considering the effect, we have obtained the endemic Role of spatial heterogeneity in the evolution of local and global infections condition of a virus strain in the SIS model on a lattice space using pair approximation. For the extreme cases (completely global or completely local strain), the endemic condition is consistent with previous studies [18,21]. For the intermediate global infection rate, the condition shows that the strain using both global and local infection can survive with the parameter sets with which the extreme strains cannot. For example, when the global density of target cells is too low to persist for completely global infection, high degree of clustering of target cells promotes the survival of a strain using local infection.
In fact, the invasibility analysis by pair approximation and the Monte-Carlo simulation show that evolutionarily stable strategies or evolutionary outcomes are strongly dependent on the global density of target cells and their degree of clustering. As target cells become aggregated, higher proportion of local infection is selected. This tendency can be seen in the negative slopes in the relationship between ESS G and the spatial aggregation measure p CC =x 2 C in Fig 2. This effect seems stronger at high global density of target cells because the slope of high x C is steeper than that of low x C (Fig 2C). When the local density of target cells is high, spread through local infection is easier for finding the next susceptible cell compared to global infection. Thus, a virus that uses local infection may predominantly be efficient in spreading and that trait is selected. However, there is also a difference between invasibilty analysis and simulation, whereby the former predicts the completely local strain can be ESS in some range but the latter does not. This discrepancy may be attributed to the limitation of pair approximation. In the simulation of evolutionary dynamics, we can generate different spatial structures that have the same first and second moment (global and pair densities) but have different higher order moments. In spite of these variations in spatial structure, the evolutionary outcomes were very similar (Figs 2, 3 and 4). Therefore, we conclude that global and pair densities are enough to explain outcome of evolution of global or local infection, and that parameters which describe higher order configuration do not affect results.
We also checked the dependence of ESS on the parameters that define infection dynamics. When the ratio β G /β L is changed with fixing p CC =x 2 C , the ESS is adjusted to use more efficient pathway (for example, an ESS G for β G /β L > 1 is larger than that for β G = β L ; Fig 2D). The analytical result show that the magnitude of β G /β L relative to p 2 CC =x C also affects whether or not the completely global strain can be an ESS. When the degree of clustering p CC =x 2 C is smaller than β G /β L , the completely global infection becomes ESS. If β G = β L , the boundary is at complete spatial randomness (CSR). If β G < β L , using some local infection can be evolutionarily stable even if targets cells are more uniformly distributed than CSR. The effect of increasing recover rate α is against to favor local infection. This result is apparently counter-intuitive because high recovery rate would mitigate the self-shading effect [14] and hence local infection would be more efficient. However, high recovery rate increases the density of disease-free clusters, which makes global infection far beneficial in accessing the isolated clusters.
The evolution of global vs local spread is also addressed in some studies in ecology and epidemiology [14][15][16][17]. In the cited studies, spatial structure is either not considered at all or limited only to some patchy distribution. Self-organized spatial structures are not definitely controlled or parameterized during simulation. Thus, our model defines spatial structure at first, and examines the dependences of evolutionary outcomes on spatial parameters. The other difference is the effect of "self-shading". If infected individuals are clustered, the number of susceptible individuals available for local infection decreases. Kamo and Boots [14] assumed that infected individuals will die and vacant sites will be covered only by host local reproduction, which has a strong effect against local infection. In contrast, our model assumes SIS model in which an infected cell changes its state to susceptible and the spatial structure is mainly caused by definition. Therefore, self-shading effect is weaker than the previous studies and our analytical results show the suitable parameter region where the completely local strain can be an ESS. In terms of dynamics on artificially organized spatial structure, Hiebeler [18,19] modeled the competition between species using different proportion of global and local reproduction. In such models, the pairwise invasibility between resident and invader (or mutant) was checked, but the evolution of the trait after invasion was not considered. Here we used similar model for the within-host viral evolution and applied adaptive dynamics framework. It suggests that invasion has the possibility of replacement by mutants and this replacement drives the evolution of traits in a population. Applying adaptive dynamics is justified by the assumption that we focus on viruses that cause persistent infection. Such persistent viruses may have sufficient time for the within-host evolution and the adaptive dynamics can be applied. There should be conflicts between the aims for the short-term increase within a host and that for the long-term spreads between hosts. However, the combined effects of these conflicting selection processes is beyond the scope of the present paper.
As in the previous models [14][15][16][17], we assumed a trade-off between global and local infection. The cost of local infection in our model is superinfection caused by the retention of virions on the infected cell surface, which is also assumed in [15,16]. The importance of superinfection is inferred by the fact that various viruses have a mechanism to avoid superinfection promoted by local infection [24][25][26][27][28][29].

Ecological implication
From an ecological point of view, our model can be applied for considering the evolution of long and short dispersal. Target cells correspond to habitats for animals or plants, with S sites and I sites corresponding to unoccupied and occupied sites, and non-target cells representing unsuitable sites to settle. According to the endemic condition (6), the persistence of a species with some local colonization rate depends on the spatial structure. Our results, if applied to a conservation biological setting, suggest that, even if we conserve the abundance of habitats for an endangered species, extinction might occur only due to the change in spatial arrangement of habitats. In terms of the evolution of dispersal distance, previous studies add other settings such as a trade-off between survivability and dispersal range [30], the population dynamics in local patches [31], kin selection [32], the existence of pests [33], or the disturbance structure [34,35]. However, the effect of spatial structure per se has not been clarified yet. We therefore checked how spatial structure affects the evolution of shot and long dispersal. Our results suggest that the ESS proportion of short dispersal depends not only on the degree of clustering but also on the density of habitats per se. Short dispersal is selected when habitats are clustered, and this tendency is strong when the abundance of habitats is high.
Our model is also regarded as a model of dispersal rate evolution that considers how much an offspring should disperse beyond its natal patch if we arrange habitats like patches. In this interpretation, local infection denotes staying in a natal patch and global infection represents outgoing from a natal patch. The infection rates β G and β L correspond to survival rate of outgoing and staying individuals, respectively. The conventional result in this situation suggests that even if the survival rate of an outgoing individual (β G ) is much lower than that of a staying individual (β L ), at least half of all offsprings should go out [36]. Pair approximation in the present study does not predict this result but this is clearly due to the limit of pair approximation. Actually, Monte-Carlo simulation shows that the completely local strategy is never selected. In addition to the previous result, we suggest that even if the survival rate of outgoing individual is higher than that of staying one (β G > β L ), there is a parameter region in which some offspring should be left due to high p CC =x 2 C . Although a similar result has been suggested [37], the reasoning differs to our study. The previous result is due to the variation of patch quality, which means that an individual born in a good patch should leave its offspring in a natal patch. In our model, high clustering of habitats prefers leaving offsprings in a natal patch without assuming differences in patch quality. When the habitats are highly clustered (p CC =x 2 C > 1), the probability of finding new habitat is higher for local dispersal (q C/C = p CC /x C ) than for global dispersal (x C ). Therefore, the advantage of finding a new habitat can outweigh the disadvantage of low survival rate.
Our model predicts the evolutionarily stable dispersal strategy only but some previous models suggest evolutional branching of dispersal rate [38,39] or evolutionary bistability in which the evolutionary outcome differs in initial state [40]. In general, branching or bistability may occur when the fitness of a phenotype depends on the frequencies of other existing phenotypes and possible phenotypes have a proper trade-off [41,42]. In our model, there is a trade-off between global and local infection. The reason why we only observe ESS is that the linear trade-off is not suitable for causing evolutionary branching or bistability. In fact, if we assume nonlinear trade-off between global and local infection is assumed, we can observe the evolutionary bistability (Fig 5).

Virological implication
Our model suggests that spatial structure has an important role, while it is not commonly considered in the field of virology. In in vitro experimental cases, two-dimensional cell culture is commonly used and this condition is similar to our model. Hence, there is a possibility that the evolution of global and local infection can occur. In an example of culture of measles viruses (MVs), a higher level of local infection was selected for in continuing passages [43], indicating the emergence of mutant viruses with a high ability to induce membrane fusion in vitro. When the evolution of global or local infection occur, other viral traits like virulence will evolve as shown in Boots and Sasaki [8]. They analytically showed that a lower virulence is predicted as infection becomes more local. The importance of local infection in the evolution of influenza viruses is shown experimentally [3,4]: cell-to-cell transmission promotes a faster expansion of the diversity of virus quasispecies and may facilitate viral evolution and adaptation when influenza viruses' neuraminidases are inhibited and virus release from infected cells is suppressed [3,4]. Therefore, we emphasize the relationship between the spatial distribution of target cells and the evolution of viral infection mode when studying in vitro infectious dynamics. If we manipulate cell density and the efficiency of viral transmission by antibodies, viruses that have favorable level of global and local infection may be obtained.
In a host body, it is rare that cells distribute like two-dimensional cell culture systems except for epithelia. Epithelial cells form a continuous sheet and they may be different in susceptibility because of surface molecule expression, response to interferons and other immune cell activities. Therefore, epithelia can be a place where evolution of cell-to-cell infection can occur, and in fact, some viruses are known to have an ability for cell-to-cell infection in epithelia like MVs [44] (for other examples, see Table 1 in [5]). In contrast, influenza viruses can also infect epithelial cells but the evolution of local infection is not known. The reason may be the short length of infection period; influenza infection period can be as short as a week but some viruses like HIV survive in the host body for a long time. Such persistent viruses may have sufficient time for the within-host evolution and the adaptive dynamics framework can be applied. Since MVs also have an ability to establish persistent infection, these viruses may also have a chance to evolve efficient cell-to-cell infection in the host body.
It has been shown that cell-to-cell viral transmission through virological synapse occurs in retroviruses such as human T-lymphotropic virus type 1 (HTLV-1) [45] and HIV [46]. This process is thought to have important role because the contribution of cell-to-cell infection on HIV spread is estimated to be equal to or more that of cell-free infection by comparing static and shaking culture conditions [47,48]. Since infected cells can move in the lymphoid tissue and find a connection to susceptible cells, the spatial viscosity of the infected target cells in the lymph nodes should be weakened in these viruses. However, there would remain some non-random correlation of uninfected target cells (p CC =x 2 C > 1) because of the locality of T cells in the lymph nodes, and this could favor cell-to-cell transmission over cell-free transmission in retroviruses too. Therefore, we suggest new conditions that evolutionarily promote cell-to-cell infection in those viruses: highly localized distribution of target cells in the lymphoid tissues (p CC =x 2 C > 1). These points have not been suggested in the previous theoretical study [16].
Our model may explain the emergence of mutant MVs that are isolated from patients of subacute sclerosing panencephalitis (SSPE). These viruses have mutations that provide high ability of cell fusion (i.e. high level of local infection) [49] and can infect central nervous system cells, while wild type MVs cannot [43,49]. Since MVs can spread in a cell-to-cell manner between epithelial cells, epithelia is a candidate place in which the evolution of local infection occurs. Another possibility is lymph nodes because the main target cell of MVs is SLAM (signaling lymphocytic activation molecule, also known as CD150) expressing immune cells such as T and B cells etc. [50,51]. As discussed in the case of HIV, concentrating target cells in the lymph nodes satisfies the condition under which local infection is selected. Consequently, MVs are prone to evolve local infection in a host body and to gain the ability to infect cells of central nervous system but how MVs reach the central nervous system remains unknown.
In conclusion, our results suggest that the mode of viral spread, global or local infection, may undergo adaptive evolutionary change in vitro and in vivo. In the future, we can consider more realistic situation such as evolutionary dynamics in three-dimensional space or the repulsion of superinfecting virions that attenuates self-shading effect [29]. In order to examine the emergence of mutant virus or the evolution of virulence, we need to take into account the fact that the mode of infection itself is subject to selection.