Network analysis reveals differential metabolic functionality in antibiotic-resistant Pseudomonas aeruginosa

Metabolic adaptations accompanying the development of antibiotic resistance in bacteria remain poorly understood. To interrogate this relationship, we profiled the growth of lab-evolved antibiotic-resistant lineages of the opportunistic pathogen Pseudomonas aeruginosa across 190 unique carbon sources. We semi-automatically calculated growth dynamics (maximum growth density, growth rate, and time to mid-exponential phase) of over 2,800 growth curves. These data revealed that the evolution of antibiotic resistance resulted in systems-level changes to growth dynamics and metabolic phenotype. Drug-resistant lineages predominantly displayed decreased growth relative to the ancestral lineage; however, resistant lineages occasionally displayed enhanced growth on certain carbon sources, indicating that adaption to drug can provide a growth advantage in certain environments. A genome-scale metabolic network reconstruction (GENRE) of P. aeruginosa strain UCBPP-PA14 was paired with whole-genome sequencing data of one of the drug-evolved lineages to predict genes contributing to observed changes in metabolism. Finally, we experimentally validated in silico predictions to identify genes mutated in resistant P. aeruginosa affecting loss of catabolic function. Our results build upon previous mechanistic knowledge of drug-induced metabolic adaptation and provide a framework for the identification of metabolic limitations in antibiotic-resistant pathogens. Robust drug-driven changes in bacterial metabolism have the potential to be exploited to select against antibiotic-resistant populations in chronic infections.


