The “Disease Index” (DI) is a commonly used measure to quantify disease phenotypes. However, the DI is not formally defined, and different definitions are used in the literature. I will use a commonly used definition of the DI, where wilting symptoms are quantified regularly over a defined time, based on a scale of 0–4. Here, one whole number corresponds to a 25% interval of total wilted leaves per plant.

Formally, this DI is defined as:

where w is the number of wilted leaves, and t is the number of total leaves of a single plant. This is multiplied by 4 and rounded with a precision of either 1 or 0.5. One DI is recorded per individual and time point. Independently of the infection method used and the precise research question, one score is commonly recorded per plant per day. Recovery from the infection (a decrease in DI from one time point to the next) is not typical for R. solanacearum infections and therefore not considered in modeling. In experimental inoculations, the total observation time is typically between 10 and 30 days.

Depending on the scientific question one aims to answer in the analysis, the DI can either be used directly as a response variable, or a more suitable response variable can be derived from it.

This section will explain how the disease index and time can be used to analyze different aspects of disease. Figure 1 shows an example disease index data set (Figure 1A) analyzed using the three different methods (Figures 1B–D). These three methods are the area under the disease progression curve, analysis of disease indexes over time and survival analysis. Each of these measures has specific biological implications, outlined in the subsequent paragraphs and each analysis based on these different responses should be interpreted differently, to be able to make conclusions on the underlying biological phenomena. Conducting multiple analyses using distinct measures and comparing their respective outcomes can be a great aid in arriving at biologically meaningful conclusions.

The three response variables, and respective analyses are outlined below.

(I) Area under the disease progression curve (AUDPC): AUDPC is a historically used and well-established measure of disease in plant pathology and bacterial wilt research (Strange and Scott, 2005; Wydra and Beri, 2006; Hadiwiyono et al., 2007; N’Guessan et al., 2012; Lebeau et al., 2013; Meng et al., 2015). Briefly, the disease progression is drawn, and then the area defined by that curve and the x-axis is calculated. To draw the disease progression curve, one takes the mean DI for each time point, and then connects the dots. Once the disease progression curve was drawn, the area under that curve can be calculated. As the AUDPC will increase with time, reasonable comparisons of AUDPC can only be made for experiments of the same observation time (Figure 1B). As the AUDPC provides a linked measure of disease incidence and time, one can consider AUDPC a measure of disease severity, meaning it can be used to summarize the disease progression over time in a single value.

(II) Disease indices over time: in this approach, the DI is scored daily, based on the relative number of wilted leaves per plant, as described above. Subsequently, the DI is used as a response variable for the analysis. Treatment and treatment in relation to the time that has passed since inoculation are used as predictors (Figure 1C) (Franks et al., 2008; Colburn-Clifford et al., 2010; Plener et al., 2010; Jacobs et al., 2012; Monteiro et al., 2012; Ailloud et al., 2016; Mori et al., 2016). An analysis of disease indices over time usually aims to describe the disease progression curve itself as opposed to an analysis of the area encompassed by that curve outlined above. Two different approaches may be taken, either one can attempt to analyze the speed of disease progression, while disease is still actively developing, this can be done in a simple linear framework. However, an alternative could be the application of logistic regression models (Ben et al., 2013), which can accommodate the overall shape of the disease progression curve. Here, I will use linear regression to analyze the speed of disease progression, and more sophisticated, generalized linear models are employed for survival analysis (see below).

(III) Survival: methods from Survival analysis aim to describe the incidence of a certain event within a population over time. In the eponymous situation, that event is “death,” and one is interested in how long it takes until a certain fraction of the population has died, for example to assess efficacy of a treatment. However, any binary event of interest, such as symptomatic vs. asymptomatic hosts or infected vs. uninfected hosts can be used.

For the purpose of survival analysis, disease severity is transformed into a binary scoring, by defining a specific DI value as the event of interest and then recording when each individual reaches this DI value (Plener et al., 2010; Remigi et al., 2011; Poueymiro et al., 2014; Wang et al., 2015). For simplicity, I will continue here using host death as the event of interest, however, in Supplementary Material IIB it is exemplified how these methods can be used to analyze symptom onset.

