Estimating the epidemiology of emerging
 Xylella fastidiosa
 outbreaks in olives

Emerging infectious diseases and invasive species are of significant global concern for economies and public health, and are on the rise due to socioeconomic, ecological, and environmental factors, such as human population density and growth, drug resistance, and host richness (Jones et al., 2008). In plant health, diseases such as European ash dieback (Gross et al., 2014) and sudden oak death Received: 1 May 2020 | Accepted: 11 June 2020 DOI: 10.1111/ppa.13238

have spread to vast parts of the globe, causing death to millions of plants in the process. However, at the start of an outbreak, decision-making processes are less likely to be underpinned by epidemiological information and there is considerable uncertainty around the impact and relative advantage of different disease control measures (Morgan, 2019). Therefore, there is a clear need to develop rapid epidemiological understanding for decision makers.
Xylella fastidiosa is a bacterial plant pathogen with a wide host range of over 500 species (European Food Safety Authority, 2020).
It has had a significant impact on agricultural and horticulture trade around the globe and is economically (Tumber et al., 2014) and socially (Almeida, 2016) costly. X. fastidiosa is a xylem-limited gram-negative bacterium and the recognized agent of a number of severe and economically important diseases, including Pierce's disease of grapevines, citrus variegated chlorosis (CVC), almond leaf scorch, and other disorders of perennial crops and landscape plants (Purcell and Hopkins, 1996). In California alone, it has been estimated that Pierce's disease directly costs the grape-growing industry $104 million per year in reduced yield, management costs, and regulatory costs (Tumber et al., 2014). Once restricted to the Americas, a new invasive strain, X. fastidiosa subsp. pauca ST 53, the causal agent of olive quick decline syndrome (Saponari et al., 2013), was discovered near Gallipoli, southern Italy, in October 2013 (Loconsole et al., 2014). Since the initial outbreak, the disease has spread through a large proportion of the olive trees (Olea europaea) in southern Puglia (Martelli, 2016). A recent estimate suggests there are 4 million unproductive heavily infected or dead olive trees, representing a 10% loss of Italian olive oil, which equates to an economic loss of approximately €390 million (Italia Olivicola, 2019). If left unchecked, it has been estimated that olive-producing countries could lose billions of euros to the disease (Schneider et al., 2020). Other strains and outbreaks of the bacterium have since been discovered in other host plants in Spain, France, other regions of Italy (Tuscany), and Portugal (EFSA Panel on Plant Health, 2019a). These detected outbreaks demonstrate multiple and independent introductions because the strains are distinct (Sicard et al., 2018;Saponari et al., 2019).
The bacterium is generally transmitted by various widespread species of xylem-feeding insect vectors (Homoptera, Auchenorrhyncha) (Cornara et al., 2017). Specifically, in olives in Puglia, Italy, X. fastidiosa subsp. pauca is mainly vectored by the froghopper Philaenus spumarius (Saponari et al., 2014b). Currently, there is no known cure for this deadly disease of olives once the tree is infected (EFSA Panel on Plant Health, 2019b), and the only approaches for control are to destroy the host trees and limit spread by creating buffer zones around them, or managing the insect vector population by insecticides or removal of their weed habitats (EFSA Panel on Plant Health, 2019a). The outbreak in Puglia, Italy, is characterized by extensive leaf scorch and dieback, mainly in olive trees, although other cultivated and wild hosts such as oleander (Nerium oleander), almond (Prunus dulcis), myrtle-leaf milkwort (Polygala myrtifolia), and coastal rosemary (Westringia fruticosa) can also harbour the bacterium (Saponari et al., 2013(Saponari et al., , 2014a.
A particular challenge posed by X. fastidiosa subsp. pauca in olives is its apparently long symptomless period before the appearance of leaf browning, dieback, and eventual death of the plant.
Under laboratory conditions (e.g., artificial inoculation and constant greenhouse environment) symptoms can start to appear after a period of 14 months (c.1.2 years), although this depends on environment and olive cultivar  Health, 2019a). Therefore, ascertaining the symptomless lag under field conditions will help to better inform effective policy.
Host to host transmission is one of the key determinants of the speed at which X. fastidiosa spreads locally and in the landscape (White et al., 2017), but quantification of this process has so far been lacking in olives. It has been shown that within-host bacterial colonization and build-up to infectious levels can take some time (up to a year; Saponari et al., 2017), which is likely to affect vector acquisition over the infection period, with symptomless plants being less infective. However, extrapolation from laboratory studies of between-host transmission to the field is not straightforward because the vectors are polyphagous and show distinct seasonal dynamics and behaviour (Cornara et al., 2017). An additional complexity is that as the symptoms progress, during which the bacterium propagates and greater numbers of the olive branches die back, the olive tree may become less attractive to the insect vector and thus the propensity for between-host infection may drop. Furthermore, the rate at which host plants become desiccated has also not been formally established. In general, the effects of these complexities for predicting X. fastidiosa outbreaks in olives are unknown (but see Daugherty and Almeida, 2019, for outbreaks in grapevines), and simple models that do not allow for this variation may fail to accurately estimate disease dynamics.
In this paper, we aim to address these gaps in our knowledge by fitting epidemiological mathematical models to observed field data, and performing model selection to balance best model fit against model complexity. In doing so, we aim to quantify the underlying epidemiological parameters in the X. fastidiosa subsp. pauca Puglian outbreak that are essential for accurate pest risk assessment and effective control strategy planning. We also reveal the sensitivity of the models to each of these parameters and how they contribute to the overall epidemic.

| Epidemiological model
Our aim is to model the epidemiological dynamics of X. fastidiosa subsp. pauca in Puglian olive groves. Previous models of X. fastidiosa dynamics (White et al., 2017;Abboud et al., 2019;Daugherty and Almeida, 2019) have generally assumed simple tree to tree infection growth dynamics (e.g., logistic-or Gompertz-type growth models of the number or proportion of infected host plants), or used a simplified compartmental Susceptible-Infected (SI) model that omits critical features of the X. fastidiosa epidemiology, such as the symptomless period and host mortality rate (Soubeyrand et al., 2018). Our aim is to build a more realistic model that captures additional details of the infection process using a determin- We model the numbers of susceptible and infected host plants in the olive groves. Model compartments were chosen to represent the major stages of disease progression in olive trees. A susceptible compartment (S t ) included the healthy uninfected trees at time t (years). To reflect the potential for differential infectiveness of olive trees at different stages of infection, the infected (I t ) compartment from the standard framework (Allen, 1994) was subdivided into three sequential compartments representing the numbers of infected trees without symptoms ('symptomless', I A,t ), with symptoms ('symptomatic', I S,t ), and with the majority or all of the branches dead ('desiccated', I D,t ). We expected that, compared to trees with symptoms, there would be lower infectivity from trees without symptoms because they have a low bacterial load  and are probably less likely to transmit X. fastidiosa subsp. pauca to vector insects. Similarly, we expected desiccated trees to have lower infectivity because, although the bacterial load is likely to be high, much of the foliage is withered and brown, which probably makes the trees less attractive to vector insects and may support lower vector densities. Thus, we included in the model two proportion parameters, b A and b D , representing the relative infectiveness of infected symptomless and desiccated trees, respectively, as a proportion of the infectiveness of symptomatic trees.
A further modification to the standard compartmental model was made to characterize observed disease progression. In standard compartmental models (Allen, 1994), individuals may pass to the subsequent compartment immediately after entering their current compartment. In this case, transition times are exponentially distributed for continuous time models and geometrically distributed for discrete time models. However, the disease progression data indicated a delay between appearance of initial symptoms and onset of extensive symptom expression and rapid die back. This probably represents a minimum time period for bacterial load in the trees to build up to lethal levels. Therefore, we extended the basic compartmental model so that infected olive trees would not start to enter the desiccated I D compartment until after a sojourned symptom expression period of τ years had passed (subsequently referred to as desiccation delay parameter). We denote the τ th compartment of symptomatic sojournment by I S,t . The total number of trees with symptoms is given by The epidemiological dynamics at population level are given by capturing the following individual annual transition probabilities for infection of susceptible trees (α t ), symptom appearance in symptomless infected trees (σ), and mortality of infected trees with symptoms (γ), which are given by Here, N t = S t + I A,t + I S,t + I D,t = N is the total number of olive trees in the population. The contact rate β is the expected number of Xylella-transmitting contacts per year from each infective tree with symptoms (Allen, 1994) and, as explained above, proportions b A and b D scale the relative infectiveness of I A,t and I D,t . Parameter T A is the mean duration of the symptomless period in years and T D is the mean time to desiccation after τ years.
It can be seen from this that our model is a generalization of standard compartmental models in the Susceptible-Exposed-Infected-Removed (SEIR) framework (Allen, 1994) whereby the exposed and removed classes are able to infect the susceptible class. In the case that there is no delay (τ = 0), the model collapses to a standard discrete-time compartmental model (Data S1). Setting b A = 0 means that I A is equivalent to an uninfectious exposed E compartment, while setting b D = 0 means that I D is equivalent to an uninfectious R compartment. When b A = 1 or b D = 1, I A and I D are, respectively, equally infectious as I S and thus these compartments can be thought of as a single infected I compartment. It can be seen from this that

| Data to estimate the model parameters
Models were fitted to data from 2 to 3 yearly censuses of symptom

| Model fitting and model comparison
The generalized SI A I S I D model and its simplifications from the SEIR family (SI-like, SIR-like, SEI-like, and SEIR), all with a fixed range of values of the desiccation delay parameter, τ, were fitted to these data in a Bayesian framework. Then the fits of the candidate models were compared using the deviance information criterion (DIC) (Gelman et al., 2013) in order to identify the most parsimonious model structure and value of τ consistent with the data.
For any set of model parameters, we calculated deterministic disease dynamics for 12 years starting from an initial proportion of infected trees without symptoms, I A,0 , which was treated as an additional parameter to be estimated because it was not known when the trees became infected. Twelve years was considered an upper limit on the plausible length of time the plots could have been infected Parameter estimation from these likelihoods was achieved using an adaptive Metropolis-Hastings MCMC sampler in which the variance-covariance matrix of the multivariate normal proposal distribution was tuned to an acceptance rate of 30% (Scheidegger, 2018).
Each model was fitted with a chain of 10 6 iterations after a burn-in of The adaptive MCMC sampler uses a multivariate normal proposal distribution for the parameters, which is not well suited to discrete parameter distributions. Therefore, we ran separate chains for values of the desiccation delay parameter, τ, fixed at integer values between 0 and 6 years. For each model specification (SI A I S I D , SI-like, SIR-like, SEI-like, and SEIR) the optimal value of τ was identified as that giving the minimal DIC, calculated using the variance of deviances in the MCMC chain (Gelman et al., 2013). Following this, the model specifications were then compared based on their DIC values.

| Sensitivity analysis
To assess the relative contribution of the model parameters to X. fastidiosa epidemics, we performed a sensitivity analysis on each fit-

| Model comparison
For all five epidemiological model specifications, MCMC chains had the lowest DIC when the desiccation delay parameter τ (the number of years with symptoms before desiccation begins) was set to 3 years (Table S1). In all cases, setting τ to a different value resulted in large increases in DIC of at least 7 units (Table S1) indicating very low support compared to models with = 3 years. With = 3 years, comparison of the different epidemiological model specifications showed that the SEIR and SEI-like models gave the best fits to the data and were essentially indistinguishable from each other (Table 1). The full SI A I S I D model was slightly less parsimonious and was penalized for estimating two extra parameters that did little to differentiate the generalized model from the simpler SEIlike and SEIR models (Table S1). For example, the SI A I S I D model posterior median estimate of b D (relative infectivity of desiccated trees) recapitulated its uniform prior distribution (Table 2, Figure 1), showing that this parameter has no effect on model fit and supporting the inability to distinguish between SEI-like and SEIR fits (Table 1).
Furthermore, the SI A I S I D model posterior median estimate of b A (relative infectivity of trees without symptoms) was very close to zero (Table 2) causing the trees without symptoms in the SI A I S I D model to function equivalently to the exposed (E) compartment of SEI-like and SEIR models. Supporting this, models in which trees without symptoms were fully infective (SI-like and SIR-like, with b A = 1) gave the poorest fits to the data on the basis of DIC (Table 2).

| Parameter estimates and model dynamics
Posterior parameter estimates from the SEI-like and SEIR model were extremely similar (Table 2) and gave a very good fit to the observed disease progression in the monitored plots (Figures 2 and 3; Figure   S1). The posterior for β suggested that an olive tree with symptoms is able to infect an average of 19.5 other trees in the local population each year (95% credible range 14-26; Table 2, Figure 1). The posterior distribution for T A suggested a mean period without symptoms of around 1.2 years (14 months, 95% credible range 1.0-1.3 years;  Figure 1). Accounting for the 3-year delay before desiccation, the estimated mean time from symptom establishment to desiccation (τ+T D ) was 4.3 years (95% credible range 4.0-4.6 years;   Figure 1). Posterior estimates for the initially infected proportion (I A,0 ) were very low (Table 2, Figure 1) and equate to approximately one tree being initially infected in the median-sized plot.

TA B L E 1 Comparison between the generalized full epidemiological model (SI
There was appreciable posterior correlation between estimates for some parameters, most notably the strong negative correlation be-

| Sensitivity analysis
The results of the sensitivity analysis for the full SI A I I I D model is given in Figure 5 where we compare −50%, −10%, +10%, and +50% changes in individual median posterior parameter values given in  Figure S2) found to be uninformative in the Bayesian model fitting (Table 2).
Interestingly, other model parameters show sensitivity in particular infection compartments, but not in the overall infection as seen in For example, changes in T A cause large effects on the symptomless compartment (I A ), but less on the overall disease dynamics (S).

| D ISCUSS I ON
X. fastidiosa is a devastating plant pathogen that has the potential to infect large numbers and species of cultivated and wild plants One of the most critical epidemiological parameters to measure is the transmission rate, the rate at which susceptible hosts become infected, as this is a measure of disease spread, and strongly influences how difficult achieving control may be (Allen, 1994). We estimated that one infected symptomatic olive infected approximately 19 susceptible hosts per year, suggesting a fast rate of transmission. During the symptomless period, the bacterium begins to replicate within the host, colonizing large parts of the plant xylem . Over this period the host plant is potentially infectious as it is possible for xylem-feeding insect vectors to take up the bacterium (Sicard et al., 2018  Furthermore, given the substantial time to desiccation and that the speed of spread is often determined by the dynamics at the invasion wavefront (Kot and Schaffer, 1986)  C C C C desiccation events. Saponari et al. (2017) found that full desiccation in laboratory studies using saplings occurred after approximately 2 years. However, it is likely that time to desiccation in mature trees in the field is highly dependent on environmental factors, which may induce stress (e.g., water stress), and other factors such as the number and timings of inoculations from multiple vector bites, size, age, and health of the host plant, season, location of inoculum and prevalence of vectors, to name a few. Despite these factors, our time to desiccation estimate had a relatively narrow 95% confidence interval. This may be due to the fact that the data were collected from adjacent locations over the same period.
In Puglia, P. spumarius is the main vector of X. fastidiosa subsp.
pauca (Saponari et al., 2014b;Cornara et al., 2017). Vector abundance is difficult to predict as their distribution is patchy (Morente et al., 2018), and female oviposition is not correlated with olive proximity (Latini et al., 2019). However, it has been shown that P. spumarius occurrence and abundance does increase with olive grove land cover, although evidence of other biotic and abiotic factors affecting insect populations is less well supported (Santoiemma et al., 2019). Based on this, and a lack of data on the abundance of insect vectors, we chose not to explicitly model vector dynamics (but see EFSA Panel on Plant Health, 2019a, for a seasonal vector model). One further caveat with our approach is that our data only relates to the symptoms expressed rather than a quantitative measure of bacterial load. This is especially important because symptomless plants may be either uninfected or infected without symptoms. We have attempted to counteract this in our fitting method by combining the two classes in the model when comparing it to the data with a severity score of zero.
X. fastidiosa has been predicted to be able to establish over wide geographical areas, with the most suitable environment in the ical parameters, such as the transmission rate at different stages of the disease, symptomless period, and time to death, help to achieve informed policy. These will be important for predicting the spread and control of this disease using spatially explicit landscape-scale models (White et al., 2017;Soubeyrand et al., 2018;Abboud et al., 2019; EFSA Panel on Plant Health, 2019a). Furthermore, these epidemiological parameters are difficult and costly to measure experimentally, but models fitted to monitoring data can provide estimates and uncertainty distributions, giving valuable knowledge on X. fastidiosa disease dynamics and, potentially, other plant diseases. Parnell and Andrea Maiorano. We also thank Maria Saponari and Pierfederico La Notte for numerous useful discussions. Finally, we also thank the editor and two anonymous referees for their valuable comments. All authors confirm no conflict of interests.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are available from https://github.com/Quant alab/Xf-NPlan ts-2018/tree/maste r/data or the corresponding author upon request.