The metabolome as a link in the genotype-phenotype map for peroxide resistance in the fruit fly, Drosophila melanogaster

Background Genetic association studies that seek to explain the inheritance of complex traits typically fail to explain a majority of the heritability of the trait under study. Thus, we are left with a gap in the map from genotype to phenotype. Several approaches have been used to fill this gap, including those that attempt to map endophenotype such as the transcriptome, proteome or metabolome, that underlie complex traits. Here we used metabolomics to explore the nature of genetic variation for hydrogen peroxide (H2O2) resistance in the sequenced inbred Drosophila Genetic Reference Panel (DGRP). Results We first studied genetic variation for H2O2 resistance in 179 DGRP lines and along with identifying the insulin signaling modulator u-shaped and several regulators of feeding behavior, we estimate that a substantial amount of phenotypic variation can be explained by a polygenic model of genetic variation. We then profiled a portion of the aqueous metabolome in subsets of eight ‘high resistance’ lines and eight ‘low resistance’ lines. We used these lines to represent collections of genotypes that were either resistant or sensitive to the stressor, effectively modeling a discrete trait. Across the range of genotypes in both populations, flies exhibited surprising consistency in their metabolomic signature of resistance. Importantly, the resistance phenotype of these flies was more easily distinguished by their metabolome profiles than by their genotypes. Furthermore, we found a metabolic response to H2O2 in sensitive, but not in resistant genotypes. Metabolomic data further implicated at least two pathways, glycogen and folate metabolism, as determinants of sensitivity to H2O2. We also discovered a confounding effect of feeding behavior on assays involving supplemented food. Conclusions This work suggests that the metabolome can be a point of convergence for genetic variation influencing complex traits, and can efficiently elucidate mechanisms underlying trait variation.