To be suitable for survival analysis, the disease index dataset needs to be turned into a survival dataset. In a survival dataset, if an individual dies, the day of death is recorded, together with a status indicating “dead.” If an individual survives until the end of the trial, this is recorded as status “alive,” and the last day of observation is recorded as date. Based on the number of subjects in a cohort alive at a given time point, a Kaplan-Meier survival estimate can be calculated and further, the survival over time can be fitted to a specific distribution to proceed with parametric testing (Figure 1D). Survival analysis provides a way to analyze survival of populations upon bacterial challenge. Additionally, survival analysis offers methods to estimate and compare the hazard, which is the risk of dying at a given time point, different populations are exposed to.

Statistical analysis is based on table calculations. Proper table formatting is crucial to be able to properly interface with the R framework of analysis and syntax. In the Supplementary R scripts, formatting is done within R to generate data that conform with the concepts of tidy data (Wickham, 2012).

The R language and environment is maintained by the R foundation and available from R-project.org (R Core Team, 2014). RStudio (RStudio Team, 2015) is an integrated development environment for R, free for academic researchers.

Many of the functions used in the Supplementary scripts are not part of the R base installation. Instead, these are made available through packages and development of new packages and methods is an active area of research. R packages can be retrieved from the Comprehensive R Archive Network (CRAN).

To be able to use the Supplementary R scripts the following packages need to be installed.

For data manipulation and visualization: tidyr (Wickham, 2016c), dplyr (Wickham and Francois, 2016), ggplot2 (Wickham, 2009), modelr (Wickham, 2016a) and stringr (Wickham, 2016b) [these and other helpful packages are part of the tidyverse (Wickham, 2016d)].

Functions from the following packages are used for analysis: MESS (Ekstrøm, 2016), survival (Therneau and Grambsch, 2009; Therneau, 2015a), survcomp (Haibe-Kains et al., 2008; Schroder et al., 2011), rms (Harrel, 2015, 2016), coxme (Therneau, 2015b), lme4 (Bates et al., 2014), lmerTest (Kuznetsova et al., 2016), multcomp (Hothorn et al., 2008), and rcompanion (Mangiafico, 2017).

To able to work with the Supplementary “R Markdown” (.Rmd) files, the packages rmarkdown (Allaire et al., 2015) and knitr (Xie, 2014, 2015, 2016) are required. Some additional table output formatting is done using stargazer (Hlavac, 2015).

The dataset used for Figure generation in this manuscript (S1) is based on a modified version of the dataset published in Wang et al. (2015) (S2) and both the modified and the original datasets are available in the supplement. Specifically, Strain3 and Strain4 were edited to exhibit an unusual behavior in terms of disease development. While this may not be a biologically probable behavior it is helpful in understanding the outcome of different analysis methods. Datasets used in S3 were published in Lowe et al. (2015) and Lowe-Power et al. (2016). The dataset from Ravelomanantsoa & Prior (S4) has not been previously published and is therefore only available in a pseudonymized form. A summary of the dataset used in the main text is given in Table 1.

Many popular types of statistical analysis are based around linear regression. As implied by the name, linear models assume that y and x exhibit a linear relationship.

When performing a linear regression analysis, one sets out to solve a linear function. A simple linear function, with one response variable (y) and one predictor (x) can be written as:

Here, “a” is the intercept, “e” is an error term and “b” is the slope. In all linear models discussed and employed here, the part that one aims to estimate and subsequently compare, corresponds to a or b. The term fitting is used to mean “optimally solving the formula for a, b and e given the values of x and y recorded.” Assuming that y is a single summary measure of disease, and x is used to denote treatments, we aim to estimate b for each individual treatment.

One can estimate the value for b and a, that best fit to the observed data. This “best fit” is optimized to exhibit the least distances to the recorded data. These distances are known as residuals. When a linear model for a single predictor (x) is solved regarding a, one essentially performs a pairwise comparison of y and x. Each distinct y is recorded paired to a single x value, and these pairs are compared. This can be understood visually if x and y are both continuous variables, one will be able to draw a line that determines y based on x. If x is not a continuous variable, but instead, for example different treatments, this becomes harder to visualize as a line, instead this can be thought of individual means that will be obtained for y depending on the value x.