Introduction
With the threat of a 'post-antibiotic era' looming, there is a critical need to develop new strategies to treat bacterial infections (1,2). The design of such approaches could be guided by a better understanding of the relationship between antibiotic resistance and bacterial metabolism. Bacterial metabolism has been shown to be an important factor in the efficacy of certain classes of antibiotics (3)(4)(5). For example, mutations to the electron transport chain have been shown to reduce proton motive force (PMF) and limit PMF-dependent influx of aminoglycosides (5,6). Such mutations are commonly observed in aminoglycoside-resistant clinical isolates (7,8). Conversely, quinolone efflux is often PMF-dependent and decreased PMF can result in increased drug susceptibility (9). Links between metabolism and antibiotic resistance have resulted in a variety of proposed treatments to clear infections including drug cycling approaches, which aim to continuously re-sensitize resistant populations, and metabolite supplementation strategies, which jumpstart metabolism to restore drug susceptibility in antibiotic tolerant cells (10)(11)(12). These tactics do not require new antibiotics, rather they help to prolong the efficacy of existing drugs through the manipulation of metabolism in resistant or tolerant populations.
A potential complementary approach to prolong the efficacy of existing antibiotics is to promote the metabolism and growth of antibiotic-sensitive populations over antibiotic-resistant populations, analogous to how prebiotics can support growth of beneficial microbial communities (13). However, application of this method relies on the discovery of robust resistance-specific metabolic limitations with known genetic mechanisms. Studies of the fitness and metabolic phenotypes of antibiotic-resistant bacteria focused around clinical isolates have yielded conflicting results (14)(15)(16). Resistance mutations that reduce antibiotic susceptibility have been shown to exhibit positive, negative, and null effects on bacterial fitness (17,18). The impact of individual mutations on fitness are further obscured by the presence of compensatory mutations (19). Reports of metabolic phenotypes of clinical isolates vary across studies and even within the same patient (14,16,20). The direct impact of sustained antibiotic pressure on metabolic adaptation is confounded by in vivo pressures, including nutrient stress, oxidative stress, host inflammation, and competition with co-infecting pathogens (21)(22)(23)(24). There has recently been success using adaptive laboratory evolution (ALE) experiments to study the evolution of antibiotic resistance in a controlled in vitro environment (8,(25)(26)(27)(28)(29)(30). Sequencing and expression profiling of lab-evolved resistant bacteria have revealed genetic mutations responsible for resistance phenotypes (31,32); however, connecting specific mutations to metabolic phenotypes remains a significant challenge (33,34).
Genome-scale metabolic network reconstructions (GENREs) can provide the framework to contextualize genetic and metabolic changes accompanying the development of antibiotic resistance (35,36). A GENRE is a quantitative formalism that captures all known metabolic reactions in an organism (37). Within a GENRE, gene-protein-reaction (GPR) rules link annotated metabolic genes to the reactions that their gene products catalyze (38,39). Among other functions, GENREs can be used to mechanistically evaluate the impact of mutations in single genes on bacterial growth across many environmental conditions (40). This analysis allows for the systematic prediction of the metabolic consequences of individual mutations identified in sequenced drug-evolved lineages. If they can be identified, robust resistant-specific changes in metabolism could be exploited to select for antibiotic-sensitive bacterial populations in chronic infections.
To better understand the relationship between antibiotic resistance and bacterial metabolism, we profiled the metabolic phenotypes of previously published lab-evolved antibiotic-resistant lineages of the opportunistic pathogen Pseudomonas aeruginosa. We evaluated growth of piperacillin-, tobramycin-, and ciprofloxacin-resistant P. aeruginosa on 190 unique carbon sources. We also evaluated growth of the starting ancestral lineage and a media-evolved lineage. This effort resulted in the generation of over 2,800 individual growth curves from which we captured resistance-specific changes in growth dynamics. Phenotypic data was integrated with previously collected genomic information on each lineage to probe the importance of reported resistance mutations on the observed metabolic phenotypes (8). Finally, phenotypic and sequencing data were used in tandem with a recently published GENRE of P. aeruginosa strain UCBPP-PA14 to predict the individual impact of 343 gene deletions in the piperacillinresistant lineage on loss of catabolic function (41). To the best of our knowledge, this study is the first to use this combined experimental and computational approach to predict and validate genetic mutations driving metabolic adaptation during the evolution of antibiotic resistance.
Altogether, we report that in vitro adaptation to antibiotics resulted in systems-level changes to metabolic function and growth dynamics in P. aeruginosa. We contextualized our data with a computational model to identify important genotype-phenotype relationships in the drug-evolved lineages and experimentally interrogated model-driven predictions. By improving our mechanistic understanding of the metabolic adaptations accompanying the evolution of antibiotic resistance, we aim to help guide the development of novel treatment strategies against bacterial pathogens.

Results
Profiling growth phenotypes of antibiotic-resistant Pseudomonas aeruginosa P. aeruginosa was previously evolved to lysogeny broth (LB) media through serial passaging for 20 days (8). In parallel, the same starting ancestral lineage P. aeruginosa was also evolved to each of three antibiotics: ciprofloxacin, piperacillin, and tobramycin ( Fig 1A). As previously reported, the minimal inhibitory concentration (MIC) of each drug measured in its respective drug-evolved lineage increased at least 32-fold while the LB-evolved control lineage had no increase in MIC to any of the drugs relative to the ancestral lineage (8). Phenotypic data and previously collected whole-genome sequencing data were paired with a genome-scale metabolic network reconstruction of P. aeruginosa UCBPP-PA14 to predict mutations impacting loss of catabolic function in the piperacillin-evolved lineage. In brief, the model was used to predict the individual impact of each gene that was deleted from the piperacillin-evolved lineage on the ability to grow on a given carbon source. If a simulated gene knockout resulted in a loss of model growth on a carbon source where an experimental loss of catabolic function was also observed, then the prediction was experimentally validated with a transposon mutant. This process was repeated for the piperacillin-evolved lineage on 41 carbon sources that could be both experimentally and computationally evaluated.
To characterize phenotypic changes in metabolism that arise with resistance, we evaluated the growth of each lineage on 190 unique carbon sources in the absence of drug-pressure in triplicate ( Fig 1A). These growth phenotyping experiments resulted in the generation of 2,850 individual growth curves (S1 Data). Curves for each lineage on each carbon source were averaged across three biological replicates and these 950 averaged growth curves are referred to for the remainder of the analysis (S1 Fig  To better understand growth differences across lineages, we measured and summarized the maximum growth densities across growth-supporting carbon sources for each lineage (Fig 2C,   S2 Fig, S2 Data). The maximum growth density was defined as the maximum optical density measured at 600nm (OD600) on each averaged 48-hour growth curve after background subtraction. Interestingly, all evolved lineages, including the LB-evolved control, exhibited a significantly decreased average maximum growth density relative to the ancestral lineage (Pvalue < 0.001). Each lineage was found to have a different number of growth-supporting carbon sources ( Fig 2D). The ancestral lineage grew on 33.2% (63/190) carbon sources.
Comparatively, ciprofloxacin-and piperacillin-evolved lineages were the most metabolically limited, only growing on 18.4% (35/190) and 24.7% (47/190) of the tested carbon sources, respectively. Somewhat unexpectedly, the LB-evolved control was only able to grow on 28.4% (54/190) of carbon sources, while tobramycin-evolved lineage grew on 32.6% (62/190) of sources. Overall, adapted lineages were more metabolically limited and had decreased growth relative to the unevolved ancestral lineage.

Development of antibiotic resistance results in altered growth dynamics in P. aeruginosa
To more rigorously examine metabolic differences of our antibiotic-resistant lineages, we quantified the growth dynamics of each individual growth curve across carbon sources that supported growth of all five lineages ( Fig 1A, S3 Data). For each of the 20 growth-supporting carbon sources, three parameters of bacterial growth were calculated: maximum growth density, growth rate, and time to reach mid-exponential phase (a proxy for the duration of lag phase) (Fig 3). The maximum growth density varied by carbon source and by lineage ( Fig 3A).
All evolved P. aeruginosa lineages except for the tobramycin-evolved lineage exhibited significant decreases in maximum growth density relative to the ancestral lineage (P-value < 0.001) (Fig 3D). The LB-evolved control had the largest decrease in growth density, indicating that media condition impacted the maximum level of growth perhaps more so than antibiotic pressure.
For each growth curve, we also calculated the growth rate and the time taken to reach midexponential phase (Fig 1A). The irregular shape of many of the growth curves made it difficult for existing growth curve analysis software to determine the growth rate. To better analyze our data, we adapted a previously published sliding-window algorithm (42). The modified algorithm determined the growth rate to be the maximum slope of the natural log of the growth curve. The slope was calculated across a minimum of eighty minutes (eight time points). The time to reach mid-exponential phase was defined as the time at which the culture reached its maximum growth rate (S1 Code, Methods). Across the 20 universal growth-supporting carbon sources, all antibiotic-evolved lineages grew at significantly slower rates than the ancestral lineage, with piperacillin-evolved and ciprofloxacin-evolved lineages exhibiting the largest decreases in growth rates (P < 0.01) (Fig 3B, E). The average growth rate of the LB-evolved control also decreased; however, this change was not significant (P > 0.01). The time to mid-exponential phase was defined as the beginning of the window when the growth rate was calculated. The ancestral, LB-evolved, and tobramycin-evolved lineages all reached mid-exponential phase in under 5 hours. Conversely, the piperacillin-evolved lineage was characterized by significantly longer lag phases, taking up to 29 hours to reach mid-exponential growth (Fig 3C, F).
Ciprofloxacin-resistant P. aeruginosa growth curves were the most irregular in shape and therefore lag times were more varied, with growth on some carbon sources beginning earlier than ancestor and growth on other carbon sources delayed more than the piperacillin-evolved lineage. The cause of this variability, particularly whether it is biology-driven or noise-driven, remains to be determined. These findings demonstrate that adaptation to different antibiotics can heterogeneously impact the growth dynamics of P. aeruginosa across many growthsupporting conditions.

Tobramycin-resistant P. aeruginosa shows enhanced growth on N-acetyl-D-glucosamine
On multiple carbon sources, we observed that an antibiotic-evolved lineage exhibited enhanced growth relative to the ancestral lineage. This phenotype is potentially of greater concern in a clinical setting since a drug treatment that selects for both resistant and metabolically advantaged isolates could have negative clinical outcomes.
One example of enhanced growth occurred in the tobramycin-evolved lineage, which grew to a higher growth density than the ancestral lineage on the carbon source N-acetyl-D-glucosamine (Fig 4). N-acetyl-D-glucosamine is a component of human mucin as well as the cell wall of Gram-positive organisms, including S. aureus (43,44). Consistent with phenotypic profiling data, all four biological replicates of tobramycin-evolved P. aeruginosa from a previous study (8) showed enhanced growth on 20mM of N-acetyl-D-glucosamine relative to the ancestral lineage (S3 Fig). One replicate in particular reproducibly showed enhanced growth beyond what was observed in the initial phenotypic screening. This replicate also contained a single-base deletion causing a frameshift in the nuoL gene encoding NADH Dehydrogenase I. Mutations in the nuo complex were seen in multiple drug-evolved lineages previously described (8) and have been previously associated with aminoglycoside tolerance (7,8,45,46). Mutations to the nuoL gene as well as other disruptions to the electron transport chain and oxidative phosphorylation have been shown to reduce proton pumping, likely preventing aminoglycoside influx into the cell and lessening drug efficacy (6,47,48). To determine the impact of the nuoL mutation of the tobramycin-evolved lineage on N-acetyl-D glucosamine utilization, we selected a nuoL transposon mutant from the non-redundant PA14 library (49). Similar to the tobramycin-evolved lineage, the nuoL transposon mutant showed enhanced growth on N-acetyl-D-glucosamine relative to the ancestral lineage (Fig 4).
The nuoL transposon mutant also exhibited an increased MIC to tobramycin relative to wild type P. aeruginosa, although based on defined clinical breakpoints, the MIC was not elevated enough for the mutant to be considered resistant (Table 1) (50). Inhibition of the electron transport chain by the inhibitor carbonyl cyanide m-chlorophenyl hydrazone (CCCP) did not impact wild-type growth on N-acetyl-D-glucosamine, suggesting that oxidative phosphorylation was not primarily responsible for this phenotype (S4 Fig). The exact mechanism by which nuoL has enhanced growth on N-acetyl-D-glucosamine remains to be determined. However, these results demonstrate one example of how resistance-associated genetic mutations can unexpectedly impact metabolism in P. aeruginosa. Second, the P. aeruginosa metabolic network reconstruction had to contain both an exchange reaction and a transport reaction for the carbon source. An exchange reaction can be thought of as a way to include a metabolite in the simulated media while a transport reaction can be viewed as a way for P. aeruginosa to uptake or secrete a particular metabolite. Finally, for a carbon source to be included in the analysis it had to support the production of a non-zero biomass in the complete model, where biomass can be considered a proxy for bacterial growth.
There were 41 carbon sources that met these requirements. On each of these carbon sources, we simulated single-gene knockouts to evaluate the contribution of each gene on model growth and carbon source catabolism. Genes that when knocked out in the model resulted in a loss of biomass production were predicted to be essential for growth on the tested carbon source.
To predict which genes in the large deletion of the piperacillin-evolved lineage could have impacted observed changes in metabolic phenotype, we looked for overlap between deleted genes identified by sequencing and model-predicted essential genes (S5 Fig). In total, there were 17 genes shown to be deleted in the piperacillin-evolved lineage that were predicted to be essential for growth on at least one of the 41 experimentally measured carbon sources (Fig 5).
For it to be possible to validate an essentiality prediction, the ancestral lineage needed to be able to grow on the carbon source of interest. Otherwise, we could not determine the effect of knocking out a single gene on catabolic function. We were unable to evaluate the effect of six of the 17 genes because they were predicted to be essential for growth on three carbon sources that were unable to support growth of the ancestral lineage: D-ribose, D-serine, and L-serine.
Two more genes, fahA and hmgA, were excluded from further validation because ancestral growth on L-phenylalanine was just above the defined AUC cutoff value. Another two of the 17 genes, bacA and glgA, were predicted to be essential across all environmental conditions. While neither of these genes were essential for growth in vitro, they were associated with the ability of the model to produce biomass. Predictions of the knockout of these genes, while perhaps not directly informative for our study of antibiotic resistance, provide useful new information for model curation. Of the remaining seven predicted essential genes, five genes deleted in the piperacillin-evolved lineage were predicted to be essential for L-leucine utilization: gnyA, gnyB, gnyD, gnyH, gnyL. A cluster of PAO1 genes orthologous to the gnyABDHL cluster has been previously reported to be involved in L-leucine catabolism, specifically in the downstream catabolism of isovaleryl-CoA into the citric acid cycle substrate acetyl-CoA (Fig 6A) (51). The piperacillin-evolved lineage was found to lack the ability to grow on L-leucine while the ancestral strain was able to grow on Lleucine ( Fig 6C). Based on the predictions of the genome-scale metabolic network reconstruction, we hypothesized that this loss of catabolic function in piperacillin-evolved P.
aeruginosa was due to the deletion of these five genes in the gny operon. To test the validity of these predictions, we measured growth of a gnyA transposon mutant in minimal media with 40mM of L-leucine. We found that like the piperacillin-evolved lineage, optical density of the gnyA mutant in culture remained constant over 48 hours (Fig 6C). These results further suggest that gnyA is necessary for L-leucine utilization. Our model predicts that deletion of any one of the remaining gny genes would have a similar effect on the ability to utilize L-leucine. The gnyA gene was also predicted to be essential for utilization of L-isoleucine; however, the piperacillin-evolved lineage was able to grow on this carbon source (Fig 5). Because this lineage contained a large deletion with over 300 genes, we could not rule out that gnyA was independently essential but that added mutations had compounding effects which restored growth. To elucidate the role of gnyA on L-isoleucine utilization, we tested growth of the gnyA transposon mutant on 20mM of L-isoleucine. P. aeruginosa was able to grow on this carbon source without a functioning gnyA gene, indicating that this gene was a false positive essential gene in our model (S6 Fig). The exact reason that gnyA was erroneously predicted as essential remains to be determined and indicates that the current understanding of L-isoleucine utilization in P. aeruginosa is incomplete. Finally, our model identified two genes deleted in the piperacillin-evolved lineage as essential for utilization of 4-hydroxybenzoic acid (Fig 5). These genes, scoA and scoB, have been reported to encode for subunits A and B of a CoA-transferase (52). Unlike the gnyABDHL cluster, which has been previously associated with L-leucine catabolism, it was not immediately obvious how scoA and scoB impacted 4-hydroxybenzoic acid utilization. According to the P.
aeruginosa GENRE, scoA and scoB use a downstream product of 4-hydroxybenzoic acid degradation to catalyze the conversion of succinyl-CoA to succinate in the citric acid cycle ( Fig   6B). To validate the prediction that the piperacillin-evolved lineage was unable to utilize 4-hydroxybenzoic acid due to deletion of scoA and scoB, we attempted to grow mutants with transposon insertions in each gene on 20mM of 4-hydroxybenzoic acid. We found that both mutants were able to grow on 4-hydroxybenzoic acid. Growth curves of both mutants matched the growth of the ancestral lineage while the piperacillin-evolved lineage was unable to grow ( Fig 6D). From this result we conclude that scoA and scoB are not individually necessary for growth on 4-hydroxybenzoic acid. Further evaluation is required to determine why these genes were falsely predicted to be essential.
Altogether, we paired whole-genome sequencing with a genome-scale metabolic network reconstruction to identify five genes associated with loss of metabolic function in piperacillinresistant P. aeruginosa. Three additional genes predicted to be essential by the model were experimentally determined to be non-essential, highlighting gaps in our current knowledge of metabolism that need to be further investigated. While some of our predictions may be retrospectively obvious, they were all guided by the model, without which we would have been unable to reconcile direct genotype-phenotype relationships in such a complex mutational landscape. A similar approach could be applied to predict metabolic deficits in any bacterial strain with an annotated genome.

Discussion
We have shown that P. aeruginosa experiences systems-level changes to metabolism when exposed to sustained antibiotic pressure. Through the semi-automated analysis of over 2,800 individual growth curves on 190 unique carbon sources, we determined that antibiotic-resistant P. aeruginosa exhibit differential growth dynamics from antibiotic-sensitive as well as other and Gram-positive cell walls and is assumed to be present in sputum of CF patients (43,44).
Once taken up by P. aeruginosa, N-acetyl-D-glucosamine can be catabolized to produce energy for growth or used to signal the production of the virulence factor pyocyanin (43,55). While we were unable to identify the specific mechanism driving enhanced growth of a nuoL mutant on Nacetyl-D-glucosamine, we speculate that this mutation disrupts pyocyanin production and consequently has compensatory effects on catabolism and growth. In other words, there may exist a tradeoff between growth and virulence factor production that is sensitive to disruption of the nuoL gene. This hypothesis is supported by the finding that mutation and downregulation of genes in the nuo operon have been associated with decreased pyocyanin production (56,57) and more generally that mutations to the electron transport chain result in reduced virulence (58).
Adaptation resulted in significantly decreased average growth rates across all antibiotic-evolved lineages ( Fig 3B). While growth rate was not significantly decreased in the LB-evolved control lineage, we did observe a significant decrease in maximum growth density relative to the ancestral lineage, indicating that evolving to media and antibiotic at the same time may have had confounding effects on metabolism. Media condition has been shown to impact evolution of resistance as well as growth phenotypes in the absence of drug (28,59). For future adaptive laboratory evolution studies, we recommend that the ancestral lineage be first adapted to the base media and then subsequently adapted to the antibiotic of interest to control for the effect of media (60).
Genome-scale metabolic network reconstructions are becoming increasingly prevalent in the study of antibiotic resistance (25,28,61,62). Although the P. aeruginosa metabolic network reconstruction does not account for many genes associated with resistance (e.g. regulation of efflux, membrane permeability), the model allowed for the rapid identification of potential mutations impacting growth phenotypes. Erroneous predictions provided useful starting points for future model curation and increased understanding of metabolic functionalities. Notably, the majority of the overlap between model genes and sequenced mutations occurred in the large deletion of the piperacillin-resistant lineage. Due to the limited number of mutations to metabolic genes across the other evolved lineages, we were unable to find specific genetic mechanisms for many of the experimentally observed changes in metabolic phenotype. We speculate that metabolic adaptations that could not be explained by our model simulations may have been affected by resistance mutations that impact regulation. We would expect the integration of a transcriptional regulatory network or transcriptomic profiling data with the metabolic network reconstruction to uncover more genotype-phenotype relationships in the resistant lineages (61,(63)(64)(65)(66).
While this study shows that there are systems-level metabolic changes following sustained antibiotic pressure, ultimately, evolution of resistance is a stochastic process that is sensitive to a wide variety of environmental and host factors not accounted for in our experimental design.
As such, observed mutations may have varied with a different selection of antibiotics or media conditions. Nevertheless, the mutations we observed were consistent with those in sequenced clinical P. aeruginosa isolates, indicating that the general trends we observed may be robust to nutrient differences between LB media and CF sputum (14,47). While we have chosen to focus on three antibiotics and one media condition in a single bacterial species, our combined experimental and computational approach can be applied to a large variety of organisms and environmental pressures. Moving forward, we have laid the groundwork for the interrogation of broad metabolic consequences of antibiotic resistance.
There are many potential advantages to investigating how bacterial metabolism is impacted by antibiotic treatment. Clinically, it may be of value to understand how a prescribed antibiotic impacts the growth of an infecting pathogen. For example, a clinician may decide to administer a drug that limits the metabolic flexibility of a resistant pathogen over another drug that would promote a metabolic advantage. From an engineering perspective, we propose that robust metabolic changes can be leveraged in the design of new antibiotic treatment strategies. For instance, supplementation with metabolites that cannot be catabolized by resistant pathogens could help to restore antibiotic sensitivity in chronic infections. Future work should focus on exploiting antibiotic-driven metabolic adaptations to mitigate the development of antibiotic resistance.

Bacterial strains and culture conditions:
Pseudomonas aeruginosa strain UCBPP-PA14 was previously evolved to one of three antibiotics (piperacillin, tobramycin, or ciprofloxacin) in lysogeny broth (LB) media for 20 days (8). As a control, UCBPP-PA14 was also adapted to LB media in the absence of drug for 20 days. One replicate of each 20-day evolved lineage described previously (8)

Semi-automation of calculation of growth dynamics:
Three growth parameters (growth rate, maximum optical density, and time to mid-exponential phase) were calculated semi-automatically with custom scripts in MATLAB (Mathworks, R2016b). Growth rate was calculated using a sliding window algorithm modified from (42). The algorithm was modified to identify the maximum slope of ln(OD600) vs. time prior to reaching lag phase. A window size of 80 minutes (eight time points) was used. Prior to calculating the maximum optical density, the average of the negative control at each time point was subtracted from the average of three replicate growth curves for each lineage on each carbon source. The maximum optical density was then reported as the maximum of this mean backgroundsubtracted curve. The time to mid-exponential phase was defined as the earliest time point included in the calculation of growth rate; it was assumed that mid-exponential phase begins when the culture reaches its maximum growth rate. The value was averaged across three replicates. All code used to calculate growth dynamics is available in the supplemental materials (S1 Code).

Genome-scale metabolic network reconstruction of Pseudomonas aeruginosa:
The recently published genome-scale metabolic network reconstruction of wild-type P.
aeruginosa UCBPP-PA14 was used for this study (41). The model accounts for the functions of 1129 genes and 1495 reactions. Previous analysis revealed that the model was 81% accurate at predicting growth phenotypes of 91 experimentally measured carbon sources (41).

Gene essentiality predictions:
69 carbon sources of the 190 measured carbon sources had exchange reactions in the genome-scale metabolic network reconstruction, iPau1129. An exchange reaction represents the ability to import a carbon source into the cell from the media. We temporarily added exchange reactions for another ten carbon sources that had transport reactions but lacked exchange reactions. The base model supported growth (non-zero biomass) on 41 of these 79 carbon sources. We used iPau11290 to simulate growth of P. aeruginosa in each of these carbon sources with minimal media. Using the singleGeneDeletion function in the Cobra Toolbox (67), we simulated single-gene knockouts of each gene in the presence of each unique carbon source. From these simulations, we predicted genes for each media condition that when knocked out, resulted in the model no longer being able to support growth. We refer to genes predicted to be required for growth as essential genes for a given environment. We looked for overlap between essential genes predicted by the model and mutated genes in our evolved resistant strains. This overlap was used to predict resistance-specific essential genes as well as to identify incorrect model predictions.

Statistical tests:
Statistical comparisons between the wild type ancestor and all evolved lineages were made with two-tailed Wilcoxon rank sum tests (p-value < 0.01) (S1 Text). Samples were assumed to be independent.

Data availability:
MATLAB code used to calculate the growth rate and time to mid-exponential phase of each growth curve can be found in the supplemental materials (S1 Code). Code to reproduce gene essentiality predictions (Fig 5, S5 Fig, and S4 Data) are also included (S3 Code). The R script and all data used to generate figures and supplemental data (excluding Fig 1 and S7 Fig) are additionally provided (S2 Code). Complete growth data, binary growth calls, growth dynamics, The negative control contains media that was not inoculated (n = 6 colonies of gnyA mutant, error = standard deviation). had a product between 1-2kb.