Measuring splash-dispersal of a major wheat pathogen in the field

Capacity for dispersal is a fundamental fitness component of plant pathogens. Empirical characterization of plant pathogen dispersal is of prime importance for understanding how plant pathogen populations change in time and space. We measured dispersal of Zymoseptoria tritici in natural environment. Primary disease gradients were produced by rain-splash driven dispersal and subsequent transmission via asexual pycnidiospores from infected source. To achieve this, we inoculated field plots of wheat (Triticum aestivum) with two distinct Z. tritici strains and a 50/50 mixture of the two strains. We measured effective dispersal of the Z. tritici population based on pycnidia counts using automated image analysis. The data were analyzed using a spatially-explicit mathematical model that takes into account the spatial extent of the source. We employed robust bootstrapping methods for statistical testing and adopted a two-dimensional hypotheses test based on the kernel density estimation of the bootstrap distribution of parameter values. Genotyping of re-isolated pathogen strains with strain-specific PCR-reaction further confirmed the conclusions drawn from the phenotypic data. The methodology presented here can be applied to other plant pathosystems. We achieved the first estimates of the dispersal kernel of the pathogen in field conditions. The characteristic spatial scale of dispersal is tens of centimeters – consistent with previous studies in controlled conditions. Our estimation of the dispersal kernel can be used to parameterize epidemiological models that describe spatial-temporal disease dynamics within individual wheat fields. The results have the potential to inform spatially targeted control of crop diseases in the context of precision agriculture.


Introduction 44
Ability to spread within host populations is a fundamental requirement for plant pathogens.
Spatial spread directly inuences the number of new hosts plants that a pathogen can 46 potentially invade and aects the spatial distribution of the pathogen population. For polycyclic diseases, even small dierences in the spread during one infection cycle can 48 result in considerable dierences in epidemic outcomes after several cycles. Understanding the mechanisms and scale of the spread improves our ability to predict and control 50 plant disease epidemics.
Although the epidemiological importance of dispersal has long been recognized (Heald, 52 1913), measurements in eld conditions are rare, for numerous reasons: (i) eld experiments on dispersal are dicult to design and to conduct (McCartney et al., 2006).

