The effect of aggregation of pathogen and biocontrol microbe propagules on biocontrol potential: a simple modelling study

Effective use of biocontrol agents (BCAs) is a potentially important component of sustainable agriculture. The ecological processes determining the success of biocontrol are complex, which may partly explain the limited success of biocontrol against plant diseases in field crops. Understanding the ecological characteristics of BCAs in addition to biocontrol mechanisms and direct biocontrol efficacy, including their survival and dispersal under heterogeneous conditions, is critically important to improve biocontrol efficacy. In this simulation study, we focused on the effects of the spatial aggregation of initial pathogen and biocontrol propagules (inocula) under spatially homogeneous and heterogeneous conditions on biocontrol potential. The simulation showed that, as expected, increasing the biocontrol propagule density led to increased biocontrol potential. Under a given inoculum density, increasing spatial aggregation of BCAs is not only likely to reduce biocontrol potential but also to increase variabilities in biocontrol outcomes. A spatially random distribution of BCA propagules is most likely to result in best biocontrol outcomes. Spatial aggregation of BCAs is more important than spatial heterogeneity in influencing biocontrol potential. Thus, the present simulation study illustrates the importance of ensuring a close-tohomogeneous distribution of BCA propagules for maximising biocontrol potential. By the same reasoning, a faster host growth rate will reduce biocontrol potential if BCA cannot keep up with host growth in terms of coverage: increasing BCA aggregation essentially leads to decreased effective coverage.


Background
The biological control of plant pathogens involves the use of non-pathogenic (or weakly pathogenic) microbes to control pathogens. Several mechanisms are possible for biocontrol, including direct parasitism, competition for resources (nutrients and host sites), antibiosis and induced resistance (Whipps 1992;Elad 2003). In biocontrol of plant diseases, the focus is usually on the augmented introduction of antagonists to control diseases.
Despite extensive research and development in biocontrol of plant diseases, success in biocontrol of plant diseases in field crops has been limited (Barratt et al. 2018;Syed Ab Rahman et al. 2018), most successes being achieved under more controlled conditions, e.g. greenhouse crops and produce in post-harvest stores. The complex ecological processes involved in biocontrol of crop diseases have been often cited to explain the lack of biocontrol success as well as the variable biocontrol efficacy achieved. Despite the recognition of the importance of ecological knowledge in biocontrol agents (BCAs) for predicting and optimising biocontrol (Xu et al. 2010;Juroszek and von Tiedemann 2011;Jeger 2013a, 2013b), there is limited information on the fate of biocontrol organisms under natural conditions. For instance, theoretical modelling studies of biocontrol of both foliar and soil borne diseases indicated the importance for BCAs to be able to colonise effectively susceptible healthy host tissues (Jeger et al. 2009;Cunniffe and Gilligan 2011). Yet, this specific ecological feature of BCA is often neglected.
One important ecological trait affecting microbial population dynamics is related to the spatial distributional properties of their propagules and their abilities of being spatially dispersed. The spatial distribution and dispersal of propagules have been shown to affect biocontrol efficacies for specific pathosystems, including Ulocladium atrum against Botrytis cinerea (Kessel et al. 2005), biocontrol of fungal pathogens with mycoviruses (Liu et al. 2000), and mycoparasites against soilborne pathogens (Jeger et al. 2004). Bacteria can form biofilms in response to external environmental conditions, which may affect biocontrol outcomes in terms of microbial survival and the extent of contact between antagonists and pathogens on the host tissue surface (Morris and Monier 2003;Nongkhlaw and Joshi 2014). The survival of Pseudomonas syringae cells in aggregates promoted a highly clustered spatial distribution of bacteria on leaf surfaces (Monier and Lindow 2003). Biofilm formation explains the potential importance of cooperative interactions of epiphytes among both homogeneous and heterogeneous populations, influencing the development of microbial communities.
Biocontrol organisms need to share similar ecological niche requirement as the target pathogens in order to survive and even increase their population size in the target environments to achieve biocontrol. This is particularly true if biocontrol is primarily achieved through competition and direct parasitism; but it can also be important for localised antibiosis or induced resistance as the response intensities likely depend on the concentration (hence diffusion) of active antibiotics or chemical inducers (Costet et al. 1999). Several bacterial antagonists occupied different sites from pathogenic P. syringae strains on bean plants (Wilson et al. 1999), reducing the contact between antagonists and pathogens and hence potential biocontrol. When several pairs of common epiphyte bacterial species (P. syringae, P. fluorescens, and Pantoea agglomerans) were pairwise co-inoculated onto leaf surfaces, they formed mixed aggregates, with the degree of cell segregation in mixed aggregates differing between specific pairs (Monier and Lindow 2005). Only a small fraction of the cells of different species was in contact even on heavily colonised leaves.
In commercial practice, microbial BCAs are usually applied as conventional fungicides through augmented application of formulated BCA products. The initial distribution of introduced BCA propagules can largely influence their distance to pathogen propagules, affecting subsequent biocontrol outcomes. This is supported by a limited number of modelling studies, including a spatially explicit simulation study of Ulocladium atrum against Botrytis cinerea (Kessel et al. 2005), spatial host heterogeneity in affecting biocontrol of fungal pathogens with mycoviruses (Liu et al. 2000), and spatial distribution of mycoparasite propagules in soil on controlling soilborne pathogens (Jeger et al. 2004). Furthermore, environmental heterogeneity, on both spatial and temporal scales, is the rule rather than an exception when foliar pathogens are concerned. Such spatio-temporal heterogeneity can affect pathogens and BCAs differentially, thus affecting biocontrol outcomes.
Leaf surface is a nutrient-limited environment and the microclimate on the leaf surface could change rapidly both spatially and temporally depending on external conditions and canopy/leaf architecture. There have been no studies on the competitive interactions between microbes on leaf surfaces in relation to the spatial heterogeneity and distributional properties of microbial propagules. In this study, we conducted simulation studies to assess the joint effects of initial aggregation of microbial (pathogen and BCA) propagules and spatial heterogeneity on biocontrol potential. Specifically, we investigated the effects of the abundance and spatial aggregation of BCA/pathogen propagules, and the extent of environmental heterogeneity on biocontrol potential. Following a previous study (Jeger et al. 2004), we used the distance between a pathogen propagule and its nearest BCA propagule to represent biocontrol potential. Such a distance is critically important for biocontrol outcomes, particularly for biocontrol mechanisms of competition and direct parasitism and, to a lesser extent, for localised antibiosis or induced resistance.