Linear models can be extended to include multiple predictor variables. This leads to an introduction of additional “x” predictors in the formula. Each predictor has its specific intercept. In Figure 2B, estimated “a” for each strain (corresponding to x) based on the AUDPC value (corresponding to y) is shown.

Linear regression models will attempt to estimate “a” for each known x. Statistical analysis of linear regression models can be used to assess if the estimated values for a and more importantly differences between them, are well-supported by the data. Inherently, linear models assume a normal distribution of response variables. All calculations in a linear model to estimate the true mean are performed based on this assumption. An example of this is the use of the mean to assess differences. The mean is not necessarily a useful measure for the center of a distribution, if that distribution is not normal.

Statistical analyses that assume a defined distribution and aim to estimate a certain parameter of that specific distribution, for example the mean in case of a normal distribution, are commonly referred to as “parametric.”

After model fitting the model can be explored regarding its goodness of fit, effectively assessing how close the model is to the data. Briefly, a model that has the least amount of assumptions (predictors) is preferred over one that does not provide a significantly better explanation but includes more assumptions. Model testing and selection, however, are beyond the scope of this manuscript.

In a linear model with a single predictor the estimated value of a for each level of the predictor is the mean of y across all observations for that predictor, relative to a baseline (see next section) (Table 2).

The term “estimate” often causes confusion. One reason why one has to estimate, as opposed to precisely calculate the true mean, is that one should assume that the model is incomplete and the data one has in hand is a random sample of the true population. In every experiment and in every replicate, despite best efforts to control as many factors as possible, there are things that are beyond the experimenters control. The fact that these are not controllable, does not exclude these have an influence. Hence, statistical analysis will always aim to provide a measure of certainty that what one measures is due to a change of known treatment, and not the consequence of other, uncontrolled, factors. Since it is not possible to know the influence of unknown factors, one should estimate the influence of the known treatment, with a specific level of certainty. This is related to significance testing, and is explored below in more detail.

Analysis of variance (ANOVA) is a specific case of linear models, where the means of individual treatments, are compared to the “grand mean.” The grand mean is calculated across all treatments. Subsequently, all treatments are compared to the grand mean. In the language of linear models, this is termed effect coding of the predictor variable, as it compares the effect of each individual predictor relative to the grand mean.

Linear models allow for other types of coding of the predictor variable(s), which may be easier to interpret. I will make use of linear models, where the influence of each strain is analyzed in relation to a reference strain (treatment coding).

Similar to classical linear models, linear mixed effect models (LMM) attempt to fit the observed data to a linear function.

In contrast to classic linear models, LMMs account for two types of predictors, also called effects, those of interest (fixed effects) and those that are not of direct interest (random effects).

One can use an LMM to analyze disease by using a measure of disease as the response variable. Depending on the experimental design and research hypothesis possible fixed effects could for example be: bacterial strain, plant genotype, soil type, competitor strains, or fertilization status.

The accompanying script and dataset provides examples, where LMMs with the AUDPC or repeated measurements of the disease index as response variables are used to assess the influence of different strains on disease development.

The term survival analysis unites a range of methods which aim to characterize time-to-event data across different populations. Classically the event of interest is death of an individual, or failure of a product.

The choice of statistical analysis used should be made based on the underlying research hypothesis. If one is interested in the steepness of the disease progression curve, e.g., because one assumes that treatments will change the speed of disease development, linear regression of repeated DI recordings may be a useful approach. If however, one is interested in the fraction of survivors per timepoint, for example when comparing different plant cultivars in field trials, survival analysis may prove more powerful and relevant. Throughout this manuscript and the supplement I will largely employ generalized linear hypothesis testing, while adjusting for multiple comparisons using Tukey’s method, to assess statistically significant differences (Hothorn et al., 2008).