54
(ii) They are inherently multidisciplinary spanning the interface between biology and physics. (iii) It is challenging to measure the primary dispersal gradient (i.e., resulting 56 from a single round of pathogen reproduction) independently of secondary dispersal and environmental gradients (Gregory, 1968). (iv) Characterization of ner scale dispersal 58 occurring on sub-meter scales can be dicult in the eld. Existing literature contains many more studies characterizing dispersal of airborne dispersal of plant pathogens com-60 pared to splash dispersal (McCartney et al., 2006): Fitt et al. (1987) found 305 datasets on air-borne dispersal in the eld, only 10 datasets on splash-dispersal in controlled 62 conditions, and no studies on splash-dispersal in the eld. A possible explanation for this disparity is that commercially available spore traps or disease assessments readily 64 capture spatial scales of meters or longer, characteristic of air-borne dispersal, but those are more challenging to apply over sub-meter scales typical for splash-dispersal. 66 Zymoseptoria tritici (formerly Mycosphaerella graminicola) is a major fungal pathogen of wheat (Jørgensen et al., 2014;Dean et al., 2012) that causes septoria tritici blotch borne pycnidiospores (asexual). Pycnidiospores are the main driver of polycyclic epidemics during the wheat growing season. Ascospores are the main source of primary inoculum in the beginning of an epidemic and contribute to epidemic development dur-72 ing the season (Zhan et al., 1998(Zhan et al., , 2000Duvivier et al., 2013). According to some studies, initial inoculum via air-borne ascospores is uniform across the horizontal spatial extent 74 of wheat elds and does not represent a limiting factor for epidemic development (Morais et al., 2016), implying that vertical dispersal from lower to higher leaf layers is more rele-76 vant epidemiologically than horizontal dispersal across the spatial extent of wheat elds.
However, it is not clear how general these conclusions are: it is possible that they hold 78 only in regions with temperate wet climate and intense wheat production dominated by STB-susceptible wheat cultivars, where the studies have been conducted. Yet most of 80 research on splash dispersal of Z. tritici has focused on vertical dispersal of spores from initial infection on basal leaves of seedlings to emerging leaf layers (Shaw, 1987;Lovell 82 et al., 1997;Bannon and Cooke, 1998;Lovell et al., 2004b;Vidal et al., 2018). In this context, the interaction between the pathogen and its host plant has been described as a 84 race, where the pathogen population needs to "climb" up to the next leaf layer before current layer becomes senescent and resources are depleted (Robert et al., 2018). 86 We argue that horizontal dispersal is as important as vertical dispersal, because it greatly inuences the ability of a specic pathogen genotype to take over a eld and 88 therefore it can play a major role in the dynamics of emerging pathogen genotypes adapted to control measures such as fungicides or disease resistance genes in host plant.

90
However, because of insucient empirical knowledge on rain-splash driven horizontal dispersal, it has not been considered in modeling the emergence of new pathogen strains 92 (e.g., (Mikaberidze et al., 2017;Willocquet et al., 2020)). Horizontal expansion of disease foci (Zadoks and van den Bosch, 1994) corresponds to the growth of a genetically uni-94 form pathogen population which is well-adapted control measures applied in a spatially uniform manner (e.g., fungicides or disease resistance genes). For this reason, spatially heterogeneous control strategies have been suggested such as varietal mixtures in individual elds or landscape mosaics (Mundt and Browning, 1985;Brophy and Mundt, 98 1991;Newton et al., 2009;Sapoukhina et al., 2010;Newton and Guy, 2011;Mikaberidze et al., 2014;Djidjou-Demasse et al., 2017;Ben M'Barek et al., 2020). Spatial scale of 100 such strategies can be optimized based on the knowledge on spatial scales of horizontal dispersal of the pathogen. This knowledge is captured by a dispersal kernel function 102 that describes the probability of an individual dispersal event to end up at a certain location relative to the source (Nathan et al., 2012).

104
Dispersal via asexual spores of Z. tritici and of the wheat pathogen Parastagonospora nodorum (formerly Septoria nodorum) has been studied in controlled conditions using 106 inoculations via infected straw or spore suspension combined with articial rain (Brennan et al., 1985;Saint-Jean et al., 2004;Vidal et al., 2017). Bannon and Cooke (1998) 108 studied the eect of wheat-clover intercrop on dispersal from plates by articial rain and detected a reduction of dispersed spores at 15 cm distance from the source. No 110 experiment has been conducted within a host plant canopy in the eld that would allow to parameterize the dispersal kernels associated with splash dispersal.

112
We conducted an experiment to measure splash-dispersal of Z. tritici in eld conditions. We carried out localized articial inoculations with two Z. tritici strains as well 114 as their mixture that contained equal proportions of each strain. Conducive weather conditions allowed us to measure primary dispersal gradients resulting from single dis-116 persal events following the initial articial inoculation. We used automated digital image analysis (Karisto et al., 2018) to estimate sizes of pathogen populations on wheat leaves 118 based on detection of pycnidia (fungal fruiting bodies). In addition, we genotyped the re-isolated strains using strain-specic primers to dierentiate the inoculated pathogen populations from the background natural population. Using a spatially-explicit math-ematical model that took into account the spatial extent of the source area (Karisto 122 et al., 2019b) with the distance of 12.5 cm between the rows. Plots were randomly assigned to four treatments with ve replicates of each treatment as shown in Fig. 1A. The four treat-ments were: (i) inoculation with the single strain ST99CH_1A5 (short identier 1A5), (ii) inoculation with the single strain ST99CH_3D7 (short identier 3D7), (iii) inocula-146 tion with both strains 1A5 and 3D7 (mixed inoculation) and (iv) control (no inoculation).

152
In each plot, a 40 cm-wide area across the middle of the plot was inoculated. Disease levels were measured in the middle of the inoculated area (x 0 = 0 cm) and at eight 154 locations outside of the inoculated area, four on each side, at distances x ±1 = 40 cm, x ±2 = 60 cm, x ±3 = 80 cm and x ±4 = 120 cm from the center of the inoculated area (see

Z. tritici inoculation
Inoculum was prepared by growing the fungus for seven days in yeast-sucrose-broth 162 (https://dx.doi.org/10.17504/protocols.io.mctc2wn). The liquid culture was then ltered, the blastospores were pelleted in centrifuge and re-suspended into sterile water.

164
The washed spore suspension was diluted with water to achieve the concentration of 10 6 spores/ml. For mixed inoculation, the spore suspension was obtained by mixing the 166 same volume of each single strain suspension so that the nal spore concentration was 10 6 spores/ml and each strain was present with the concentration of 5 × 10 5 spores/ml.

168
Finally, we added 0.1 % (v/v) of Tween20 and kept the inoculum suspension on ice until spraying on the same days.
Inoculation was performed by spraying 300 ml of spore suspension onto the inocula-  202 At t 0 , STB incidence was measured at the leaf scale in the following manner. Thirty to forty stems were inspected on each measurement line. The highest diseased leaf layer 204 was recorded for each stem. The leaves lower than that were assumed to be diseased as STB is usually more prevalent in the lower leaf layers (Lovell et al., 2004a ). Additionally, 206 naturally senescent leaf layers at the bottom of the canopy were recorded. In this way, incidence was estimated for all non-senescent leaf layers. After estimating the incidence, 208 eight infected leaves were collected from up to two consecutive leaf layers that had incidence higher than 20% (to avoid removing too much of the inoculum). The collected 210 leaves were then mounted on paper sheets and scanned with 1200 dpi resolution. The resulting images were analyzed using automated image analysis method to measure two 212 aspects of conditional STB intensity that represent the host damage and the pathogen reproduction, as described in Karisto et al. (2018). Host damage was measured as the 214 percentage of leaf area covered by lesions and pathogen reproduction was measured as the number of pycnidia per leaf. The sampled leaf layers at t 0 were the ag leaf layer 216 (F) and the layer immediately below the ag leaf (F-1).
At t 1 , the wheat leaves were already mostly chlorotic and hence the incidence mea-218 surement was not possible in the eld. Instead, we collected at random about 24 leaves from each measurement line. The leaves were taken into the lab and inspected for the presence of pycnidia. Incidence was recorded based on the presence of pycnidia, and the leaves with pycnidia were scanned as described above to quantify the conditional 222 STB intensity. Due to vast chlorosis, the measurement of host damage was considered unreliable and only pathogen reproduction was used in the subsequent analysis. Thus,

224
we measured the conditional disease intensity using only the numbers of pycnidia per leaf.

226
We estimated the number of asexual reproduction and dispersal events between t 0 and t 1 using the following arguments. First, based on the data from Shaw (1990, as revisited 228 in Fig. A1 of Karisto et al. (2018)), the latent period after inoculation was estimated to be longer than 20 days (average daily temperature during rst 19 days was 19 • C).

230
Thus, there was likely no spread from the inoculation areas during the rainy period only 1317 days after inoculation (dai) (Fig. 1F). This was further conrmed by a visual 232 assessment on 8 June (22 dai), when we observed only few tiny lesions and mostly no pycnidia on plants inside the inoculation areas, concluding that substantial spread did 234 not yet occur. Second, at t 0 (14 June, 28 dai) we observed substantial disease ( Fig. 2A) in the inoculation areas. There was no rain during the week preceding t 0 , two strong 236 showers during the night about 36 hours after t 0 and no more rain during the following week. We conclude that there was most likely only one asexual spread event right after 238 t 0 , which caused the disease gradients that we observed at t 1 (48 dai) outside of the inoculation areas.

240
Based on the arguments above, we conclude that the second inoculation was successful: it caused high disease levels observed in the inoculation areas at t 0 , after a latent period population. Following the spatially-explicit method for parameterizing a dispersal kernel (Karisto et al., 2019b), the dispersal process can be described mathematically using two 262 area integrals: one over the source area and the other over the destination area. We used the exponential dispersal kernel that ts well for splash dispersal (Fitt et al., 1987;264 Saint-Jean et al., 2004).
The disease intensity at t 1 in a measurement line at distance x * from the inoculation 266 area is given by 268 B x -1 x +1 x 0 x +4 x +3 x +2 x -2 x This general approach can take into account arbitrary shapes of the two areas of integration D and S. In the specic case of our experimental design, the two areas correspond 270 to the following rectangles: is the destination area that corresponds to a measurement line (the dark rectangles in Fig. 1B),

284
Note that 10 cm width of the measurement lines was practically the smallest possible width that could be achieved in the eld measurements, because the foliage of even a 286 single straw spans more than 10 cm, limiting the spatial resolution of our measurements.
For this reason, we simplied the model by neglecting the width of the measurement 288 lines and assigning the disease intensity values of each measurement line to the middle of the line. With this simplication, the disease gradient is calculated as As implied by Madden et al. (2007) and the analysis of Fitt et al. (1987), dispersal is 292 often modeled as a one-dimensional process. This simplication may lead to incorrect estimates of the dispersal kernel (Karisto et al., 2019b). However, to enable comparisons between our outcomes and the earlier estimates in controlled conditions, and to highlight the distorting eect of this "standard" simplication, we constructed the one-dimensional 296 function describing the dispersal process that follows an exponential kernel: 298 In Eq.
(3) the integral takes sum over the inoculated area along the length of the plot.
We conducted the integration in Eq.
(3) and obtained a simpler expression: where x * takes the values outside of the inoculation area.

302
Estimation of dispersal and transmission parameters. The data used in the analysis were obtained in the following way. First, we conducted STB incidence measurements 304 and we acquired conditional STB intensity measurements from the digital image analysis. Second, the STB intensity values in each measurement line were multiplied by 306 the corresponding STB incidence to obtain the full or unconditional intensity values, thereafter called intensity for the sake of simplicity. Third, we calculated the average 308 full intensity in each measurement line to obtain ve data points (from ve replicates) for in each treatment (at distances x 0 , x ±1 , x ±2 , x ±3 and x ±4 ). These average values over 310 each measurement line were used for tting the model functions. The dispersal gradient functions (Eqs. (1), (2)) were tted to the data to estimate α and I 0 β.

312
To compare the treatments and the two directions, we used the bootstrapping approach. We re-sampled the collected samples with replacement to create a large set of 314 bootstrap samples. Variation in the bootstrap samples reects the variation that we expect to observe if the actual experiment was repeated several times (see for example Davison and Hinkley, 1997). Bootstrapping allowed us to model explicitly the variation related to the incidence counts and the variation related to the leaf collection, inde-318 pendently of each other. This approach also allowed us to assess uncertainties in the parameter estimates without making any assumptions about the distributions of the 320 data or the parameter values.
We created 100 000 bootstrap samples for each measurement line in each replicate in 322 the following manner. First, we simulated the incidence counts on the measurement lines to create a distribution of incidence values. We assumed a population of 82 stems within 324 each measurement line (based on the observed stem density) and simulated a random leaf sampling 100 000 times with all possible incidence values. We recorded the real 326 incidence value each time the simulation produced the observed value. That created a distribution of real incidence values corresponding to our observation. Second, we 328 sampled with replacement the analyzed leaves to generate 100 000 new samples with the original sample size. Third, the mean disease intensity of each bootstrap leaf set was 330 multiplied by an incidence value drawn from the corresponding incidence distribution to obtain the mean intensity for each measurement line. Finally, we grouped these 332 unconditional means over measurement lines into sets representing the ve replicates.
As a result, we obtained 100 000 bootstrap replicates of the entire experiment.

334
The one-dimensional disease gradient in equation (4) was tted to each of the 100 000 bootstrap replicates. The disease gradient function in equation (2) Additionally, we tested dierences between conditions in the parameters αand βI 0 350 simultaneously using a two-dimensional hypothesis test based on the joint distribution of dierences in α and βI 0 (the equidensity test, analogous to Johansson et al., 2014).

352
A kernel density estimate of the joint distribution was obtained to dene the degree of extremity of a point in the two-dimensional parameter space. The point reecting 354 the observed dierence was compared to the joint-distribution of dierences between bootstrap replicates. The observed dierence was considered signicantly dierent from 356 zero, if it was located in a suciently sparse area, such that less than 5% of the bootstrap estimates were located in regions with equal or lower density, analogous to being in the 358 thin tail of a one-dimensional distribution.
We present 95% condence intervals for the parameters derived from the distribution 360 of bootstrap results, i.e. the limits of 2.5th and 97.5th percentile of the distribution.
Additionally, dierences in disease levels between treatments 1A5, 3D7 and mixed inoc-362 ulation were tested at t 0 , x 0 and t 1 , x ±1 with the Kruskal-Wallis test and the pairwise Dunn's posthoc comparison with the Bonferroni correction.

364
Statistics implementation. All data analysis was implemented in Python (versions 3.5.2 and 3.6.0) and the code is provided together with the data (supplement for review).

Results
The second inoculations with strains 1A5 and 3D7, and their mixture were successful: 374 at t 0 we observed increased disease levels in the inoculation areas of all three treatments compared to controls (Fig. 2 A). At time point t 1 , there was a gradient of disease 376 intensity from higher levels at x ±1 to lower levels at x ±4 (Fig. 2 B). Genotyping the re-isolated strains conrmed the successful spread (Appendix A). In total, 4190 plants 378 were inspected in the course of incidence measurements, 2527 leaves were collected and analyzed using the digital image analysis, and 153 isolates were genotyped The entire 380 dataset including raw data, bootstrap replicates, best tting parameter estimates and weather data is available in DATADRYAD (TBA after acceptance in the journal).

382
Estimates of pathogen dispersal parameters Fitting the equation (2) to the observed disease gradients allowed us to estimate pa-384 rameters α (dispersal parameter) and βI 0 (transmission rate × initial intensity at the source). For 1A5 treatment, the estimates of α were very low and estimates of βI 0 were 386 very high, compared to 3D7 and mixed inoculation ( Table 1). The results are biologically unrealistic (as discussed below in subsection Estimates of disease transmission 388 rates) and were likely due to an insucient disease intensity within the inoculation area and consequently a shallow gradient outside the inoculation area ( Fig. 2 A and B). Less successful spread of 1A5 than 3D7 was conrmed by comparing disease intensities between the treatments at t 0 , x 0 and at t 1 , x ±1 . At t 0 , x 0 , the disease intensity was (5) for parameter α: p α = 0.17, parameter βI 0 : p βI 0 = 0.13, Fig. 3 A), nor in mixed 404 inoculation treatment (p 2D = 0.74, p α = 0.60, p βI 0 = 0.95). This similarity between directions suggests isotropic dispersal during the experiment.

406
Considering this isotropy, we estimated the parameters using the dataset that combined the two directions. We established a signicant dierence between 3D7 and mixed 408 inoculation treatments (p 2D = 0.014, p α = 0.020, p βI 0 = 0.018, Fig. 3 B). The dispersal parameter α was higher in mixed inoculation while βI 0 was higher in 3D7 treatment 410 (Table 1). We thus highlight that following the mixed inoculation, the pathogen has dispersed further but caused lower disease transmission compared to 3D7 treatment.

412
Estimates of disease transmission rates In addition to estimating the pathogen dispersal parameter, we estimated the trans-414 mission rate of the disease. The tting yielded estimates of βI 0 , from which we extracted the transmission rate β by dividing βI 0 by estimates of I 0 . We estimated 416 I 0 = 249 pycnidia/leaf for 3D7 and I 0 = 227 pycnidia/leaf for mixed inoculation measured as the mean intensity in ag leaves. Based on these estimates we calculated 418 β = 7.7 (unitless) for 3D7 and β = 5.6 for mixed inoculation treatment.
Parameter estimates for the strain 1A5 were not realistic. It is not biologically plausi-420 ble that spores of the strain 1A5 would disperse by only 3 cm, while spores of the strain 3D7 disperse some 14 cm, because pycnidiospores of the two dierent strains are ex-422 pected to have similar physical properties. Likewise, the estimates of βI 0 for strain 1A5 were unrealistically high (almost 100 000, Table 1). However, we still inferred the param-424 eter βI 0 for the strain 1A5 assuming that the physical process of spore dispersal via rain droplets is the same for the two strains. Under this additional assumption (i.e., by setting 426 α 1A5 = α 3D7 = 13.5 cm), we found no dierence in the spread of the strain 1A5 between the two directions (p = 0.38, parameter βI 0 ). When we combined the data from the two 428 directions, we estimated βI 0 = 349 pycnidia/leaf and β = 3.0 (I 0 = 118 pycnidia/leaf ).
The transmission rate estimates of the strain 1A5 and the mixed inoculation were lower  The proportion of identied 1A5 and 3D7 isolates was lower further away from the inoculation area in the mixed inoculation treatment. We can thus conclude that the 442 strain 1A5 has spread less successfully than the strain 3D7: the proportion of the 3D7 isolates in the 3D7 treatment was higher than the proportion of 1A5 isolates in the 1A5 444 treatment and the same eect was observed in the mixed inoculation treatment.
How good are the simplications? 446 We compared the three dierent models (Eqs. 1, 2, 3) with regard to (i) the accuracy of the estimates and (ii) the computational time.

448
(i) The parameter estimates are virtually the same when using the two two-dimensional models (line destination vs. explicit destination area), which justies our use of the 450 slightly simplied model in Eq. (2). In contrast, the one-dimensional model resulted in substantially higher estimates of α (Table 2).

452
Additionally, the relationship between the population spread and dispersal parameter α is dierent between one-and two-dimensional models. This dierence becomes clear 454 for example when dispersal is described based on mean dispersal distance which is α for one-dimensional exponential kernel but 2α for two-dimensional. For the strain 3D7, reproduction is more relevant for pathogen ecology and evolution. One example of this relevance is that pycnidia numbers are more powerful than PLACL for predicting the 482 future PLACL (Karisto et al., 2018).
The dispersal kernels estimated here correspond to eective dispersal of the pathogen, 484 not to the dispersal of all spores. A dierence between those may arise from densitydependent post-dispersal mortality (Nathan et al., 2012;Klein et al., 2013 from Brennan et al. (1985), but they used a one-dimensional model of dispersal. Our one-dimensional estimates (Table 2) correspond to half-distances of 500 10.5 cm and 17.0 cm for 3D7 and mixed inoculation, respectively (half-distance ≈ 0.69α in 1D). Thus, experiments in controlled conditions appear to translate well into eld 502 but details of the models need to be taken into account for a sound interpretation of results. We computed the one-dimensional estimates solely to compare them to results of 504 previous studies, while we consider the two-dimensional estimates to be more accurate.
While the horizontal dispersal measurements made in previous experimental studies 506 were useful, they are dicult to extrapolate to real eld conditions for two reasons.
First reason is that the modeling approaches did not consider the spatial extent of the 508 source area as we did here but assumed a unique point source (Karisto et al., 2019b).
Second reason is that plant canopy, acting as a barrier that limits the eective dispersal 510 distance, was not always included in the design (e.g. Brennan et al., 1985).
The estimated transmission rates were β 1A5 = 3.0, β 3D7 = 7.7, and β mix = 5.6. The 512 intermediate transmission of the mixed inoculation likely resulted from the contrasting reproductive capacities of the two strains. Strain 1A5 is known to produce fewer and 514 smaller pycnidia than strain 3D7 on cultivar Runal in greenhouse (Stewart et al., 2018), which is further conrmed in eld conditions by our results.

516
STB development depends on weather conditions (Henze et al., 2007). Also, the disease levels within the source and along the disease gradient were measured only on 518 upper leaf layers, but dispersal occurred likely to, and possibly also from, the lower leaf layers as well, which were not included in our estimates of transmission rates. For 520 all these reasons, the relative dierences between transmission rate estimates in our experiment are informative, while their comparison with transmission rates measured in 522 other studies is of limited value.
Dierences in dispersal. According to our analysis, the mixture of the two strains 524 has spread further than the single strain 3D7. The dierence may arise from various sources. Density-dependent eects may have attened the disease gradient in the mixed 526 inoculation treatment, for example if the mixed pathogen population is more sensitive to saturation at high population densities close to the inoculated area due to cross-528 suppression between the two pathogen strains. In this case, the beginning of the gradient, where the population density is highest, would experience stronger saturation after mixed 530 inoculation compared to single strain inoculations. That would result in lower overall transmission but also atter gradient (i.e., longer dispersal), as the tail of the gradient 532 would be relatively stronger after the mixed inoculation compared to the single strain inoculations due to relaxed cross-suppression. Interestingly, these expected patterns 534 matched with our observations in the present study.
Also, sexual reproduction in the mixed inoculation treatment could result in a shal-536 lower gradient due to wind-dispersed ascospores which are expected to have a much higher dispersal distance. The average latent period following infections by ascospores 538 was shown to be longer than the latent period associated with asexual pycnidiospores (Morais et al., 2015), but the distribution of ascospore latent periods might still be wide 540 enough to contribute to our observed gradient (Suert and Thompson, 2018). Ultimate causes of the dierence remain unknown.

542
Sampling distances. The measurement lines were at closest 20±5 cm from the edge of the inoculated area. Measuring the disease intensity closer to the source and even inside the source might have improved the estimates, as the dierences between the gradients are more pronounced close to the source. However, closer to source the reliability of 546 measurements may suer from saturation and also from dispersal via direct contact (Fitt et al., 1989).
excluded from tting, because our aim was to analyze only the newly spread infection 550 and capture the primary disease gradients. Increase of the disease intensity in the source area from t 0 to t 1 was not only due to secondary infections but likely also from latent 552 infections that can become symptomatic after a long time (Karisto et al., 2019a ). An additional reason for excluding data at x 0 is that saturation eects are expected to be 554 strongest there.
Spatially-explicit modeling. The spatially-explicit modeling that accounted for the 556 spatial extent of the source allowed us to parameterize the actual dispersal kernel, while overcoming the practical issues related to a point source. We were able to create suf-558 ciently strong, measurable disease gradients in two out of three treatments. Further extension of the source might have created a suciently strong gradient for strain 1A5 560 too. Overall, using an extended source provides a major improvement for the purpose of measuring dispersal over small scales, where using a point-source is practically impossi-562 ble. The spatially-explicit modeling also allowed to estimate transmission rates β in a biologically meaningful manner. Benets of the two-dimensional modeling compared to 564 the one-dimensional tting are hence undeniable.
Statistical advances. The bootstrapping methods allowed us to test for dierences 566 in our parameter estimates and the associated uncertainties in a robust manner without making any assumptions about underlying distributions. Computationally inten-568 sive bootstrapping with a large number of replicates (100 000) is possible using modern computing resources. Non-parametric bootstrapping is a useful alternative to standard 570 parametric tests, which are often used in biology even when their assumptions are violated. Moreover, we adopted a two-dimensional hypotheses test based on kernel density 572 estimate of the bootstrap parameter distribution. This allowed us conduct joint testing for dierences in two parameters simultaneously, which is likely to be more robust than the combination of two one-dimensional tests for each of the parameters separately. This is demonstrated in the example where a combination of 1D-tests would fail to reject the 576 null-hypothesis (Fig. 3A). The provided source code and raw data will facilitate the application of these methods to other contexts.

578
Implications for epidemiology and control. According to our estimates, the median dispersal distance of new infection from a focal source was around 30 cm. However, the 580 limit of 95% of new infections would extend up to one meter (Karisto et al., 2019b, Eq. A3). In the second dispersal step, the infection will not only spread further outwards, 582 but the secondary spread within the area covered by the rst step would likely cause signicant damage (Shaw and Royle, 1993) and this is when a clearly visible focus of the 584 disease can be seen (Zadoks and van den Bosch, 1994). Thus, when one observes a visible focus in the canopy, the spread has likely already happened at least one step further, Appendix A: Genotyping the strains with PCR Primer design. We designed four primer pairs targeted at each of the two strains.

612
The primers were rst aimed to be fully specic to the target isolates 1A5 and 3D7 within the set of four commonly used lab strains 1A5, 3D7 1E4 and 3D1 (ST99CH_1E4, 614 ST99CH_3D1, Zhan et al., 2002) and second, as specic to the target strain as possible in the eld. Specicity here means that the primers designed for 1A5 should produce 616 an amplicon in PCR only with 1A5 genome and not with other strains. Strain specic primers would allow for a convenient detection of the focal sub-population after the 618 experiment as in a mark-recapture experiment.
To design the primers, we used presence-absence data of predicted genes from Hart-620 mann and Croll (2017). We chose the target regions that were present in the target strain (either 1A5 or 3D7) and absent in the other three isolates (1E4, 3D1, and ei-622 ther 3D7 or 1A5). From those potential targets, we selected ten regions that were least frequently present in the 27 Swiss isolates analyzed by .  Table   A1.

630
Validation of primer specicity. First validation of the primers was done with qPCR among the four strains 3D7, 1A5, 3D1 and 1E4 (Tables A2 and A3, Figures A1 and A2).

632
Successful amplication of the target DNA and no amplication on non-target DNA suggested that each of the eight primer pairs was specic to their target strain among 634 the four strains, indicating successful primer design based on the genomes.  Primers' specicity was then validated in a natural population using multiplex-PCR 636 (Table A4, Table A5) combining each strain-specic primer pair with a primer pair that is specic to Z. tritici generally (Zt_gen primers, Duvivier et al., 2013). Zt_gen provided 638 a positive control for success of the PCR. Primers were tested against 37 natural strains isolated from the control plots of the experiment. Reaction with primers 1A5.9 did not    Thus, it is possible that they were the actual strain 3D7 either left on the eld from previous years of eld experiments or it was a spill-over from the current treatments.