Nonlinear disease tolerance curves reveal distinct components of host responses to viral infection

The ability to tolerate infection is a key component of host defence and offers potential novel therapeutic approaches for infectious diseases. To yield successful targets for therapeutic intervention, it is important that the analytical tools employed to measure disease tolerance are able to capture distinct host responses to infection. Here, we show that commonly used methods that estimate tolerance as a linear relationship should be complemented with more flexible, nonlinear estimates of this relationship which may reveal variation in distinct components such as host vigour, sensitivity to increases in pathogen loads, and the severity of the infection. To illustrate this, we measured the survival of Drosophila melanogaster carrying either a functional or non-functional regulator of the JAK-STAT immune pathway (G9a) when challenged with a range of concentrations of Drosophila C virus (DCV). While classical linear model analyses indicated that G9a affected tolerance only in females, a more powerful nonlinear logistic model showed that G9a mediates viral tolerance to different extents in both sexes. This analysis also revealed that G9a acts by changing the sensitivity to increasing pathogen burdens, but does not reduce the ultimate severity of disease. These results indicate that fitting nonlinear models to host health–pathogen burden relationships may offer better and more detailed estimates of disease tolerance.


Introduction
Disease tolerance is broadly defined as the host's ability to limit damage and maintain health when faced with increasing pathogen burdens, and is a general feature of host responses to infection [1][2][3][4][5]. Identifying mechanisms that underlie host variation in disease tolerance may therefore offer potentially novel therapeutic targets to treat infections [4,[6][7][8], an approach already , the sensitivity (c), slope (s), or severity of the dose-response curve (a) [7,15].
being explored in the context of sepsis [9], HIV [10], influenza [11] and malaria [12]. The key to understanding tolerance is that it cannot be measured by considering host health or pathogen growth separately, but is instead defined by their relationship. This idea is embedded in the original statistical framework of tolerance [13], where it is analysed as a linear reaction norm of host health measured over a range of increasing infectious doses or pathogen burdens. Steep negative slopes for this linear relationship describe groups of hosts that experience a loss in health with increasing loads, while hosts with shallow slopes are able to maintain relatively higher levels of health even as pathogen loads increase, and are therefore relatively more tolerant of infection [1][2][3][4]. While the linear reaction norm approach is intuitive and has been useful in advancing the study of infection tolerance (reviewed in [14]), in some cases it may be hindering our ability to achieve a greater mechanistic understanding of the processes underlying host tolerance of infections. For instance, there is no reason to expect the relationship between host health and pathogen burdens to be linear [15], and assuming so may be misleading. To circumvent this problem, a quadratic term in addition to the linear model has occasionally been used to resolve the nature of this nonlinear relationship [10,16]. Following on from quadratic models, analytical approaches that allow more flexible, nonlinear sigmoid relationships between host health and pathogen burdens-or 'tolerance curves' (figure 1)-have been proposed [4,7,15]. One advantage of a sigmoid tolerance curve is that in addition to the rate of health decline with increasing infection loads (the slope), it also allows other health parameters to be estimated, such as host vigour, host sensitivity to increases in pathogen load, and the ultimate severity of infection, which determines how sick a host can get during infection (figure 1). For example, a recent study of disease tolerance fitted a 4-parameter logistic model to the median survival of Drosophila infected with the bacterial pathogen Listeria monocytogenes, and this allowed to disentangle changes in fly health during infection that arose due to bacterial pathogenesis and host immunopathology [15], which would not have been possible using classical linear analyses.
Despite these analytical advances, many studies continue to infer tolerance phenotypes from separate measures of host health and pathogen burdens measured at a single infectious dose (for example, [11,12,17]). This approach may provide a general indication that groups of hosts differ in their ability to tolerate infection-for example, when differences in survival are not accompanied by changes in pathogen loads-but they are less useful at describing the rate of health loss with increasing pathogen burdens (the very definition of tolerance), and are also not informative about tolerance at varying infectious doses. This multitude of analytical methods also makes it difficult to draw general conclusions across studies in different species about the ability of hosts to tolerate infection, and how this may vary with genotypic and environmental variation [5,18].
Using both linear and nonlinear analyses, here we examine the tolerance response to a systemic viral infection in Drosophila melanogaster, which has been used extensively as a model host to dissect the mechanisms underlying tolerance of bacterial and viral infections [15,17,[19][20][21]. Drosophila melanogaster infected systemically with Drosophila C Virus (DCV) develops pathology in the reproductive and digestive organs, severe abdominal swelling due to enlargement of the crop and eventually death [22][23][24]. An epigenetic regulator of the JAK-STAT pathway, G9a, was previously identified as a mediator of tolerance to RNA virus infection in D. melanogaster [17]. When exposed to a single lethal dose of DCV, fly mutants with a dysfunctional G9a showed higher mortality than those with a functional G9a, even though there was no difference in the viral loads of the two lines measured at a single time point. This work was notable in providing one of the first examples of immune-mediated tolerance of RNA virus infection, but because it focused on the functional basis of hypersensitivity to DCV, it was designed in a way that did not allow a comprehensive assessment of tolerance. First, infection tolerance was extrapolated from separate analyses of host mortality and viral loads, providing limited information on the role of G9a on host vigour, sensitivity to viral growth or the severity of infection that flies can withstand. Further, tolerance was only measured at a single viral dose, making it unclear if the observed tolerance phenotype is dose-specific. Finally, this work only assessed the tolerance phenotype of female flies challenged with a viral infection. Given the prevalence of sexual dimorphism in immunity [25][26][27], and its expected epidemiological and evolutionary consequences [28,29] it is important to test if infection tolerance may also vary between sexes.
We employed systemic infections in both males and females of two Drosophila lines with identical genetic backgrounds, differing only in having a functional or non-functional G9a (G9a +/+ and G9 −/− ). We challenged these flies with a range of DCV doses and then quantified their tolerance responses using both the slope of linear reaction norm and a nonlinear sigmoid model. This allowed us to measure infection tolerance in the most comprehensive way, identifying which components of infection tolerance are affected by a single regulator of fly immunity (G9a), while also providing a useful comparison of current methodology to estimate components of infection tolerance.

The magnitude of G9a-mediated antiviral protection is dose dependent
Following infection with a range of doses of DCV, we found that overall G9a −/− flies showed significantly higher mortality compared with G9a +/+ (figure 2a,b and table 1) in line with previously reported effects of this gene on fly survival [17]. However, we found that G9a +/+ and G9a −/− responded differently to each viral dose, and that the magnitude of the survival benefit of having a functional G9a varied with the infectious dose of DCV (table 1, fly line × dose interaction). Notably, mortality was similar between G9a +/+ and G9a −/− when challenged with the highest dose (10 9 ), indicating that the protective effect of a functional G9a is no longer observed when flies were challenged with very high doses of DCV.

G9a +/+ and G9a −/− flies exhibit similar Drosophila C virus viral loads at all infection doses
Overall, female flies achieved higher viral loads measured 5 days post-infection compared to males (table 2, sex effect). Viral load increased in a dose-dependent manner in both lines (table 2, viral dose  effect) and did not differ between the two fly lines across all infection doses (figure 2c,d, table 2, line × dose interaction).
2.3. The slope of the linear reaction norm suggests that G9a-mediated tolerance is sex-specific Given that we found a significant positive relationship between DCV infection dose and the viral titre (males: p < 0.0001, r 2 = 0.56; females: p < 0.0001, r 2 = 0.64; figure 2c,d). Differences in tolerance between G9a +/+ and G9a −/− are indicated by a significant interaction between the viral dose and the fly line for survival, which reflects that the rate at which survival changes with increasing viral doses (tolerance) varies between fly lines. We detected a significant interaction in females (table 3 and figure 3b), where G9a +/+ females showed higher tolerance (Slope: −1.2 ± 0.4) compared to G9a −/− females (Slope: −1.8 ± 0.3). In males, however, no significant difference in slopes was detected between G9a +/+ and G9a −/− (table 3), suggesting that per cent decline in health of G9a +/+ males was comparable to that of G9a −/− (figure 2a).

Nonlinear models reveal that G9a affects the sensitivity to infection but not its severity
We then analysed the relationship between fly survival and viral dose using a nonlinear 4-parameter logistic model [7,15,30]. Using this nonlinear model allowed us to assess how G9a affected the sensitivity of flies during infection and the ultimate severity of DCV infection in both functional and deficient versions of this regulator of the JAK-STAT pathway (figure 1). A comparison of the overall fit of the curves showed that G9a +/+ and G9a −/− genotypes have distinct tolerance profiles during DCV infection (p < 0.0001, males: F 2,196 = 14.8, females: F 2,196 = 50.71). We further used this model and equation   extract parameters-'c' (inflection point: sensitivity) and 'a' (disease severity) (figure 1). In contrast to the analysis assuming a linear relationship, we found that both male and female G9a −/− flies differed from G9a +/+ in their ability to tolerate DCV (figure 3c,d and table 4). This nonlinear analysis showed that G9a +/+ flies have a significantly higher inflection point compared to G9a −/− flies, suggesting that G9a +/+ flies are less sensitive (or more tolerant) to increasing viral doses. The sensitivity to infection differed by 70-fold between G9a +/+ and G9a −/− in females (p < 0.0001) while we detected a 10-fold difference in males (p < 0.0001) (figure 3e). We did not find any difference between G9a +/+ and G9a −/− in the severity of infection (

Lack of G9a leads to sex-specific expression of the JAK-STAT ligand upd3
The differences between G9a −/− males and females in the sensitivity to increasing viral loads (figure 3e) prompted us to investigate if this may be due to sex-specific regulation of fly immunity by G9a. As G9a is known to be an epigenetic regulator of the JAK-STAT pathway, we measured the expression of JAK-STAT pathway genes in flies receiving the highest viral concentration, 5 days following DCV infection (the same day on which viral loads were quantified). We chose this dose since we detected high viral titres at 5 DPI, and were therefore more likely to detect changes in immune gene expression. We measured the expression of an extracellular ligand upd3 which is released upon stress or infection, its transmembrane receptor domeless, a suppressor of JAK-STAT pathway, i.e. socs36E, and a stress response element totA. Both socs36E and totA are downstream regulators of JAK-STAT pathway. Compared to flies with a functional G9a, G9a −/− mutants showed a significant increase in the expression of the JAK-STAT ligand upd3, and this effect was stronger in male flies (figure 4 and table 5). Males showed generally higher expression of the JAK-STAT receptor domeless, although this effect was independent of G9a status (figure 4 and table 5). G9a −/− flies also showed significantly higher expression of the negative regulator of JAK-STAT, socs36E, but this effect did not differ between males and females. Finally, we found no effect of either sex or G9a on the expression or turandotA (totA), which is commonly expressed in response to stress.

Discussion
Targeting mechanisms that promote greater tolerance of infection is a promising addition to our current arsenal of strategies to fight infection [3,[6][7][8]. However, if this is to be a successful undertaking, it is crucial that the analytical tools we use to measure tolerance are able to capture distinct host responses during infection. A major aim of this study was to evaluate two analytical methods to measure tolerance (linear reaction norms and nonlinear curves), in order to assess the benefit and limitations of each approach. Both linear and nonlinear models consistently showed that G9a +/+ females have higher tolerance than G9a −/− females, when measured across a range of DCV doses. However, the two models presented  different results in the case of males. Additionally, we observed similar mortality between G9a +/+ and G9a −/− at the highest dose, highlighting the importance of studying infection at a range of doses, as artificially high infectious doses can mask these protective effects. Differences in infection tolerance between groups of hosts are commonly extrapolated from single infection doses and by performing independent analysis for host survival and pathogen burdens [11,12,17]. These experiments are useful in detecting the effect of specific mechanisms underlying host tolerance, but the limitations of single-dose tolerance experiments arise because a host's ability to limit the damage, and therefore tolerate infection, is not necessarily independent of the within-host pathogen burden it suffers. The approach arising from the evolutionary ecology of infection, which measures tolerance as a linear reaction norm, provides additional information, such as an estimate of general vigour from the intercept, while its slope gives an estimate of the rate of decline in the health [1,14]. While useful, there are at least two limitations to interpreting infection tolerance as a linear reaction norm. First, there is no biological requisite for the relationship between pathology and pathogen burdens to be linear       In this regard, using nonlinear sigmoid models offers greater flexibility to fit a range of possible relationships (including a linear approximation) between host health and pathogen loads [15]. Moreover, nonlinear tolerance curves provide additional information on several components of host responses to infection, such as the sensitivity and severity of infection. Differences between groups of hosts in any of the parameters extracted from a nonlinear model (slope, sensitivity, severity) may therefore reflect distinct underlying mechanisms that either promote greater damage prevention or increase damage repair during infection [15]. Employing such methodology to a range of host-pathogen systems may therefore yield useful targets for therapeutic interventions that increase host disease tolerance [7].
We chose to use viral dose instead of viral load as the measure of pathogen burden in our analysis of tolerance. This approach has the advantage that the covariate (dose) is measured without error, which is one of the assumptions of ANCOVA, and avoids common problems with underestimation of slopes in ANCOVAS when there is experimental variance in the covariate [33]. Using dose as a covariate is also useful in the cases where load is not normally distributed or when the relationship between pathogen dose and pathogen load is condition-dependent, which could lead to biased tolerance estimates [20]. In the situations when the correlation between the viral dose and viral load is weak, it would be ideal to use pathogen load in estimating tolerance.
Despite all the stated advantages, one potential drawback of using nonlinear sigmoid models to study tolerance is that the slope of health decline is measured over a very small range of pathogen doses. Accurate estimates of that slope would therefore require many estimates of host health within a very narrow range of pathogen burdens, which may be experimentally challenging. For this reason, in the current analysis we fixed the slope to −1, which also increased the statistical power to estimate the two parameters of interest (sensitivity and severity). In this respect, we propose that using the linear reaction norm approach could be useful if the main variable of interest is the rate at which hosts lose health, which is estimated along the full range of pathogen burdens. Nonlinear approaches are more informative if the parameter of interest is the dose that causes the greatest shift in host health, or if the question relates to how hosts may differ in the ultimate severity of an infection.
A subsidiary interest of this work was to quantify G9a-mediated tolerance of DCV in both male and female Drosophila. Our results support previous work using a single dose of DCV, which found that the lower survival of G9a −/− flies was not attributed to differences in their viral load [17]. Applying the classical linear reaction norm approach, the linear slope showed no difference in the rate of survival with increasing DCV doses, while the nonlinear fit indicated significant differences between G9a +/+ and G9a −/− males in their sensitivity to infection. Further analysis of the sensitivity and severity of infection parameters extracted from the nonlinear model also showed that the effect of G9a in mediating host sensitivity to DCV was greater in females. Our results therefore show that dysfunctional G9a accelerates the onset of infection-associated mortality (affecting sensitivity), without causing substantial effects on the disease severity ultimately experienced by infected flies. It is important to note that while previous work also found G9a −/− females were hypersensitive to DCV infection, this was only assessed by measuring survival at one dose and infection tolerance was not measured in male flies [17].
While males and females are generally susceptible to the same pathogens, sexual dimorphism in immunity is present in a wide range of species [25,34,35], and sex differences in infection tolerance are documented for all classes of viral, bacterial, fungal and parasitic infections (see [28] for a review). Differences between males and females in the ability to tolerate infection will directly impact on the pathogen loads within hosts, and as a consequence, also affect host shedding of pathogen transmission stages [36,37]. Sexually dimorphic tolerance is therefore predicted to generate potentially important heterogeneity in pathogen spread and evolution [28].
The mechanisms underlying sex-specific G9a-mediated effects on infection tolerance are not clear. G9a is a histone methyltransferase [38,39], and the protective effect of G9a during viral infection has been previously shown to be driven by the regulation of the JAK-STAT pathway [17]. Specifically, G9a is known to alter the methylation state of the positive regulator domeless, the negative regulator socs36E and downstream pathway components (totA, vir-1) of the JAK-STAT pathway [39]. Fly mutants without a functional G9a show increased mortality due to immunopathology caused by excessive expression of these downstream genes [17]. Notably, G9a, domeless and its ligands upd1, upd2 and upd3 are all found on the X-chromosome in Drosophila [40,41]. We hypothesized that these X-linked regulators of fly innate immunity could underlie the sexually dimorphic tolerance response we observed. However, while we detected sex differences in the expression of domeless (figure 4), these effects were not a consequence of G9a function. Instead we found the largest sex-specific effect on the expression of the upd3 ( figure 4). Previous work had shown that upd3 expression is higher in G9a −/− flies, resulting in immunopathology [17]. Our work shows that these effects differ between sexes and are especially strong in males. In summary, we show that G9a mediates tolerance during infection with DCV over a large range of viral doses and that this response differs between males and females. Our study therefore places emphasis on the importance of incorporating both males and females in studies of immunity. Our results also stress that conclusions about disease tolerance will vary according to the method used to estimate the relationship between host health and pathogen burdens. We suggest that a combination of linear and nonlinear models is ideal to achieve estimates of the rate of decline in health following infection but also of subtler components of hosts' responses, such as the sensitivity to increases in pathogen loads and the ultimate severity of disease experienced during infection. Further, these methods are applicable not only to the declines in health experienced during infection, but could in principle be applied to other diseases, such as cancer [42]. In general, our understanding of host responses to disease requires a more complete assessment of resistance and tolerance mechanisms across a range of genetic and environmental contexts.

Fly stocks
Experiments were conducted on Drosophila melanogaster mutants G9a +/+ and G9a −/− (also known as G9a DD2 ), kindly provided by R. van Rij (Radboud University, Nijmegen, NL). G9a −/− was originally constructed by excision of P-element KG01242 from the 5 UTR of the gene on the X-chromosome [17,39]. As a control for the mutant line, a P-element excision line was generated in the same genetic background with a functional phenotype of G9a, referred to as G9a +/+ . Both fly lines were maintained on standard Lewis Cornmeal medium under standard laboratory conditions at 25°C, 12 h : 12 h light : dark cycle.

Generation of experimental flies
To set up the experiment, we collected eggs from 15 males and 15 females of each line, kept in vials containing 6 ml Lewis medium supplemented with dry yeast to encourage egg laying. We set up 10 replicate vials per genotype. Flies were left to oviposit in the vials for 24 h before being removed and allowing eggs to develop under standard rearing conditions. To control the larval density of these flies, egg density was maintained between 80-100 eggs per vial by removing excess eggs when required. Flies emerging from these vials were challenged with DCV and their mortality post-infection and viral load were measured.

Virus preparation and titration
Drosophila C virus (DCV) is a ssRNA virus of the family Dicistroviridae. We obtained a viral stock by amplifying DCV in Drosophila S2 cells as described previously [43]. Cell homogenate containing DCV was passed through a sucrose cushion, and the resulting pellet was suspended in 10 mM Tris-HCl (pH 7.3). The virus stock was stored in small aliquots at −80°C until further use. To estimate the viral infection dose, we measured the absolute quantity of virus in the stock culture using quantitative Real-Time PCR (qRT-PCR). Briefly, an aliquot of this virus stock was taken to extract RNA using TRI reagent. Total RNA was extracted from virus stock using Tri Reagent (Ambion), using a Direct-zol RNA miniprep kit, which includes a DNAse step (Zymo Research), reverse-transcribed with M-MLV reverse transcriptase (Promega) and random hexamer primers, and then diluted 1 : 2 with nuclease-free water. The cDNA synthesis was carried out at 37°C for one hour following initial ligation of 5 min at 70°C. This cDNA was then serially diluted 10-fold until dilution 10 −10 . qRT-PCR was performed on an Applied Biosystems StepOnePlus system using Fast SYBR Green Master Mix (Applied Biosystems) and DCV primers containing 5 AT-rich flaps DCV_Forward: 5 AATAAATCATAAGCCACTGTGATTGATACAACAGAC 3 , DCV_Reverse: 5 AATAAATCATAAGAAGCACGATACTTCTTCCAAACC 3 ). The reaction conditions used were: 95°C for 2 min followed by 40 cycles of 95°C for 10 s and 60°C for 30 s. The first dilution where no viral cDNA was detected was taken as zero, and stock viral quantity was back-calculated from this reference point. Using this method, viral copy number in the stock was estimated to be around 10 9 DCV IU ml −1 .

Virus infection
We exposed experimental flies to 5 viral concentrations-0 (control), 10 6 , 10 7 , 10 8 and 10 9 DCV IU ml intra-thoracic pricking with a needle immersed in DCV suspension under light CO 2 anaesthesia. The effect of injury caused by pricking was controlled by including sham-infections performed using a needle dipped in sterile 10 mM Tris-HCl (pH 7.3).

Experimental set-up 4.5.1. Post-infection survival
We systemically infected 3-4-day-old adults by intra-thoracic pricking with DCV. For each fly line × sex × dose combination, we infected 20 replicate individual flies (400 flies in total). Following infection, flies were housed individually in a vial and flies were monitored daily for mortality. Flies were transferred onto fresh medium once a week. Survivorship was followed for 25 days post-infection and flies that were still alive at the end were censored in the analysis.

Viral load
An additional five individuals for each genotype × sex × dose combination were infected as described above to quantify differences in viral growth following infection, using the expression of DCV RNA. Flies were individually placed in TRI reagent (Ambion) following five days of infection (5 DPI) and stored at −80°C. Total RNA was extracted from flies homogenized in Tri Reagent and exactly the same protocol was followed as described previously in §4.3 (Virus preparation and titration). The expression of DCV transcripts was normalized to transcript levels of the housekeeping gene rp49 (Dmel_rp49 Forward: 5 ATGCTAAGCTGTCGCACAAATG 3 ; Dmel_rp49 Reverse: 5 GTTCGATCCGTAACCGATGT 3 ) and expressed as fold change relative to the control flies, calculated as 2 − Ct [44].

Gene expression
Using the same cDNA samples described above, we measured the expression of immune genes likely to be affected by a G9a deletion, to test for genotype-specific and/or sex-specific differences upon infection. Given that G9a is a regulator of the JAK-STAT pathway, we focused on measuring the expression of the extracellular ligand upd3, its transmembrane receptor domeless, and downstream regulators socs36E, and totA. Gene expression was quantified in flies exposed to the highest dose treatment (10 9 ) relative to their expression in uninfected controls, as we expected the strongest effects on gene expression to be detected under elevated viral challenge.

Statistical analysis
To analyse differences in survival of each host line according to the virus dose they were challenged with, we analysed post-infection survival data using a Cox proportional hazards model with 'fly line', 'dose' and 'sex' and their interactions, fitted as fixed effects. To assess differences in viral titres following infection with infection doses, we analysed the Log 10 viral titre measured at 5 DPI using ANCOVA and fitted 'fly line' and 'sex' as categorical fixed effects, 'dose' as a continuous covariate and their interactions as fixed effects. All the non-significant interactions were removed, and the minimal model was used to estimate correlation coefficient between viral load and viral dose.
We tested for differences in infection tolerance between fly lines in two ways. First, we used the classical approach and analysed tolerance as a linear reaction norm [20,45,46]. A general linear model was fitted separately for males and females to study the effect of 'viral dose' and 'fly line' on fly health, i.e. fly survival (DPI). In this analysis 'viral dose' was treated as a continuous covariate and 'fly line' as a categorical factor. Given that infection tolerance refers to the rate at which hosts lose health with increasing viral loads, the analysis of interest is how the slope of the two fly lines differs, that is, if survival is significantly affected by the 'fly line × viral load' interaction.
As a second approach to quantify infection tolerance, we fitted a nonlinear 4-parameter logistic model (equation (4.1)) which is commonly used to assess dose-response curves [15,30,47] Table 6. Summary of AIC results for models fitting the data for estimating tolerance. For each model fitted to the dose response data, we report the AICc as a measure of the goodness of fit to the data. The AICc weight can be interpreted as the probability that model is the best model among the fitted models. AICc were calculated as AICc = AICc(i) − AICc min. . R 2 shows the proportion of the variance response that is explained by each model. For ease of comparison between males and females, we opted to use the 4 parameter logistic model, as this was the best model for female data, and in males the difference between the AICc for the 3-parameter and 4 parameter models was negligible and the R 2 identical. Using GraphPad Prism v. 6.0, we therefore fitted a 4-parameter logistic model described by the following equation: where a is the level of disease severity, b is a host's general vigour, c is the sensitivity to increases in pathogen dose (in this case, the pathogen dose that results in the level of health halfway between the min and max health, similar to EC 50 ), s is the slope of the logistic model, and x is the pathogen dose.
Since the experiment was terminated 25 days post-infection, we constrained the upper limit of the model (host vigour) to 25 (DPI) (see also [15]). To test if these tolerance curves differed between fly lines, we compared the overall fit of the curves using extra sum-of-squares F-test [48]. Using these fitted models, we also examined if lines differed in their sensitivity and disease severity, after fixing the slope to −1. To test for differences in the severity and sensitivity to infection between male/female groups, or between G9a lines, we carried out two-way ANOVA with line and sex as fixed effects, using Tukey's HSD as a post hoc analysis of pairwise comparisons. Unless otherwise stated, all analyses were carried out in JMP 12 (SAS).
Data accessibility. All data are available at the Dryad Digital Repository: http://dx.doi.org/10.5061/dryad.7b2p6 [49]. Authors' contributions. V.G. carried out the experimental work; V.G. and P.F.V. conceived the study, designed the experiment, carried out statistical analyses and wrote the manuscript. Both authors gave final approval for publication.
Competing interests. We declare we have no competing interests. Funding. This work was supported by a strategic award from the Wellcome Trust for the Centre for Immunity, Infection and Evolution (http://ciie.bio.ed.ac.uk; grant reference no. 095831), and by a Society in Science-Branco Weiss fellowship (http://www.society-in-science.org), both awarded to P.F.V. P.F.V. is also supported by a Chancellor's Fellowship (School of Biological Sciences, University of Edinburgh).