Useful and informative statistical analysis requires a clear hypothesis that describes the expected outcome. Usually the research hypothesis is that a (specific) change of treatments will lead to a (specific) change of outcome. Statistical testing attempts to lend credence to the research hypothesis via falsification of the null hypothesis. A null-hypothesis matching to the research hypothesis above would be: a change of treatment will lead to no change of outcome. One indicator whether the null hypothesis is true, is the p-value. For the purpose of significance testing one needs to define alpha, the significance threshold. Commonly used is an alpha of 0.05 (5%). If a p-value below alpha is obtained, this is taken as an indication that the null hypothesis is wrong and usually the research hypothesis is accepted instead. But, how does one get a p-value?

To arrive at a p-value, one assumes the null hypothesis to be true. Then, one estimates the true mean for each treatment. Next, one compares the estimated mean to whatever is stated in the null hypothesis. Usually the null hypothesis states either that the true mean is a specific value, or alternatively that the difference between two means obtained for two different treatments is zero. A difference in means is also called an effect.

Based on the difference in means (or the difference from a single mean to a defined value, corresponding to the mean under the null hypothesis), the degrees of freedom and test-specific calculations one arrives at a certain value, known as the test statistic. This test statistic is compared to the distribution of the test statistic. The p-value describes the region of the test distribution, where the obtained test statistic is located. For example, a p-value of 0.01 indicates that 1% of area of the test distribution are further away from the center of the test distribution than the calculated test statistic. In other words, the p-value describes the probability of observing an event as extreme or more extreme than the experimental outcome assuming the null hypothesis and all test specific assumptions were fulfilled.

If the obtained p-value is below alpha, one can reject the null hypothesis. Usually, the research hypothesis is accepted. p-Values are sensitive to a number of factors for example, larger sample sizes will usually decrease the calculated p-value, even if the observed effect is the same, because it is assumed that with more observations, an estimation of the true mean is more precise. Alpha and sample size should be defined before conducting the experiment to minimize the chance of wrongly rejecting the null hypothesis (type I error) and to maximize the chance of finding a true effect (i.e., minimizing the type II error). Power analysis is beyond the scope of this manuscript, but a starting point is the R package “pwr” (Champely, 2016).

Interpretation of the p-value is not necessarily easy and may be confounded by test-specific assumptions. A low p-value does not necessarily indicate a large difference in means; it should be understood to imply that one is unlikely to observe that outcome assuming the null hypothesis was true. Often this is not of direct interest, but instead what is more important is to know the difference in means. A significant test result for a comparatively small effect does not necessarily imply biological significance, but significance within the model used for analysis.

To combine the display of (difference in) means while providing a way to visually assess certainty of the estimate, I will make use of estimated means, and 95% confidence intervals (CIs). In such a display, if the mean of one treatment is not within the 95% CI of another, those two treatments are different assuming a significance level of 0.05 (Krzywinski and Altman, 2013).

Additionally, I will make use of “compact letter displays” to indicate groups of statistical significance. Treatments that are assigned to the same “letter” group are not significantly different. Treatments that are assigned different letters using this method, display a significance difference. For example, a treatment in group “ab” is different from one in group “c.” A treatment in group “ab” is not different from those in group “a” or “b.” However, treatments that are assigned only group “a” are different from those in group “b.” For the purpose of inferring treatment specific differences from a model that includes an interaction and covariates, the treatment coefficients are averaged.

The two types of analysis discussed previously are both based on the disease development of an individual plant. However, if a sufficiently large number of individuals are analyzed as part of an experiment, these may be viewed as a population and disease can be analyzed per population.

Before being able to start with survival analysis, the disease index scorings need to be converted to a survival table. To generate a survival table one needs to check which of the repeated observations for one individual is the first to cross a threshold that defines an event. This time point is recorded, together with a binary status indicator set to the state of “dead.” If an individual never passes the threshold, the last day of observation is recorded together with the status “alive.” Subjects that leave the study before the last day of observation can be recorded as alive on that day, known as right censoring. All variables that specify a fixed or random effect such as Strain, Plant or Batch, should be retained in the survival table. Based on such a table one can now analyze how a population survives over time, and, e.g., compare the impact of different strain treatments on the survival.

The threshold of event generation has to be set by the user, and should not be done without inspecting the data first. As exemplified in Figure 5 changed survival threshold can have dramatic impact on the resulting survival estimates.