Results
This study considered the following simulation factors (variables): abundance of BCA/pathogen propagules, the extent of their spatial aggregation, and the extent of environmental heterogeneity (grain [patch] size and patchiness [defined as the ratio of BCA survival to the pathogen survival in a given patch]). Contributions of these simulation factors and their interactions to various measures (size and variability) of D min (the distance of a pathogen propagule to its nearest BCA propagule) were very similar. Thus, subsequent presentation focused on the average of D min ( D min ) over all pathogen propagules for a given simulation run.

Homogeneous areas
Under this scenario, both BCA and pathogen propagules were assumed not to suffer any mortality across the entire simulation area.

Random distribution of pathogen and BCA propagules
BCA abundance accounted for 96.76% of the total variability in D min , with the residuals accounting for 3.21%. D min as well as its absolute variability (variance) decreased as BCA abundance increased from 50 to 400. Average D min was 0.074, 0.052, 0.037 and 0.025 for the BCA abundance of 50, 100, 200 and 400, respectively. D min correlated positively with its standard deviation across all simulation runs, r = 0.97.

Clustered distribution of pathogen and biocontrol propagules
The main effect of BCA abundance and aggregation, and their interaction accounted for 31.47%, 47.44% and 1.18% of the total variability in D min , respectively; the residuals accounted for 19.64% of the total variability. D min decreased with increasing abundance and with decreasing aggregation (i.e. increasing sigma value) ( D min correlated positively with its standard deviation across all simulation runs, r = 0.90. Although the absolute variability in D min decreased with increasing BCA abundance and sigma values ( Fig. 1), the relative variability (measured as CV) was greater for intermediate aggregation values (Fig. 2).

Heterogeneous areas
Under this scenario, BCA and pathogen propagules may experience differential mortality (defined as patchiness) in specific niches.

Random distribution of pathogen and BCA propagules
Patchiness accounted 82.30% of the total variability in D min , with the residuals accounting for 17.68%; the effect of grain size and its interaction with patchiness was virtually non-existentindeed statistically insignificant. D min as well as its absolute variability decreased as patchiness increased (Fig. 3). However, most of the Fig. 1 The density of the mean distance between a single pathogen propagule and its nearest BCA propagule in a given simulation run for each combination of biocontrol abundance (50,100,200,400) and aggregation (0.01, 0.02, 0.04, 0.08, 0.16, 0.32) in a homogeneous simulation area. The dashed line indicates the median value patchiness effects (> 95%) was due to the differences in the number of pathogen/BCA propagules retained in each simulation run; increasing patchiness led to reduced number of pathogen and BCA propagules retained. Average distance was 0.034, 0.029 and 0.027 for the respective patchiness of 1:0.25, 1:0.5, and 1:0.75; the corresponding values of the number of pathogen (BCA) propagules retained were 303.1 (303.7), 336.4 (336.4), and 368.7 (368.2).
D min correlated positively with its standard deviation across all simulation runs, r = 0.94.

Clustered distribution of pathogen and BCA propagules
As for the Poisson (random) process, the effect of grain size was minimum; all terms involving this factor only accounted for a total 0.10% of the variability in D min . The single most important factor was BCA aggregation, accounting for 84.98% of the total variability in D min . The residuals accounted for nearly all the remaining variability (14.26%); including numbers of pathogen and BCA propagules retained did not lead to noticeable changes. Increasing aggregation of BCA propagules (i.e. decreasing sigma values) led to increasing D min (Fig. 4). Average D min was 0.101, 0.087, 0.065, 0.042, 0.032 and 0.030 for the sigma value of 0.01, 0.02, 0.04, 0.08, 0.16 and 0.32, respectively. D min correlated positively with its standard deviation across all simulation runs, r = 0.79.

Discussion
We used a simple simulation approach to illustrate the importance of the spatial distribution properties of BCA propagules on potential biocontrol outcomes. The distance of each pathogen propagule to its nearest BCA propagule was used to represent biocontrol potential. As expected, increasing BCA abundance led to increased biocontrol potential. Under both homogeneous and heterogeneous conditions, spatial distributions of BCA propagules was the most important factor affecting biocontrol potential. The more spatially aggregated BCA propagules were the more variable and lower the biocontrol potential was. Unlike previous studies of biocontrol under spatially heterogeneous conditions (Liu et al. 2000;Jeger et al. 2004;Kessel et al. 2005), the effects of initial aggregation of BCA and/or pathogen propagules on the biocontrol efficacy were not confounded with other dynamic aspects of BCA/pathogen development. Moreover, the present study investigated the joint effects of spatial heterogeneity on leaf surfaces and spatial distributional properties of pathogen/BCA propagules on biocontrol potential.
The distance between a pathogen and its nearest BCA propagule is likely to overestimate biocontrol potential, particularly under clustered situations. In the present study, many pathogen propagules could all share the same BCA propagule as their nearest BCA neighbour, especially when pathogen propagules were extremely aggregated. For instance, propagules within a sporulating lesion or bacterial biofilm are extremely aggregated. Bacterial biofilms can protect bacterial aggregates from adverse conditions and reduce their contact with potential antagonists (Monier and Lindow 2003;Morris and Monier 2003;Nongkhlaw and Joshi 2014). Therefore, under these circumstances, a single BCA propagule may not be able to suppress the highly aggregated pathogens, particularly if the biocontrol mechanism involved is competition or direct parasitism.
We showed that biocontrol potential was primarily determined by BCA aggregation and abundance. However, this does not mean that pathogen distributional properties are not important. The lack of effects of pathogen distributional characteristics on biocontrol potential results simply from the fact that the minimum distance between pathogen and BCA propagules was calculated from the perspective of pathogen propagules (i.e. conditioned on the spatial pathogen distribution). If we had estimated biocontrol potential from the BCA perspective, namely using the distance of each BCA propagule to its nearest pathogen propagule, pathogen characteristics would then have become most important because of Fig. 3 The density of the mean distance between a single pathogen propagule and its nearest BCA propagule in a given simulation run for each combination of heterogeneity grain size and patchiness for randomly (Poisson) distributed pathogens and BCAs. The dashed line indicates the median value the symmetry in the pathogen and BCA parameter values. We estimated biocontrol potential from the pathogen perspective because of the following two reasons. Firstly, for biocontrol of plant diseases in commercial practice a BCA is usually applied to crops via regular augmented applications. Thus, disease locations can be viewed as pre-existent and hence biocontrol efficacy is conditioned on pre-existing disease patterns. Secondly, agricultural practitioners can usually alter BCA properties via formulation and application technologies.
The extent of patchiness, in terms of both the grain size and relative BCA/pathogen survival, did not affect biocontrol potential in the present study. This apparent lack of patchiness effect results most likely from the symmetric simulation parameter values used in the present study. Both the pathogen and BCA have the same levels of aggregation at the same abundance; the overall relative survival across the whole simulation grid is the same for both the pathogen and BCA. Thus effectively even if we had exchanged the 'pathogen' and 'biocontrol' labels it would not have resulted in any differences in the simulation results. However, if the net survival across the whole simulation grid for the BCA had been set to be greater than the pathogen, then biocontrol potential would have increased; vice versa.
The present results suggest that reducing aggregation of BCA propagules whilst maintaining the overall population density is the most important first step to improve biocontrol potential in practice. One way to achieve this is to improve spray coverage. The coverage of wheat spikes by B. amyloliquefaciens post application was inadequate for protecting spikes from Fusarium head blight pathogens in the field (Crane and Bergstrom 2014). Commercial field spraying equipment for strawberry could result in not only large unprotected areas (30-70%) but also large variability in the coverage (Xu, unpublished). Application technologies could thus be Fig. 4 The density of BCA abundance in a given simulation run for each combination of biocontrol aggregation (0.01, 0.02, 0.04, 0.08, 0.16, 0.32) and patchiness (0.25, 0.50, 0.75) in a heterogeneous simulation area with a small grain size. For all simulations, the initial abundance for the BCA was 400, some of which may be discarded when they were in the patches favouring the pathogen. The dashed line indicates the median value further developed to improve coverage to take into account weather conditions and characteristics of target tissues in relation to crop canopy structure. BCA formulation could be improved to facilitate their local dispersal (i.e. adding adjuvants) and survival (i.e. adding nutrients). We need to pay attention to other factors that could affect BCA survival and distribution postapplication, such as rapid host growth that could lead to a spatial distribution of BCA propagules with increased spatial aggregation as well as with a reduced overall density. BCAs should ideally possess the ability of rapidly colonising crop surfaces to increase coverage and hence keep up with host growth.
Improving BCA movement (dispersal) within and among host tissue units are similarly important to ensure better coverage of susceptible tissues by BCA. There was limited dispersal of B. subtilis among leaves on the same strawberry plants under both open and protected environment (Wei et al. 2016). This limited dispersal/spread/colonisation of BCAs within and between host tissue units is another reason that may explain why biocontrol of post-harvest diseases is more likely to succeed when BCAs are applied post-harvest than pre-harvest (Dukare et al. 2019). For biocontrol of post-harvest diseases, host tissues usually do not expand and environmental conditions are generally very stable. Thus, good initial application coverage (often achieved by immersing in biocontrol suspensions) and ability for BCA to survive under the storage conditions will maximise the chance of successful biocontrol.

Conclusions
We may conclude that under a given inoculum density, increasing spatial aggregation of BCA propagules is not only likely to reduce biocontrol potential but also to increase variabilities in biocontrol outcomes. Thus, it is critically important to ensure a close-to-homogeneous distribution of BCA propagules in order to maximise biocontrol potential in commercial agriculture.

Overview of the study methodology
The present study investigated the effect of initial aggregation of plant pathogen and/or BCA propagules within a given area (e.g. leaf surface) on the distance between pathogen and biocontrol propagules, which would greatly affect biocontrol outcomes. The spatial coordinates of individual pathogen or biocontrol microbial propagules were simulated using simple stochastic point processes with the R package "mobsim" (May et al. 2018), either as a Poisson process, or as a Thomas process (Wiegand and Moloney 2014). For the Poisson process, individual propagules were placed randomly; in contrast for the Thomas process, individuals of the same species were spatially clustered (Wiegand and Moloney 2014). For simulating the Thomas process, a common model of intraspecific aggregation in ecology, "mobsim" needs the numbers and sizes of the clusters, as well as the number of individuals per cluster, either independently or jointly for all species. The Thomas process only considers intraspecific aggregation whilst individuals of different species are distributed randomly with respect to each other (McGill 2011).
In the present simulation study, we did not specify the number of clusters or the number of individuals per cluster. Thus, the "mobsim" programme uses the square root value of the total abundance as the number of clusters and the number of individual per cluster. Therefore, this study did not consider the effect of variable cluster sizes. Furthermore, simulation parameter values as well as the number of values for each simulation variable were chosen arbitrarily, particularly regarding the nature of a spatially heterogeneous simulation area. However, the exact values of these simulation variables did not affect main conclusions drawn from this simulation study as the focus was on the relative distance between pathogen and BCA propagules following a BCA application. We did not consider subsequent microbial dispersal that critically depend on the absolute grain size relative to the propagule size and the extent of its dispersal.

Nature of the simulated area
This simulation study considered two scenarios: spatially homogeneous and heterogeneous areas. Under the scenario of a homogeneous area, both pathogen and biocontrol propagules had an equal survival ability (= 1.0, i.e. without mortality) across the entire simulation area. To simulate a heterogeneous area, we considered two specific aspects of spatial heterogeneity: grain (patch) size and patchiness (defined as the ratio of BCA survival to the pathogen survival in a given patch). Two grain sizes were considered: large and small. For the large grain size, the simulated area was divided to 25 equal square sections with the pathogen surviving better than the BCA in eight of the 25 square sections, and vice versa. For the remaining nine sections, both the pathogen and BCA had an equal survival ability as for the homogenous case. For the small grain size, the simulated area was divided into 100 equal square sections. In 32 of these sections, the pathogen survived better than the BCA and vice versa; in the remaining 34 sections both the pathogen and BCA had an equal survival ability as for the homogenous case. The exact location of individual sections for each of the three categories (favouring pathogen, favouring BCA, equal) was randomly assigned. Examples of simulated distributions of pathogen and propagules are given in Fig. 5.

Simulation parameters
For each combination of simulation parameters, there were 100 replicate simulation runs.

Homogeneous area
There were four simulation variables: pathogen aggregation index and abundance, and BCA aggregation index and abundance. Six values of aggregation index (i.e. the sigma value in the "mobsim" programme) were used for both the pathogen and BCA: 0.01, 0.02, 0.04, 0.08, 0.16, and 0.32. The lower the value, the more aggregated the simulated points are. In addition to these 36 combinations of pathogen and BCA aggregation, a single random case was included: both pathogen and biocontrol propagules were randomly located in the simulation area (i.e., simulated by a random Poisson process). Four levels of abundance were used for both the pathogen and BCA: 50, 100, 200, and 400, giving a total 16 pathogenbiocontrol combinations. Thus, for the homogeneous scenario, there were 592 (37 × 16) combinations of simulation parameters. Each combination of these parameters was used by the "mobsim" to simulate distribution of the pathogen and biocontrol propagules.

Heterogeneous area
For each grain size, there were three simulation variables: pathogen and BCA aggregation, and patchiness. As for the homogeneous case, there were 37 combinations of pathogen and BCA aggregation. For the patchiness, three ratios of the BCA survival to the pathogen survival for the grains favouring the BCA (or the pathogen survival to the BCA survival for the grains favouring the pathogen) were used: 1:0.25, 1:0.5, and 1:0.75. Thus, there were a total 222 combinations of simulation parameters across the two grain sizes. A single abundance of 400 was used for both the pathogen and BCA as the patchiness would lead to further reductions in the number BCA and pathogen propagules. For this reason, we did not use abundance values lower than 400 in the present simulation study.
A two-step simulation was used to simulate the distribution of pathogen and BCA propagules in the grid. First, as for the homogeneous case, locations of 400 pathogen and 400 BCA propagules were simulated with appropriate aggregation indices. Then for those BCA propagules located in the square sections favouring the pathogen, a uniform random number (x) in the range of 0-1 was generated to determine whether each BCA propagule can survive (i.e., to be retained). If x was less than the patchiness parameter (0.25, 0.5 or 0.75), it was retained; otherwise it was discarded. For all other BCA propagules located in the sections of either favouring BCA or neutral), they were all retained. The same procedure was applied to the 400 simulated pathogen propagules. Therefore, the number of pathogen or BCA propagules retained was less than 400 and varied with replicate simulations.