Background
Phenotypic variation among individuals in a population arises from variation in the genotype, the environment, and the interaction between the two. Genetic variation is a major determinant of many complex traits and, while numerous genetic association studies have failed to explain substantial portions of heritable variation in a given complex trait, use of highly polygenic models have closed this gap considerably [1]. However, highly polygenic models do not easily allow us to identify gene-level associations with traits, and may not yield mechanistic insight into the pathways that shape trait variation. Genetic variation ultimately affects phenotype though the effect of genes on downstream 'endophenotypes'the epigenome, transcriptome, proteome, metabolome and microbiome [2,3]. Several authors have proposed that these endophenotypes, and the metabolome in particular, may serve as a powerful tool in mapping genotype to phenotype, as well as a source of information upon which to construct mechanistic hypotheses [3][4][5][6][7]. In this study, we explore the possibility that genetic effects on phenotypes are filtered through the profile of small molecules, the metabolome, that function downstream of genotype but upstream of phenotype.
To explore the potential of metabolomic profiling to bridge the genotype-phenotypic gap and to identify underlying mechanisms of natural variation, here we study resistance to peroxide (H 2 O 2 ) stress in a fruit fly model of genetic variation, the Drosophila Genetic Reference Panel (DGRP) [22]. This system provides an ideal model to study the ability of metabolic profiling to bridge the genotype-phenotype gap. First, H 2 O 2 resistance assays can be performed on hundreds of flies in parallel, and the resulting survival data can be analyzed within a rigorous statistical framework [23][24][25]. Second, many association studies have examined genetic variation for survival in Drosophila [26][27][28], including in response to oxidative stressors [25,29]. The DGRP, a set of genetically diverse and fully sequenced inbred lines, now enables labs around the world to quickly identify loci associated with any trait of interest [26,30]. Third, numerous studies over the past decade have shown that metabolite profiles in flies are highly sensitive to variation due to genotype and environment [9,[31][32][33][34][35][36]. And so, a metabolomic study of H 2 O 2 resistance in the DGRP could be used to simultaneously measure the metabolomic response to stress and its association with resistance phenotype within a genetically diverse population. We hypothesize that resistance to H 2 O 2 among diverse genotypes may associate with metabolic pathways that mediate the resistance phenotype. This hypothesis supposes that we would find resistance-associated metabolic pathways among a panel of genetically diverse lines. Alternatively, metabolites could have genotype-specific associations to an extent that we may fail to detect trait associations with metabolic pathways across the genotypes measured here. We note that our experimental design does not rule out genotype-specific metabolic activity as a mediator of stress resistance. Rather, we seek to investigate the potential for trait-associated genotypes to converge on the metabolome.
We present an analysis of survival time and metabolic profiles in flies from DGRP lines held on H 2 O 2 or control food. We show that survival is heritable, and that 34.2% of the total variance could be explained by additive effects of the known genetic variants in our study population. By mapping genetic variation associated with resistance, we identified at least two genes associated with lifespan on H 2 O 2 food, including NPF and u-shaped (ush). To simultaneously assess the contribution of the metabolome to phenotypic variation, we profiled the metabolome of a subset of genotypes. By comparing highly resistant with highly sensitive lines, we modeled H 2 O 2 resistance as a binary trait. Consistent with the hypothesis that genotypes could converge on common metabolomic endophenotypes, we found a consistent metabolomic signature of resistance to H 2 O 2 . Multivariate analysis of metabolome variation across these genotypes allows us to distinguish resistant from sensitive lines, even in samples of flies not exposed to H 2 O 2 food. We compare the potential of genotype and metabolome to explain trait variation among the resistant and sensitive lines. Whereas multivariate clustering based on the metabolome leads to clear distinctions between resistant and sensitive genotypes, similar methods applied to genetic variation alone fail to differentiate their resistance phenotype. These results suggest that a variety of even quite diverse genotypes that share a phenotype may do so by converging on a similar metabolomic endophenotype.
Additionally, using univariate analysis of individual metabolite features, we found glycogen and folate metabolism are associated with stress resistance and validate this analysis by showing that metabolite feeding or genetic manipulation of candidate pathways both affect survival on H 2 O 2 food. Finally, we found a strong effect of H 2 O 2 on feeding behavior, suggesting that variation in survival of flies on food supplemented with H 2 O 2 could be explained in part by variation in starvation resistance. Our results suggest that starvation or nutrient assimilation might be the underlying cause of mortality in historical assays where the stressor has been administered to Drosophila in the food.

Results
Resistance to peroxide food within the DGRP We found substantial variation among 179 DGRP lines for survival of adult females on H 2 O 2 food. We used a mixed model to study the variation in mean lifespan and, in addition to significant effects of genotype (random effect, likelihood ratio χ 2 = 36.2, df = 1, P = 1.7 × 10 − 9 , Methods), we found that weight was a significant predictor, with larger flies surviving longer on H 2 O 2 (fixed effect, β = 0.49, P = 1.7 × 10 − 4 ). We used individual fly lifespans to estimate a broad sense heritability for survival on H 2 O 2 food of H 2 = 44.8% ± 0.04 (mean ± se, Methods).
We then compared our measures of H 2 O 2 resistance in the DGRP with traits measured in the DGRP in previous studies. These studies have measured either survival or behavioral responses of the DGRP to two other oxidative stressors, paraquat and menadione [25,29]. We observed Pearson correlations of r = 0.34 (n = 179, P < 10 − 4 ) and r = 0.35 (n = 179, P < 10 − 4 ) between the mean survival time reported here and those measured on food supplemented with paraquat or menadione, respectively (Fig. S1 [25]). We found no correlation between our measures of H 2 O 2 survival and two behavioral traits, the startle response and climbing, measured by Jordan et al. [29] following chronic (13-16 day) exposure to menadione food (data not shown). We also found that the survival times of the DGRP on H 2 O 2 food correlate highly with survival measured under starvation in two different labs (Fig. S1 [26,37]) as well as data from our own lab (Fig. 1a). It is notable that the correlation between H 2 O 2 resistance and starvation was greater than the correlation between H 2 O 2 resistance and any other trait published for the DGRP, including survival on food containing paraquat or menadione (Fig. S1). These results suggest that H 2 O 2 food may affect feeding or nutrient assimilation in Drosophila.
To test the possibility that H 2 O 2 affects feeding behavior, we measured the amount of food consumed by flies in the CAFE and dye incorporation assays described in the Methods section [38,39]. Peroxide reduced feeding in each of three different genotypes tested, including strains with both relatively long and short survival time on H 2 O 2 food (Fig. 1b and c). During the 24 h feeding period in the CAFE assay, flies exposed to liquid H 2 O 2 food consumed no more than the volume lost due to evaporation in chambers without flies, suggesting that the flies consumed very little H 2 O 2 food over 24 h (Fig. 1c). This finding suggests that mortality in flies exposed to H 2 O 2 is due at least in part to starvation, in addition to any oxidative stress caused by H 2 O 2 exposure.

Genetic associations with peroxide resistance
The DGRP is both highly inbred and sequenced to high confidence for the majority of SNPs and small indels [22], therefore we sought to estimate the extent of genetic variation captured by the characterized genetic variation in the DGRP. We used restricted maximum likelihood to estimate the proportion of variance in phenotype that could be explained by the genomic relationship between the lines used in our study, the socalled SNP heritability (ĥ 2 SNP ) at 34.2% (Methods). Thus, a substantial amount of the heritable variation in H 2 O 2 resistance could potentially be explained by the characterized genetic variation in the DGRP. We then sought to identify individual genetic variants and genes that might associate with H 2 O 2 resistance.
We used a linear regression model in PLINK [40] to test for associations of H 2 O 2 resistance with approximately 1.9 million SNPs with a MAF ≥ 5% and <30% missing genotypes, while accounting for population structure and a significant effect of the major inversion In(2 L)t (Methods). We used the q value approach [41] to control the false discovery rate (FDR) and at 20% FDR, 14 variants (all were SNPs) were associated with resistance ( Fig. 2a, Table S1). Pairwise linkage disequilibrium (LD) among the 14 SNPs indicates that they associate with H 2 O 2 resistance as seven loci, or groups of SNPs in LD (r 2 > 0.5, Fig. 2b). With only 179 genotypes, we lack the power to analyze SNP-SNP interactions among these loci. However, our data indicate that H 2 O 2 resistance is polygenic within the DGRP.
To investigate these genetic associations further, we performed gene ontology (GO) analysis, looking for biological processes and signaling pathways that are overrepresented among the markers associated with survival on H 2 O 2 food. To identify gene-level associations, many studies use the minimum P-value of all variants in a gene (P min ). One might expect bias in this approach, such that genes with more variants are more likely to have a smaller (more significant) P min by chance alone. Indeed, we found that -log 10 (P min ) was positively associated with the number of variants per gene (Fig. S2), potentially biasing gene-trait associations in favor of genes with more variants [42]. To test for this bias, we compared the top 200 genes ranked by P min with the top 200 genes from ten GWAS of randomly permuted phenotypes. Out of 15,322 gene models in the Fly Base release 5.49, the null expectation for such intersection would be 2.6 genes. In contrast, we found an average of 11.6 ± 2.0 genes (mean ± se) in common across the permutations (χ 2 = 28.3, df = 3, P < 1.1 × 10 − 7 ), consistent with a bias caused by SNP density. To correct for this bias, we used a permutation approach to derive gene-level P-values, P gene , while also accounting for population structure (Methods). Unlike P min , P gene did not associate with the number of variants per gene, and it reduced the number of false-positives when compared to top genes from GWAS of randomized phenotypes (Fig. S2, χ 2 = 2.2 × 10 − 26 , df = 3, P = 1). Thus, P gene increases the accuracy of The correlations between mean lifespan of 31 DGRP lines on food containing either 2% glucose and 2% hydrogen peroxide in replicate trials (peroxide 1 and peroxide 2) or no glucose (starvation). Below the diagonal are plots of trait values (mean survival (h)). Least-squares linear regression lines are shown in red. Above the diagonal are Spearman correlation coefficients for each pair of traits. b and c Feeding assays, the mean (+SD) absorbance at 630 nm of extracts from three replicate vials of seven to eleven 1 to 5 day old mated females from the indicated DGRP line (b). Flies were exposed to peroxide (open bars) or control food (colored bars) for two hours of feeding prior to dye extraction. c The mean (+SD) volume of liquid food consumed by ten replicates of ten 2 to 4 day-old mated females was measured using the CAFE assay. CAFE food contained 2% glucose and, either water (colored bars) or 2% peroxide (open bars). Additional apparatus were set up without flies (no flies) to measure volume loss due to evaporation. Asterisks indicate p < 0.05 (Welch's t-test) gene-trait associations. Only one gene, ppk14, was significant below an FDR of 0.05, so we chose to look for biological processes or pathways that might be enriched among the 29 genes at FDR ≤ 0.5 (Table 1). This gene ontology (GO) enrichment analysis identified several biological processes and two biological pathways (FDR ≤ 0.05, Table 2). Most of the enriched processes are nested in hierarchical categories and thus are not independent. Also, the enrichment of 7 of the 9 biological processes, and the endothelin signaling pathway, is due entirely to three genes occurring in a small gene cluster, each encoding adenylyl cyclase (ACXA, ACXB and ACXE, Table 2). Adenylyl cyclase is involved in several signaling pathways, including G-protein coupled and calcium-based signaling. Separately, the platelet-derived growth factor (PDGF) signaling pathway is enriched (FDR = 0.02) due to three other genes in our dataset, Rab2, Ets21C and the c-Myc-binding protein homolog CG17202.
We validated these gene-trait associations by using RNAi to manipulate the expression of six of the 29 candidate genes, nAChRbeta3, Ets21C, ush, Nha1, Jon25Bi and Marcal1 (P gene < 0.0003 in all cases), and testing the effect of a mutation in a seventh candidate NPF (P gene < 0.0006). Several of the candidates reside near a cluster of six trait-associated SNPs within a 17.5 kb interval on chromosome 2 (Fig. 2a). This interval spans several genes, including u-shaped (ush, P gene = 7.1 × 10 − 5 ), which contains an intronic C/T SNP associated with H 2 O 2 resistance (P = 2.49 × 10 − 8 ) and several other SNPs in LD with this SNP (Fig. 2b). None of these SNPs are predicted to alter amino acid sequence or splicing, but instead may affect ush expression. Ush has roles in development and growth, including as a negative regulator of PI3K activity within the IlS pathway in the fat body [43]. To validate the effects of ush on H 2 O 2 resistance, we used the RU486-inducible GAL4 GeneSwitch driver S106 to drive RNAi targeting ush in the adult fat body [44]. RU486 treatment of flies carrying both S106 and UAS-RNAi targeting ush resulted in shorter survival times on H 2 O 2 food than the same genotype without RU486 (Fig. 3). We saw the same result with two independent RNAi lines, each targeting a different portion of the ush mRNA (Fig. 3). We saw no effect of RU486 on the survival of F 1 flies carrying the driver and the empty control P-element in the same genetic background as the UAS-RNAi flies, nor an effect of RNAi on lifespan on food lacking H 2 O 2 (see Methods). Knocking down ush in the nervous system with the elav GeneSwitch driver did not affect H 2 O 2 resistance (data not shown), and ubiquitous knockdown of ush appears to be lethal as we failed to recover Act-GAL4/ush-RNAi flies in crosses of the ush-RNAi construct to the constitutive Act-GAL4 driver. We also detected a strong interaction between RU486 treatment and an Act GeneSwitch driver in our H 2 O 2 food assay and so were unable to test the effect of ubiquitous knockdown Inset shows a locus on chromosome 2 L that includes 6 resistance-associated SNPs in LD (r ≥ 0.5) and several candidate genes, including ush. Genes shown in blue were significant in genelevel analysis (FDR < 0.5). b A pairwise LD plot for all 14 resistance-associated SNPs (FDR < 0.2), groups of SNPs in LD are indicated with bars across the top, and a color scale for r is shown. c A Q-Q plot of -log 10 P-values for the variants tested in the GWAS showing little to no inflation after adjusting for population structure (Methods) of ush in adult flies (data not shown). Similar RNAi of five other candidates did not affect H 2 O 2 resistance (data not shown).
Another candidate, NPF (P gene = 5.8 × 10 − 4 ) encodes the ligand neuropeptide F (NPF), which controls feeding, ethanol sensitivity and other behaviors [45]. The npf SK1 deletion allele reduced H 2 O 2 resistance when compared to wildtype control flies and did not affect survival on control food ( Fig. 3 and data not shown). These data demonstrate that manipulation of candidate genes from our GWAS affects the H 2 O 2 resistance phenotype.

Metabolite profiles associated with peroxide resistance
To investigate the effect of H 2 O 2 on the fly metabolome, and the potential for the metabolome to explain genetic variation in resistance to H 2 O 2 , we measured untargeted LC-MS profiles in three biological replicates for each of eight resistant (mean survival time = 106.9 h, range = 90.4-119.7 h) and eight sensitive (mean survival time = 58.9 h, range = 53.2-70.0 h) lines, chosen such that the two groups did not differ substantially in fly size (Fig. S3). These lines were subjected to another survival assay and, 24 h after being exposed to H 2 O 2 or control food, samples of flies from each line and treatment were flash frozen for aqueous metabolite extraction, while survival measurements were conducted on the remaining flies. Global metabolite profiling measured 2722 and 2691 features that passed quality control thresholds from positive and negative ionization modes, respectively and were analyzed separately (Methods). We first explored these data using principal component analysis (PCA). Samples from resistant and sensitive genotypes clearly separate along principal component one of the negative mode (PC1 neg , Fig. 4) and along PC2 of the positive mode (PC2 pos , Fig. S4). Interestingly, H 2 O 2 -treatment samples appear distinct from untreated samples among the sensitive genotypes along both PC1 neg and PC2 pos , and this separation is not apparent among the resistant genotypes ( Fig. 4 and Fig. S4). Principal component analysis thus detected between-group variation in metabolite profiles from sensitive and resistant flies, and further suggested an effect of treatment on sensitive that is not seen in resistant flies.
The separation of resistant and sensitive lines by metabolome profile is striking (Fig. 5b). However, the genotypes chosen for metabolite profiling were among the extremes of resistance to H 2 O 2 . This design raises the possibility that resistant and sensitive lines have a confounding genetic relationship. To determine if genotype could also separate lines chosen deliberately with extreme phenotypes, we carried out PCA and hierarchical clustering of the genotype data on the same lines. We analyzed the first ten PCs of genotype, which together account for 69% of the variance of > 3.2 × 10 5 genetic variants in the 16 lines used for metabolomics. These PCs failed to clearly separate the resistant and sensitive flies (for example, PC1 vs. PC2 is shown in Fig. 5). Similarly, hierarchical clustering of these 16 lines using the  . Asterisks indicate significant effect of RU486, or of the NPF SK1 mutation (P < 0.05) from the log-rank test using the survival package in R. Results are representative of at least two independent experiments, and there was no mortality observed on control food during these experiments same genetic variants also failed to separate the two phenotypes (Fig. 5). For comparison, clustering of 2691 negative mode metabolite features separated genotypes into two clusters composed of resistant and sensitive lines (Fig. 5). Thus, the distinction between resistant and sensitive flies is more obvious given an unbiased sample of LC-MS features rather than an unbiased sample of their genotype.

Metabolic pathways associated with peroxide resistance
Principal component analysis suggested that the metabolome of sensitive and resistant lines differ systematically and moreover, that resistant and sensitive genotypes differ in their metabolic response to H 2 O 2 treatment (Fig. 4). To identify the specific metabolites whose abundance was affected by treatment in a trait-dependent way we ran a linear regression model, predicting metabolite level in response to the interaction between the trait and treatment. In the negative mode data we found 105 features that were significant for a trait by treatment interaction (FDR < 0.1). Figure S4 shows the clustering of the features with a significant treatment by trait interaction term. These data come from analysis of profiles from negative mode only, as no features from positive mode reached our FDR cutoff for the interaction term. The clustering of features by z-score across samples revealed a consistent pattern with a substantially different effect of H 2 O 2 on the metabolome of sensitive genotypes compared to resistant genotypes (Fig. S5). This pattern is similar to the apparent separation of samples across the latent variables revealed in PCA (Fig. 4). We used the mummichog software package [46] to identify metabolomic pathways enriched among the features associated with H 2 O 2 resistance, or among features with significant treatment by trait interaction effects. Mummichog matches m/z ratios and retention time data for the features to those predicted to occur if those features were enriched for a given metabolic pathway [46]. It is important to point out that we use mummichog as a tool to identify metabolic pathways that could associate with variation in H 2 O 2 resistance and not as a tool for metabolite annotation. This analysis identified several pathways, including carbohydrate metabolism, amino acids and their biosynthesis, and folate metabolism (see Table S2, Additional File 2). Many of these pathways share overlapping metabolites and some pathways are nested within other pathways. We chose a subset of the identified pathways for further analysis based on a several criteria. One criterion was the significance of each pathway across the different contrasts in the linear model (see Table S2, Additional File 2). Another criterion was the strength of the identification; we gave higher priority to those features that were uniquely assigned to a particular metabolite or pathway, rather than features that were ambiguously associated with more than one metabolite, or with a metabolite that is present in several pathways.

Glycogen metabolism
Glycogen metabolism was identified by mummichog among the features that showed significant trait by treatment interaction (Fig. S6). Glycogen is a storage form source of glucose in Drosophila and variation in glycogen metabolism across the lines measured here may result in sensitivity to H 2 O 2 food. Consistent with a role for glycogen in resistance to H 2 O 2 food, flies fed supplemental maltose prior to being exposed to H 2 O 2 showed increased resistance to H 2 O 2 food in a dose-dependent manner (Fig. 6). Maltose could increase survival by providing glucose, or perhaps by some other effect as a disaccharide. Supplementing fly diet with glucose but not the disaccharide lactose extends survival time similar to supplemental maltose, which is consistent with the former hypothesis (Fig. 6).
In addition to experiments with supplemental glycogen intermediates, we also used genetic manipulation to test for a role of glycogen metabolism in resistance to H 2 O 2 food. The fat body is a site of glycogen storage in Drosophila, and RNAi targeting glycogenin (CG44244), the gene encoding the protein core of glycogen, in the fat body increased the sensitivity of flies to H 2 O 2 food (Fig. 7). We found that knocking down glycogenin using the S32 fat-body driver reduced survivorship, while the S106 fat-body driver did not (Fig. 7). Glycogen is present in the Drosophila brain, and knocking down glycogenin with the neuron-specific elav GeneSwitch driver reduced survival on H 2 O 2 substantially (Fig. 7) [47]. These data suggest that glycogen metabolism differentiates sensitive verses resistant flies.

Folate metabolism
The folate pathway was also identified by mummichog among the features associated with H 2 O 2 resistance and treatment (see Table S2, Additional File 2). The folate pathway is central to the synthesis of several amino acids, nucleotides, secondary metabolites and substrates for secondary modifications (e.g. methylation). Mummichog detected features that are consistent with metabolites both in the folate pathway as well as in peripheral pathways (e.g., S-adenosylmethionine), suggesting that  . Asterisks indicate significant effect of RU486 (P < 0.05, log-rank test). Data are representative of one or two experiments, and the mortality of flies on control food (not shown) did not exceed 15% in any condition and did not depend on RU486 exposure the activity of the folate pathway differs between sensitive and resistant flies. To investigate the potential role of the folate pathway in H 2 O 2 resistance, we tested the effect of metabolic gene knockdown on survival. Knocking down either CG8665, which encodes 10-formyltetrahydrofolate dehydrogenase (FDH), or pugilist (CG4067), which encodes tetrahydrofolate dehydrogenase (THFDH), in the fat body or in neurons reduced the survival of flies on H 2 O 2 food (Fig. 7). In parallel experiments, RU486 pretreatment failed to affect survival in flies carrying the fat body S32 driver. These data suggest that folate metabolism in the abdominal fat body and neurons is important for survival on H 2 O 2 food.
Supplemental folic acid was shown to reduce the levels of oxidative damage associated with knockdown of parkin in Drosophila neurons [48]. We also tested whether supplemental folic acid would affect survival on H 2 O 2 food and failed to see a significant effect (data not shown).

Dissecting variation for complex traits
Next-generation sequencing technology has provided geneticists with unprecedented power to identify single nucleotide variants associated with variation in complex traits in natural populations. But even with extremely large sample size (e.g., [20,49,50]), the percent of variance explained by SNPs in most studies remains small [51][52][53], and current models suggest that thousands of SNPs can contribute to any one trait [52,54]. Our modestly powered GWAS failed to detect individual SNPs that explain a large amount of the variation in H 2 O 2 resistance. However, our mixed model analysis estimates that approximately 34% of the variation in peroxide resistance could be explained by the additive effects of genetic variants within the DGRP. Similar to GWAS studies of other traits in this population, our results suggest that peroxide resistance is polygenic in the DGRP [28,55,56].
The results presented here also point to the considerable potential of metabolic profiles to distinguish genetically determined phenotypic variation (Fig. 5), and moreover, to identify novel, causal molecular pathways associated with that variation. In light of our findings, we propose a model here whereby a large number of interacting genetic loci [57] converge through a more limited number of downstream metabolic pathways, which in turn make up the functional and structural building blocks of complex traits (Fig. 8).

Metabolomic analysis
This study illustrates that untargeted metabolomic profiles, even those that include unknown chemical identities, give us tremendous power to 1) explain complex phenotypic variation; 2) identify novel pathways associated with this variation; and 3) in this particular study, suggest a novel hypothesis that resistance to stress might be caused by resistance of the metabolome to environmental perturbation.

The Metabolome as a biomarker
First, we find that the untargeted metabolome when compared to genotype associates more closely with phenotype. While a PCA of the metabolome clearly separates sensitive from resistant flies, a similar analysis of allelic variants among genotypes failed to clearly separate genotypes based on their survival on H 2 O 2 food (Fig. 5). This contrast highlights the 'proximity' of the metabolome to trait on the genotype-phenotype map [4,58].
Numerous other studies have shown metabolomic responses to diet, age, and temperature [32,[59][60][61], as well as body mass and body composition independent of diet [31], in diverse genotypes. However, these studies were not designed to test the relative association of genotype versus metabolome with phenotypic response. Moreover, studies that describe the effect of stress or environment on the metabolome often include only a single genotype and this may fail to resolve the systems-level association between genotype, environment, metabolome and phenotype [62][63][64][65]. In light of our results, future studies would benefit from a clearer characterization of the power of the metabolome relative to the genome to distinguish biologically relevant phenotypic variation.

Stress resistance pathways identified by metabolomics
Second, our metabolomic profiling suggests several possible mechanisms that underlie the phenotypic variation observed here. As we emphasize above in presenting the results, and discuss further here, in Drosophila studies, variation in the ability to survive oxidative stress could be confounded with starvation resistance. This is important to keep in mind as we discuss possible mechanisms that underlie the phenotypic variation seen here. Nonetheless, the ability of the metabolome to distinguish resistant and sensitive phenotypes is notable.
Several studies have used metabolomics to study the effects of oxidative stressors in Drosophila in a single genetic background. For example, paraquat treatment alters branched-chain amino acid, starch/sucrose, and fatty acid metabolism in the Drosophila brain [65]. We also detect significant effects of H 2 O 2 food on the first two of these three pathways (see Table S2, Additional File 2). The lack of evidence for an effect on the third pathway, fatty acid metabolism, is perhaps not surprising, given that our analysis was limited to aqueous metabolites. There are two important caveats to this earlier work as well as to the present study. First is the ambiguous nature of the untargeted metabolite data. For many features in the global metabolome profiles, we know mass/charge ratios but not chemical structure. Second, given the relationship between resistance to H 2 O 2 and starvation ( Fig. 1), we expect that the metabolites associated with exposure to H 2 O 2 or paraquat might relate to nutrient intake and/or storage, rather than oxidative stress alone. This confounding influence of H 2 O 2 on feeding could explain the results of a variety of studies that use survival on food supplemented with H 2 O 2 as a measure of oxidative stress resistance.
One way to disentangle the effects of oxidative stress from a secondary effect of a stressor administered in the diet on nutrient update is to genetically manipulate endogenous levels of reactive oxygen species. Two recent studies compared the metabolome of flies with null mutations in superoxide dismutase (sod − ) to that of flies either chronically or acutely exposed to paraquatsupplemented diet [62,63]. Each of these conditions had a distinct effect on the metabolome, with the sod − mutant metabolic profile clearly separating from the other conditions in unsupervised clustering [63]. Our analysis suggests that at least some of the differences between the sod − and paraquat-treated metabolome could be due to altered feeding or starvation in paraquat-treated flies.
Keeping these caveats in mind, we find strong evidence for two pathways associated with H 2 O 2 resistanceglycogen metabolism and folate metabolism. Indeed, the apparent decrease in glycogen content in flies sensitive to H 2 O 2 compared to control food suggests that sensitive flies are exhausting their glycogen pool in response to stress, which is consistent with previous studies of flies on food with paraquat or under starvation [36,66,67]. Maltodextrins are intermediates in the glycogen pathway, which is used to store and retrieve glucose for the glycolysis and pentose phosphate pathways (Fig. S5). We found that supplemental maltose or glucose enhances survival on H 2 O 2 food for all genotypes tested in our study and that RNAi targeting glycogenin either in the fat body or in neurons reduces survival on H 2 O 2 food (Figs. 6 and 7).
Interestingly, knockdown of glycogenin affected survival, but only with the S32 driver, which is expressed primarily in the head fat body, and not the S106 driver, which is expressed in the abdominal fat body. This suggests that glycogenin functions in the head fat body to mediate resistance to H 2 O 2 food [68,69]. However, we cannot rule out the effect of other differences between the knockdown of glycogenin by S32 compared to the S106 driver.
Mummichog also detected enrichment of the folate pathway (see Table S2, Additional File 2). Several lines of evidence from this work and previous studies indicate a role for the folate pathway in stress resistance in Drosophila [25,48,[70][71][72][73][74][75]. We show that knocking down either the folate pathway genes pugilist or CG8665 in either the fat body or neurons reduces survival on H 2 O 2 food while not affecting survival on food without H 2 O 2 ( Fig. 7 and data not shown). Several studies find that transcripts encoding enzymes of the folate pathway are up-regulated following either administration of paraquat or H 2 O 2 in food, over-expression of manganese superoxide dismutase, or in mutants with mitochondrial Fig. 8 The Metabolome as a Bridge in the Genotype to Phenotype Map. Phenotypic variation, represented by the bottom plane, is influenced by the action of a large number of genes and genetic variants (upper plane) on a smaller number of metabolites or metabolic pathways. Direct influences on phenotype are shown as dark lines whereas indirect connections are shown with grey lines. Genetic variation can influence the levels of individual metabolites or the relationship between metabolites and, as a result, affect metabolic pathways that share direct connections to phenotype; or, genetic variation can affect phenotype more directly (not depicted). Two or more genes might act epistatically to influence the level of a single metabolite, and a single gene might act pleiotropically on multiple metabolites dysfunction [72][73][74][75][76][77]. Additionally, two missense variants in pugilist associate with paraquat resistance in the DGRP [25]. Functional studies show that overexpression of Drosophila nmdmc, which encodes methylenetetrahydrofolate dehydrogenase, enhances resistance to paraquat [71], and that folic acid supplementation reduces oxidative damage to lipids and endogenous H 2 O 2 levels associated with knockdown of parkin in Drosophila neurons [48]. While our metabolomic data and previous genetic studies [70] point to a role for folate metabolism, in our own GWAS analysis, neither pugilist nor CG8665 were clearly associated with survival on H 2 O 2 (P min > 0.002 in both cases).

Adaptive shifts versus robustness of the metabolome under stress
Third, it appears that the metabolome responds to the presence of H 2 O 2 more strongly in sensitive than in resistant lines (Figs. 4 and S4). One might have expected the opposite pattern, whereby greater resistance is associated with an adaptive response of the metabolome to an environmental stressor (e.g., [78][79][80]). Instead, it appears that resistance in this study is associated with the ability to maintain the metabolome in its current state. This observation is made possible by having multiple sets of stress-resistant and stress-sensitive genotypes. We interpret these results as suggesting not only that resistance to H 2 O 2 food is explained by the metabolome, but also that the metabolome is more robust in resistant genotypes, being less likely to change when faced with a stressor. We do not propose that metabolic robustness is a universal or causal feature of stress resistance, though this is an interesting hypothesis to test, but our data suggest that metabolic robustness and/or resilience might contribute to phenotypic variation in nature.
This observation, combined with the fact that the metabolome associates with phenotype while the genotype does not, suggests the hypothesis that there are genetically diverse ways to achieve resistance to H 2 O 2 , but that these diverse genetic paths converge at a common metabolome associated with resistance. In future studies it would be worth asking if, across a broad range of phenotypes, there is a correlation between stress resistance and the 'resistance' of the metabolome to stress-induced alteration, or vice versa for traits that require metabolic adaptation.

Genomic analysis
We studied survival on H 2 O 2 food in the DGRP, a population representing a sample of natural genetic variation. We estimate broad-sense heritability for mean lifespan under H 2 O 2 at H 2 = 44.8% based on the within-genotype variance in phenotype compared to total phenotypic variation. There are significant caveats to estimating heritability in the DGRP, including the relatively small number of lines and the potential to overestimate the degree of heritability due to the low among-line environmental variance associated with repeated measures of biological replicates within inbred genotypes. With these caveats in mind, we estimateĥ 2 SNP of 34.2%. This estimate is slightly higher than the narrow sense heritability h 2 estimated for starvation resistance in other studies in the DGRP using similar methods, although estimates of h 2 in the DGRP vary widely [56,81]. Given that the DGRP is highly inbred, we have not empirically evaluated the extent of additive gene action in our study. However, our results indicate that genetic variation contributes significantly to H 2 O 2 resistance in the DGRP.
Mapping variants associated with mean survival revealed 14 SNPs that were significant and together these variants define 7 loci (Fig. 2). We limited our analysis to variants passing MAF and missing genotype filters among 179 DGRP lines. While the DGRP incorporates substantial allelic diversity from the wild, it is a population of inbred lines with very low levels of heterozygosity compared to its parent wild population and that has also been purged of deleterious alleles of strong effect [55]. Along with the missing genotypes there are also uncharacterized structural variants within this population [22]. For these reasons, we do not expect this study to identify all variants in the DGRP associated with survival under H 2 O 2 stress, nor their mode of gene action. Similar to other studies with the DGRP, we failed to find common alleles with large and highly significant effects, suggesting that variation in survival on H 2 O 2 food is influenced by many loci of small effect [55,56]. Importantly, in the relatively small population used for metabolomics, neither clustering nor multivariate analysis of genotypes associated with the discrete resistance trait whereas similar analyses of metabolome data did distinguish resistant and sensitive genotypes (Fig. 5).
Although we hoped to pinpoint specific genes that influence H 2 O 2 resistance, we face the statistical challenge that among the DGRP lines we studied, the number of variants per gene ranged from 1 to 4490 with a mean of 237. To overcome the increased risk of false positives in genes that contain a large number of polymorphisms, we used a permutation approach to measure association between genes and phenotype. Our approach, like other methods, comes with caveats, one being the imperfect annotations of variants and genes. We have limited this analysis to only those variants associated with the Fly-Base gene annotation, including 1 kb upstream and 1 kb downstream of the primary transcript [22]. Intergenic variants might affect the expression of trait-associated genes. However, we did not attempt to account for those effects in this study. Additionally, variants that are associated with candidate genes can instead exert their effect on phenotype by modifying the expression of other local genes, rather than the candidate gene, and this would result in misattributing the significance to the gene containing the variant rather than the real trait-associated gene [82].

Stress resistance pathways identified by genetic association
Our genomic analysis led to several interesting and potentially related genes associated with H 2 O 2 resistance, including ush, NPF and the pickpocket paralogs ppk7 and ppk14. Each of these genes is implicated in several processes and though we do not rule out causal associations between these processes and resistance to H 2 O 2 food, based on published studies, we argue that genetic variation in these genes influences feeding and/or metabolism to explain their effect on H 2 O 2 resistance (Fig. S7).
As we show, flies on H 2 O 2 food substantially reduce food intake, and their lifespans appear to be a function of starvation (Fig. 1). Interestingly, feeding behavior is controlled by neuronal signaling involving several of these candidate genes. Larvae require NPF to signal the intake of noxious food under starvation conditions [83], and NPF-expressing neurons in adults couple hunger to memory performance [84]. Interestingly, IlS in neurons expressing the NPF receptor represses larval feeding, suggesting that NPF controls feeding through neuronal IlS [83]. We postulate that the candidate ush influences survival on H 2 O 2 food as a peripheral regulator of IlS [43]. This appears to be a function directly attributed to the USH protein, as its human homolog FOG has been shown to directly bind and inhibit the PI3K complex [43]. Inactivation of PI3K leads to depletion of nutrient stores in the fat body and its constitutive activation reduces both nutrient stores and survival under starvation [85]. PI3K is also a negative regulator of FOXO in the IlS pathway, which has well characterized roles in the response to starvation and to oxidative stressors in food [68,[86][87][88][89]. Interestingly, ablation of insulin-like peptide-producing cells in Drosophila increases survival both on food containing paraquat and under starvation, alters whole body levels of glycogen, and leads to misregulation of metabolic genes in the glycogen pathway, which could link the effect of NPF and ush polymorphism to the variation we see in glycogen metabolism [90,91].
The effects of two other candidates, ppk7 and ppk14, might also be linked to feeding and metabolism through their potential role in nutrient signaling. ppk7 and ppk14 are members of the pickpocket gene family that are expressed in neurons and signal taste cues, modulate feeding, and may influence energy metabolism [92][93][94][95]. Interestingly, Drosophila ppk28 was recently shown to interact with glucagon-like hormone (AKH) signaling, a pathway involved in regulating glycogen metabolism in flies [91,96]. Together, the genetic analysis of H 2 O 2 resistance has revealed pathways whose role in survival could be explored in future studies.

Starvation as a confounding factor in stress assays
Our study used H 2 O 2 -supplemented food as the stressor. Oxidative stressors such as H 2 O 2 , paraquat, or menadione are often administered to Drosophila by supplementing the diet, and each of these treatments dramatically reduces survival in a dose-dependent manner [25,87]. After screening the DGRP for survival on H 2 O 2 food, we noticed that these survival times correlated closely with survival times under both paraquat and menadione exposure, and even more so with survival times measured under starvation [26].
Feeding is essential to the survival of adult Drosophila, and feeding behavior is influenced by a variety of cues, including food acidity or the presence of bitter compounds, hypoxia, and the nutrient content of food [97][98][99][100]. While some stressors affect the preference for food of particular composition [98,101], others may alter feeding behavior by affecting satiety [97]. We show that flies consume very little food containing 2% H 2 O 2 . This effect is not limited to H 2 O 2 , as many supplements have been found to reduce feeding in flies [97], and this may extend to other oxidative stressors including paraquat and menadione. While several studies have detected reduced feeding in response to paraquat-containing food [102,103], others detected no difference [67]. Contrary to the latter study, our data suggest that feeding is significantly reduced in the three genotypes tested in assays that measured feeding over either 2 h or 24 h (Fig. 1). Several differences in the experimental setting may explain the discrepancy between this study and Riahi et al. (2019) [67].
It is possible that the effect we see relates to oxidative stress and not an aversion to food supplementation, as deviation from normoxia alters feeding behavior of Drosophila larvae in a manner that appears to rely on H 2 O 2sensitive neurons [104]. Also, the UV light-avoidance of egg-laying females appears to signal through H 2 O 2 -sensitive taste receptors in Drosophila [105], and H 2 O 2 also inhibits feeding in Caenorhabditis elegans through taste receptors [106]. Additionally, two recent studies found contradictory roles for the histone methyltransferase G9a in survival in response to oxidative stressors and starvation, indicating that environment and genetic background may affect stress-response in Drosophila, but also suggesting that starvation and oxidative stress resistance may share underlying biological pathways [36,67]. These studies suggest that susceptibility to oxidative stress and starvation are partially separable. However, they do not rule out a main effect of H 2 O 2 on survival due to starvation. The relationship between lifespan under starvation and survival on H 2 O 2 food has implications for studies that draw conclusions about stress resistance in response to agents administered in food. The effects of altered feeding patterns or nutrient deprivation should be accounted for when analyzing the effects of stressors or drugs administered in the diet.

Conclusions
Genetic variation in a complex trait converges on the metabolome The sample of genotypes in this study show a consistent metabolic signature associated with their phenotype. Thus, a potentially wide degree of genotype space may converge on a smaller number of metabolic pathways to shape phenotype.

Metabolome robustness associates with stress resistance
Contrary to the metabolic change that might be expected in animals resisting stress, we find that the metabolome of resistant animals appears robust to stress treatment. This suggests that maintaining metabolism in the presence of certain stressors is a means of survival.
Glycogen and folate metabolism and several genes involved in nutrient signaling mediate resistance to peroxide food Genetic and metabolomic analysis of peroxide resistance revealed roles for glycogen and folate metabolism and genes with known roles in nutrient signaling. Future studies to understand this network may reveal novel mechanisms of stress resistance.

Starvation explains the lifespan response to peroxide food
Drosophila are sensitive to diet, including additives and food-borne treatments. We show that the response to H 2 O 2 in food can be explained by starvation. These effects may dramatically confound assays that examine responses to treatments delivered by supplementing the Drosophila diet.

Media
Standard food was made by cooking 12 g Drosophila agar (type II, Genesee Scientific, El Cajon, CA), 25 g brewers yeast (MP Biomedicals, Solon, OH), 55 g glucose monohydrate (MP Biomedicals), 30 g sucrose, 60 g corn meal, 3 g methylparaben (Genesee Scientific). 12 g 100% ethanol (Decon Labs, King of Prussia, PA), and 3 g propionic acid (Fisher Scientific, Pittsburg, PA) per liter of water. A small amount of dry active yeast was sprinkled onto standard food prior to use.
Peroxide food was made in one of two ways, for 2% food, agar was melted into 2% glucose monohydrate and 0.3% propionic acid and, after the food had cooled to less than 60°C, 30% H 2 O 2 (Fisher Scientific) was added to reach 2%, or the same volume of water was added for the control food. For the 3% H 2 O 2 food, the recipe was the same with the exception that 30% H 2 O 2 was added to reach 3% H 2 O 2 . Approximately 5 mL of food was dispensed into 25 mm wide × 95 mm tall polystyrene vials.
Starvation food was made by melting 2% agar into 0.3% propionic acid and dispensing into vials.
Food supplemented with carbohydrates was made by adding 2% of either D-(+)-maltose, β-lactose (both from Sigma, St. Louis, MO), or additional glucose to the 2% glucose control food.
RU486 food to induce the GAL4 GeneSwitch system was made by overlaying~5 mL standard food with either 50uL of 25 mg mL − 1 RU486 (mifepristone, Cayman Chemical Company, Ann Arbor, MI) dissolved in 100% ethanol or the same volume of 100% ethanol alone for the -RU control food. Ethanol was allowed to evaporate overnight at 22 to 24°C prior to using the food.
Genetic manipulation F 1 GAL4/UAS flies were collected over four days (day 0 to 3). These flies were allowed to mate for 24 h, at which time they were anesthetized and sexed, and females (ten per vial) were then allowed to feed for two days on RU486 or -RU food. After 48 h on RU486 or -RU food, flies were transferred without anesthesia to H 2 O 2 or control food to measure survival. Negative genetic controls included F 1 GAL4/attP flies which were crosses of the GAL4 driver to either the attP2 or attP40 lines from the Transgenic RNAi Project collection, where attP is the empty P-element docking site for the UAS transgenes (http://fgr.hms.harvard.edu/fly-in-vivo-rnai). Negative genetic controls were raised, induced and assayed in parallel with experimental flies.

Stress survival assays
To measure the variation in resistance to oxidative stress across a lab population, we measured the survival of mated females from 179 DGRP lines in a multi-block design on H 2 O 2 food [26]. For each block, flies were raised under low-density conditions by allowing~50 flies to lay eggs for one day on standard food in bottles. Flies for the assay were then collected over two or three days and then allowed to mate for 24 h on standard food. Four vials of flies on H 2 O 2 food and one or two vials of flies on control food without H 2 O 2 were included for each genotype in each trial. Each vial contained 10 mated females. Control vials without H 2 O 2 were included to confirm that mortality was due to H 2 O 2 . In knockdown experiments, we included 5 to 8 H 2 O 2 vials and 5 to 8 control vials to ensure that any effect of gene knockdown on survival in the absence of H 2 O 2 could be measured. To measure line weights, 24 h after beginning a lifespan assay, flies from an extra control food vial were frozen at − 80°C and were later collectively weighed on a microbalance (XS105, Mettler Toledo, Columbus, OH). In all blocks, dead flies in each vial were recorded two to four times per day using D-life software until all H 2 O 2 -treated flies had died [23]. For assays involving the Drosophila Activity Monitor System (DAMS, TriKinetics Inc., Waltham, MA), the activity of 38 to 48 individual flies per genotype was recorded simultaneously every minute over the experiment.
All calculations were performed in R [108]. Mean lifespan was estimated from H 2 O 2 assays using the restricted mean (default settings) in the Kaplan-Meier model with the survival package [109]. We ran 17 blocks with a mean of 14.5 lines (range = 4 to 35 lines) per block. We used 2% H 2 O 2 food for the first 10 blocks and 3% H 2 O 2 was used for the last seven blocks. The switch between 2 and 3% H 2 O 2 was made accidentally and was realized after the conclusion of the study. We used the following mixed model in the R package lme4 to test for an effect of these two food treatments on the log of mean lifespan: : log e li f espan ¼ f ood þ weight þ ð1jblockÞ þ ð1jlineÞ þ ε where food (2% or 3%) and weight were fixed effects, and block and genotype were both random effects along with the error term ϵ. The significance of random effects was assessed by the likelihood ratio test, and of fixed effects by ANOVA. We found no difference between 2 and 3% H 2 O 2 doses on lifespan (β = − 9.1 × 10 − 4 , P = 0.988).
To compare lifespans on H 2 O 2 food to lifespans during starvation, 2 to 5 replicates of twenty 3-to-5 day-old mated females were assayed using D-life on agar food either with or without 2% glucose (see Media). The lifespan of each line on H 2 O 2 food was measured twice in separate trials for this comparison, while the lifespan under starvation was measured once. To measure the effect of supplemental carbohydrates on lifespan ten replicates of ten 1-to-3 day-old mated females per genotype were allowed to feed on supplemented food for four days and then transferred without anesthesia to 2% H 2 O 2 or control food to assay survival.

Feeding assays
To measure feeding rate, we used both dye incorporation and CAFE assays. For both assays, flies were allowed to mate for 24 h and then separated sexes over light CO 2 anesthesia and transferred to agar-only food for 24 h of starvation. For dye incorporation, after starvation, flies were immediately transferred without anesthesia into vials that contained either H 2 O 2 or control food with 2.5% FD&C Blue Dye #1 (Spectrum Chemicals, Gardena, CA). After 2 h on dye-containing food, flies were flash frozen in liquid N 2 , homogenized in water, centrifuged at 16,000rcf for 1 min, and the absorbance of the supernatant was measured at 630 nm. The absorbance of each sample was normalized by dividing by the number of flies in the sample (n = 7 to 11 flies per sample).
For the CAFE assay, ten replicates of 10 mated females were starved for 24 h and transferred without anesthesia to assay chambers. Assay chambers were 15 mL conical bottom polystyrene tubes (Corning Inc., Corning, NY) containing water under a foam partition to maintain humidity but not allow flies to drink, and fitted with a 0.75 mm ID glass capillary (World Precision Instruments, Sarasota, FL) which had been filled with 2% glucose, 0.3% propionic acid supplemented with either 2% H 2 O 2 or water. Flies were housed in the assay chambers at 25°C on a 12/12 h D/L cycle in an incubator at 60-70% RH for 24 h before the volume of food consumed was assayed by measuring the difference in height of the top of the liquid food in the capillary and multiplying by π·0.375mm 2 .

Genetic analysis
We estimated broad-sense heritability of fly lifespan within each of the seven blocks in which at least 13 randomly chosen DGRP lines were included, treating each block as an independent measure of heritability. For each block, heritability was estimated by: σ L 2 / (σ L 2 + σ E 2 ) in an ANOVA, where σ L 2 is the among-line variance in the weight-residuals of lifespan and σ E 2 is the average within-line variance [110]. Mean heritability across blocks and its standard error were calculated from these seven estimates. We estimated SNP heritability (ĥ 2 SNP ) using genetic variance ( σ g 2 ) and the residual variance ( σ e 2 ) estimated by restricted maximum likelihood in the NAM package [54,111] using the model log e li f espan ¼ X b þ K þ ε where the lifespans were from the 9988 individuals across the study, X b is the design matrix of fixed effects, which include block, mean fly weight, Wolbachia status and the genotype at the four segregating inversions whose MAF was at least 3%: In(2 L)t, In(2R)NS, In(3R)P and In(3R)Mo. K is a genetic relationship matrix made with 712,878 LD-pruned variants (r 2 < 0.5), MAF ≥5%, genotype call rate ≥70%, in PLINK according to [112] and is modeled as a random effect, and ε is the error term.
For GWAS, we used residuals from linear regression of lifespan versus weight as our measure of H 2 O 2 resistance. To correct for block effects, weight residuals of lifespan were centered (mean = 0) and scaled to unit variance (SD = 1) by block. For lines measured in more than one block, we calculated the average blockcentered lifespans across blocks. We used stepwise regression in the MASS package to identify significant covariates among the chromosomal inversions and Wolbachia status [22]. In(2 L) t was the only significant covariate (P = 0.0277, ANOVA) and residuals from the linear regression of H 2 O 2 resistance on In(2 L) t genotype were used in a linear model of SNP-phenotype associations in PLINK [40]. Approximately 1.93 million SNPs with a MAF ≥0.05 and <30% missing genotypes were tested for association with H 2 O 2 resistance from 179 DGRP lines. To account for population structure, we used the Tracy-Windom test in the AssocTests package to evaluate eigenvalues from 20 PCs of genotype, and retained the first four PCs as covariates in the model (α = 0.05, [113]). Genome-wide significance was determined by controlling for FDR at 0.2 using the q value method [41].
Gene-level associations with H 2 O 2 resistance were derived using two rounds of permutation to reduce computational burden and to estimate more precise empirical P values. The first round tested 4810 genes that had at least one variant associated with phenotype with P ≤ 0.01, thus avoiding computation of gene-level statistics for genes with a very low chance of being significant. We refer to SNPs within 1 kb upstream and 1 kb downstream of the gene model in FB release 5.49 as associated with that gene [22,114]. In the initial analysis, 10,000 permutations of phenotype were performed, and the association between phenotype and each SNP in a gene was tested using the same linear model employed for SNP-phenotype associations described above. An initial empirical one-tailed P value for each gene (P gene ) was calculated by comparing the maximum test statistic (T max ) among the SNPs in each gene from the real GWAS to the T max for each of 10,000 permutations. A second round of selection was then run, this time choosing only genes with P gene ≤ 0.01 from the first round. The 192 genes with initial P gene ≤ 0.01 were then subjected to 1 million more permutations and the resulting P gene values were used as our measure of gene-trait association. We use q values to estimate the FDR for each P gene . Over-representation by biological process and pathways was tested using Fisher's exact test in PAN-THER (version 13.1) and the GO-slim subset of biological processes and PANTHER pathways among the 13,767 gene models in Drosophila melanogaster [115].

Metabolomic analysis
Eight of the resistant and eight of the sensitive DGRP lines were selected based on their lifespans and line weights to reduce the effect of fly size on resistance. Lifespan for these lines on H 2 O 2 food was again measured in a single block, and 24 h after exposure to the H 2 O 2 or control food, 3 replicates of 5 flies each were collected, flash frozen in liquid nitrogen, and then stored at − 80°C. Each Drosophila sample was weighed and then homogenized in 200 μL water with PBS in a microfuge tube immersed in an ice bath. Methanol (800 μL) was then added, followed by vortexing for 2 min and incubation at − 20°C for 30 min to precipitate proteins. Samples were sonicated in an ice bath for 10 min and then centrifuged at 17,000 rcf for 5 min at 4°C. From each tube, 900 μL supernatant was transferred to a new microfuge tube for drying under vacuum at 30°C (~3 h). The completely dried samples were reconstituted in 100 μL 40% water/60% HPLC-grade acetonitrile (ACN, Fisher Scientific) for liquid chromatography-mass spectroscopy (LC-MS) analysis. A pooled quality control (QC) sample was made by combining~5 μL aliquots from each reconstituted sample. The QC was analyzed once for every ten study samples to serve as a technical replicate throughout the data set to assess process reproducibility and allow for data normalization to account for any instrument drift.
LC-MS analysis was performed using an LC-QTOF-MS system (Agilent Technologies, Santa Clara, CA) consisting of an Agilent 1200 SL liquid chromatography system coupled online with an Agilent 6520 time-of-flight mass spectrometer. A 5 μL aliquot of reconstituted sample was injected onto a 2.1 × 150 mm Waters BEH-Amide 2.5 μm particle column at 35°C. The metabolites were gradient-eluted at 0.3 mL/min using mobile phase