Analysis of simulation data
In the present study, we used the distance (D min ) between a pathogen propagule and its nearest BCA propagule to represent biocontrol potential. Such a distance is critically important for biocontrol outcomes, particularly for biocontrol mechanisms of competition and direct parasitism; it is also important for localised antibiosis or induced resistance.
At the end of each simulation run, D min was calculated and stored for each pathogen propagule; average, median and variance of D min were then calculated over all Fig. 5 An example of the same simulated spatial distribution of clustered pathogen and BCA propagules in a homogenous (a), and heterogeneous with a large (b) or small grain size (c). Patches with green colour favouring the BCA, blue colour favouring the pathogen, and orange colour neutral. The same spatial patch pattern was used for all simulation runs. In this example, the ratio of survival ability between the pathogen and the biocontrol (or vice versa) is 1:0.25 pathogen propagules for each simulation run. ANOVA (analysis of variance) was then used to estimate the proportion of variance in the size (mean/median) and variability (coefficient of variation -CV) of D min that was accounted for by individual simulation variables and their interactions. As the effects of simulation variables (and their interactions) are often statistically significant in simulation studies, statistical significance was not reported specifically. For the simulated heterogeneous areas, the numbers of pathogen and BCA propagules retained were included as covariates in ANOVA; however their inclusion did not affect main ANOVA results and hence only results from ANOVA without them are presented. All simulations and subsequent statistical analyses were carried out in R version 3.4 (R Core Development Team 2019).
Abbreviations ANOVA: Analysis of variance; BCAs: Biocontrol agents; CV: Coefficient of variation; D min : The distance of a pathogen propagule to its nearest BCA neighbour; D min : The average of D min