The Tolerome: A Database of Transcriptome-Level ... - ACS Publications


The Tolerome: A Database of Transcriptome-Level...

0 downloads 66 Views 2MB Size

Subscriber access provided by Gothenburg University Library

Article

The Tolerome: A database of transcriptome-level contributions to diverse Escherichia coli resistance and tolerance phenotypes Keesha E Erickson, James D. Winkler, Danh Nyugen, Ryan T. Gill, and Anushree Chatterjee ACS Synth. Biol., Just Accepted Manuscript • DOI: 10.1021/acssynbio.7b00235 • Publication Date (Web): 10 Oct 2017 Downloaded from http://pubs.acs.org on October 11, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

ACS Synthetic Biology is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

1The Tolerome: A database of transcriptome-level contributions to diverse Escherichia coli 2resistance and tolerance phenotypes 3 4Keesha E. Erickson1,‡, James D. Winkler1,§, Danh T. Nguyen1, Ryan T. Gill1, Anushree Chatterjee1,* 51 Chemical & Biological Engineering, University of Colorado, Boulder, Colorado, 80301, USA 6 7Abstract 8

Tolerance and resistance are complex biological phenotypes that are desirable bioengineering

9goals for those seeking to design industrial strains or prevent the spread of antibiotic resistance. Over 10decades of research, a wealth of information has been generated to attempt to decode a molecular basis 11for tolerance, but to fully achieve the goal of engineering tolerance, researchers must be able to easily 12learn from a variety of data sources. To this end, we here describe a resource designed to enable scrutiny 13of diverse tolerance phenotypes. We have curated hundreds of gene expression studies exploring the 14response of Escherichia coli to chemical and environmental perturbations, from antibiotics to biofuels and 15solvents and more. Overall, our efforts give rise to a database encompassing more than 56,000 gene 16expression changes across 89 different stress conditions. This resource is designed for compatibility with 17the Resistome database, which includes more than 5,000 strains with mutations conferring resistance or 18sensitivity but no transcriptomic data. Thus, the work here results in the first combined resource 19specialized to tolerance and resistance in E. coli that supports investigations across genomic, 20transcriptomic, and phenotypic levels. We leverage the database to identify promising bioengineering 21targets by searching globally across multiple stress conditions as well as by narrowing the focus to fewer 22conditions of interest, such as biofuel stress and antibiotic stress. We discuss some of the most frequently 23differentially expressed genes, present genes commonly differentially co-expressed, and predict which 24transcription factors and sigma factors most likely contribute to gene expression profiles in a wide array 1 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

25of conditions. We also compare profiles from sensitive and resistant strains, gaining knowledge of how 26responses differ per overrepresented gene ontology terms. Finally, we search for genes that are frequently 27differentially expressed but not mutated, with the expectation that these may present interesting targets for 28future engineering efforts. The curated data presented here is publicly available, and should be 29advantageous to those studying a variety of bacterial tolerance phenotypes. 30Keywords: bacteria, transcriptome, adaptive resistance, tolerance, database 31

Bacterial stress tolerance and resistance are important considerations for applications ranging

32from strain engineering to antibiotic development. Tolerance and resistance are complex, multi-genic 33phenotypes such that efforts to engineer tolerance often harness global search or directed evolution 34techniques1–3. However, new –omics datasets evaluating tolerance to various environmental conditions 35and chemical compounds are being generated at a rapid pace. To make predictions about the emergence 36of antibiotic resistance, or to strive for rational design of tolerant strains, it is critical that the vast amounts 37of data that have been produced be assembled in a manner that facilitates analysis. 38

There are databases currently available that were conceived to facilitate query into bacterial

39resistance. These tools, including Resfams4, ARG-Annot5, CARD6, Resfinder7, and SEAR8, focus on 40DNA-sequence level information, or “genetic” information, which enables recognition of encoded 41resistance genes, collectively allowing for predictions of a strain’s drug resistance profile. They are 42designed for cases where resistance is conferred by a specialized enzyme, binding protein, or efflux 43pump, such as is common in many cases of β-lactam, aminoglycoside, and tetracycline resistance. 44Excepting ARG-Annot, these databases are not intended to predict resistance in instances where a drug is 45evaded via point mutations in the target site, as often occurs in fluoroquinolone and rifampicin 46resistance9,10. Likewise, these databases do not allow for predictions of tolerance to compounds beyond 47antibiotics, or allow for prediction of gene expression changes that might change a strain’s tolerance 48profile. Centralized repositories such as the Sequence Read Archive11 and the Gene Expression 49Omnibus12 and database like the Many Microbes Microarrays Database13 provide access to large numbers 2 ACS Paragon Plus Environment

Page 2 of 35

Page 3 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

50of transcriptome datasets, or “non-genetic” information, either from microarray or RNA sequencing 51experiments, and collectively provide an invaluable resource for analyzing phenotypes at the 52transcriptional level. These data are deposited using a plethora of organizational and normalization 53approaches, and no single resource has yet incorporated both genetic (mutations) as well as and non54genetic (gene expression) data to facilitate assessment of and across diverse tolerance phenotypes in E. 55coli. 56

Non-genetic contributions are increasingly recognized as important contributions to tolerance, as

57inherent transcriptional programs influence phenotype without the need for mutation. Biological systems 58are able to tolerate stress through the implementation of innate responses at the transcriptional and post59transcriptional level14. Non-genetic impacts on phenotype have been revealed across biological kingdoms, 60from human cancers15 to bacterial decision-making16. Non-genetic responses have the advantage of being 61able to generate diversity within a population of genetically identical cells, which can provide an 62advantage in a search for behaviors that allowing survival in a challenging environment17–19. Furthermore, 63non-genetic responses can be reversible, allowing cells to temporarily enter states, such as persistence, 64that would otherwise be costly or prohibitive for long-term survival20. Bacterial adaptive resistance is a 65non-genetic phenomenon in which slight changes in gene expression levels in response to stress promote 66transient tolerance increases21. Conversely, synthetic manipulations of bacterial gene expression have 67been shown to have significant influence on tolerance to broad classes of chemical stressors22–24. Due to 68the power of innate non-genetic responses, it is not possible to pinpoint all genes involved in bacterial 69response to stress through an analysis of the mutational landscape alone25. As such, to better understand 70and engineer resistance, it is crucial to examine how programming reflected at the gene expression level 71influences resistance phenotypes. 72

To address this gap, we here introduce a new, curated database resource that includes

73transcriptome-level information from microarray, RNA-sequencing, and quantitative PCR datasets while 74maintaining records on mutations in a variety of stress conditions. To enable related exploration of gene 3 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

75expression datasets, we built upon the platform provided by the Resistome database, which we have 76previously introduced as a means for analysis of mutation-associated resistance phenotypes in 77Escherichia coli26. The first Resistome release contained 5,511 E. coli strains with mutations in 3,639 78unique genes, delivering data on 228 distinct resistance phenotypes, from exposure to biofuels, solvents, 79antibiotic agents, and environmental perturbations. We leveraged the Resistome database to identify 80genotype to phenotype connections across multiple stress conditions, revealing trends in random genomic 81mutations regarding the propensity for nonsynonymous codon changes, and to predicting important genes 82and gene ontologies for future study of tolerance phenotypes. Here, we report significantly expanded 83capabilities beyond the analysis of resistance-conferring mutations. The resource described in this report 84contains more than 56,000 gene expression changes recorded in response to 89 different stress conditions. 85We exploit this dataset to gain unique insight into E. coli tolerance and resistance, including by 86identifying genes and pathways implicated in specific stress conditions, characterizing phenotypes by the 87sigma factors or transcription factors that most influence gene expression, and improving the selection of 88bioengineering targets by considering gene expression promiscuity. 89Methods 90Database organization 91

We have organized the data as flat key/value text files, with one file per journal article. Records

92contain information gathered from the respective paper, including species and strains used, growth 93conditions (temperature, media, oxygen), specifics on mutations (gene, codon change), observed 94resistance levels to various drugs or substances, as well as other potentially useful information. The 95database is a series of classes with specific attributes. The Paper class is at the top level, and is the parent 96of the Mutant class. Each Mutant has children representing mutations within that mutant (Mutation 97classes). A summary of the attributes in each of these classes is in Fig. 1A. To incorporate gene 98expression changes, we added a class (GEChange) at the same level as the Mutation class. The GEChange 99class has recommended attributes describing the gene that was differentially expressed and the type of 4 ACS Paragon Plus Environment

Page 4 of 35

Page 5 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

100change relative to the respective control condition (either over- or underexpressed). Additional attributes 101are provided based on their disclosure in the respective paper, including the fold change in gene 102expression, P-value and statistical test used, exposure time (in hours), growth phase, amount of stress 103applied, the accession number for the dataset in another database. All attributes that can be used to 104describe gene expression data are summarized in Table 1. 105

It should be noted that the Tolerome does not currently normalize across studies, so gene

106expression fold-change values are each with respect to whatever control or reference was utilized in the 107associated manuscript. As such, caution should be exercised when comparing fold-change values across 108disparate experiments, as the baseline gene expression state for different reference conditions is likely to 109be dissimilar. Additionally, genes expressed below the detection limit, which may or may not be 110differentially expressed, were assumed to be not differentially expressed and therefore not included in 111records. The terms “overexpressed” or “underexpressed” here indicate that genes were reported as 112significantly differentially expressed with respect to the control in that particular study. The greatest value 113of the database is being able to pinpoint trends in gene expression that persist even across different strains 114and/or in different conditions. 115

We have added one attribute within the Paper class that was not present in the initial Resistome

116release. This attribute, named GEQC, provides a measure of confidence for each gene expression dataset, 117as garnered from the information presented in the associated paper. The GEQC value ranges between 0 118and 5. A value of 0 indicates that the study has no gene expression data. A value of 1 indicates that the 119study presented the lowest level of quality control, 2 is below average, 3 is average, 4 is above average, 120and 5 is the highest level. The assignment for each paper was manually determined based on the paper’s 121analytical approach, including descriptions of experimental design, normalization techniques, sequencing 122depth, and statistical analysis. For instance, to earn an average score of 3, a microarray study needed to 123demonstrate the usage of replicates, statistical analysis, and testing of the quality of the input sequence. 124Of the records in the database, 45% earned a score of 3. To earn a 4, as did 21% of the records, a 5 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

125microarray study might do all the QC needed to earn a score of 3, then also perform and additional 126validation, such as the inclusion of qPCR data. A score of 5 would be earned if the study demonstrated 127truly exceptional validation, perhaps including very rigorous qPCR study with multiple reference genes. 128Scores of 2 were given if the study did not clearly indicate validation, as could be the case if a RNA129sequencing study did not provide depth information, or if a microarray study called differentially 130expressed genes solely by fold-change values (23% of records earned a 2). The lowest score of 1 was 131given to 8% of records (e.g. if a study only relied on qPCR data, without mention of standard controls 132such as no-template controls or no-reverse transcriptase controls). 133Literature curation 134

Entries were found via literature search as well as by searching NCBI’s Bioproject repository for

135instances with “transcriptome” data types, “SRA” project data, and with taxonomy “Escherichia coli” as 136well as by browsing NCBI’s GEO. Only SRA or GEO data with associated publications were included. 137To be incorporated in the database, studies had to include E. coli transcriptome data in the presence of an 138external or environmental stress. For instance, we did not include studies that introduced a mutation and 139measured the impact on global gene expression. Every record has at least one mutation(s) or gene 140expression change(s). 141Software 142

Software for loading and analyzing the database is written in Python 2.7. Many subpackages were

143utilized to perform the analysis described here, including scipy (v0.18.1), pandas (v0.19.2), and 144statsmodels (v0.6.1). Subpackages were installed through the Anaconda distribution (v4.3.9). We applied 145MATLAB’s bioinformatics toolbox to build heatmaps and dendrograms. Circos27 was used to generate 146circle plots. OriginPro was used to prepare other charts and graphs as well as perform statistical testing 147for linear correlation. An example of how to load the database is provided in the Supporting Text. 148Compatibility 6 ACS Paragon Plus Environment

Page 6 of 35

Page 7 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

149

ACS Synthetic Biology

Records in the initial Resistome release did not generate or include gene expression data. To

150make these records compatible with the Tolerome database described here, the only modification required 151was to add a line for each mutant indicating that the number of gene expression changes was 0. There 152were 33 records in the Resistome that did perform transcriptome analysis. Their original records only 153contained data regarding genomic changes, so they were fully updated here to incorporate transcriptome 154information. 155Data accession 156

Software and data are available at http://www.colorado.edu/UCB/chatterjeelab/Databases.html.

157Please refer to the Supporting Information t for a brief tutorial of how to load and query the database. 158Results & Discussion 159Database Description 160

Our efforts result in a resource that contains 174 Escherichia coli gene expression datasets

161published between 2000 and 2016, in addition to mutation-level data previously described26 , with sixteen 162new records added containing only mutational data and no gene expression data. In total, 56,263 gene 163expression changes in 5,049 unique genes are included. Gene expression changes were mainly measured 164using either microarray (49.8% of changes) or RNAseq (49.6%), with a small fraction of qPCR studies 165included (0.5%). As expected, earlier studies used microarray, while more recent studies migrate toward 166RNAseq (Fig. 1B). 167

Gene expression studies include 89 different stress conditions. Grouping the stress conditions by

168broad categories provides a summary of the types of conditions considered (Fig. 1C). The largest 169proportion of gene expression changes included in the database are associated with solvent/biofuel 170exposure (36% of gene expression changes), which includes strains challenged with isobutanol, ethanol, 171and n-butanol, among others. Other stress categories that comprise larger portions of the database include 172growth-associated stresses, such as growth in minimal media, at 13% of all changes, antibiotic treatment 7 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

173(12%), and exposure to organic acids (10%). The category of “other” types of stress also constitutes a 174large portion, at 12%, and includes a wide range of conditions that were not easily grouped with other 175categories, such as from E. coli grown on different food sources. As described in the previous release of 176the database26, stress categories include tags that indicate whether a mutant or strain is resistant or 177sensitive to the stress condition. 178

Gene expression data is provided from 34 different E. coli strains, including strains commonly

179used in research (MG1655, DH5α) and pathogenic isolates like O157:H7. The three strains with the 180highest representation in the database are MG1655, which accounts for 49% of gene expression changes, 181BW25113 at 18%, and W3110 with 13%. Fig. 1D depicts the proportion of gene expression changes per 182strain that are associated with each stress category. The strains with data across the highest number of 183stress categories are MG1655 (13 categories), BW25113 (6 categories), and O157:H7 EDL933 (6 184categories). 185Genes and functional classes over-represented in E. coli exposed to stress 186

The identification of genes, pathways, and functional classes that are commonly differentially

187expressed is one avenue of locating promising bioengineering targets. To locate the genes most 188commonly differently expressed, we searched across all differential gene expression records in the 189database, counting the number of occurrences for each gene and sorting the list from highest to lowest. 190Fig. 2 depicts an overview of gene expression data from all conditions, plotted on the E. coli MG1655 191chromosome. The histogram plot along the outer radius displays the number of instances in which the 192corresponding gene was differentially expressed. The average number of expression changes per gene is 19311.1±9.4, with a maximum of 74 and a minimum of 0. Out of the 4322 annotated coding genes and RNA 194genes in the MG1655 genome, 145 were not differentially expressed in any record. The seven genes most 195often differentially expressed were dps, osmY, rmf, elaB, manX, marA, and ygiW (Fig. 2). Table 2 196provides more detail on these genes and the conditions in which they were identified. Genes that are often

8 ACS Paragon Plus Environment

Page 8 of 35

Page 9 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

197differentially expressed are likely to be associated with general stress response pathways that convey the 198ability to survive in a variety of conditions. 199

To identify genes most often differentially co-expressed across all stress conditions, we generated

200a list of all sets of two genes, then searched through all records and tallied the number of strains in which 201both pairs were differentially expressed in either direction. The inner radius of Fig. 2 displays linkages 202between genes that were most often differentially expressed within the same strain. Data shown here is 203from the 2080 most frequently occurring pairwise co-occurrences, or in other words, the gene pairs 204differentially expressed together in at least 18 different gene expression datasets. Here, otsA was 205differentially co-expressed with 372 other genes (appearing differentially expressed with each of the 372 206genes in at least 18 different datasets), acnA with 322 genes, ndk with 312, otsB with 158, atpE and atpG 207with 109, atpF with 107, flhD with 89, and rpsP with 80 genes. The gene pairs most commonly co208occurring were otsA/osmY, which were differentially expressed in the same strain in 33 different datasets, 209acnA/osmY in 32 different datasets, and ndk/osmY in 31 different datasets. This analysis was performed in 210a manner naïve to whether co-expressed genes share an operon. Examining this possibility, we find that 211each of these three most prevalent co-expressed gene pairs are not co-transcribed. Indeed, none of the top 212nine most frequently co-expressed gene pairs are co-transcribed (Supporting Data File). 213

We also examined co-expression in the same direction, tallying gene pairs that were found to be

214both overexpressed in the same strain or underexpressed in the same strain. Interestingly, the six most 215prevalent overexpression pairs, which were overexpressed together in 18 to 21 gene expression datasets, 216included the inner membrane protein gene ycjF co-expressed with a molecular chaperone gene: either 217clpB, ibpB, dnaJ, htpG, groL, or dnaK. None of these genes share an operon. The most frequently tallied 218underexpression pairs often exhibited functions related to motility, with the top pair being the chemotaxis 219proteins tap and cheW, co-underexpressed in 22 datasets. These genes are in adjacent operons. Many of 220the other genes in these operons were found to be co-underexpressed with high frequency, including tar 221and cheA (in 21 cases) and tar and tap (21 cases). Motility proteins from other operons were also evident, 9 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 35

222as in the frequently tallied pair of fliZ and cheW (21 cases). Co-overexpression and co-underexpression 223counts for every pair of genes are provided in the Supporting Data File. 224

The gene most often differentially expressed across all datasets was dps. Dps is primarily

225recognized for binding to and protecting DNA from reactive oxygen species during stationary phase28, 226though it has also been linked to response to oxidative stress during exponential phase29. All 74 instances 227in the database in which dps was found to be DE were from cells in exponential phase. It was 228underexpressed in 64% of instances across 9 different stress conditions, and overexpressed in 19 unique 229stress conditions (indicated in Table 2). The direction of differential expression may be related to 230exposure time. dps was underexpressed in 100% of cases where cells were harvested after at least 72 231hours of exposure to the toxin, but for exposure times less than or equal to 8 hours, dps was 232overexpressed in 79% of cases (Fig. 3A). Four other commonly DE genes (from Table 2) also 233demonstrated this pattern of overexpression at shorter exposure times and underexpression after longer 234durations: osmY, rmf, elaB, and ygiW. 235

The mannose PTS permease gene manX was predominantly underexpressed at both short and

236long durations (Fig. 3B), and the regulator marA was overwhelmingly overexpressed, in 100% of cases 237where data had exposure duration metadata available (Fig. 3C). marA is a highly studied regulator, 238principally recognized as the key driver of the “multiple antibiotic resistance” phenotype. marA and the 239hundreds of genes in the mar regulon have been implicated in tolerance to a multitude of compounds 240beyond antibiotics30. Differences in the dynamics of expression for various genes could serve as a 241predictive tool, enabling researchers to identify genes that are consistently differentially expressed as part 242of a crucial response pathway versus those that are differentially expressed in more of a compensatory or 243reactionary manner, as has been observed31,32. 244

We sought to identify functional groups of genes that play a role in providing E. coli resistance

245and tolerance. We used the database and Panther’s statistical overrepresentation tool with Bonferroni 246correction33 to identify overrepresented biological process-associated gene ontologies (GO) in gene 10 ACS Paragon Plus Environment

Page 11 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

247expression profiles that distinguish resistant and sensitive strains (Fig. 3D). There were thirteen GO that 248were identified as over-represented in each of the four sets of genes that we examined: either 249overexpressed in resistant strains, overexpressed in sensitive strains, underexpressed in resistant strains, 250or underexpressed in sensitive strains. Each of these 13 GO can be considered generic “parent” categories 251that encompass many children, such as cellular process (GO:0009987), metabolic process (GO:0008152), 252or transport (GO:0006810). 253

More specific GO arise when searching for ontologies that are overrepresented across criteria or

254only in one set of genes. There were 17 GO terms found enriched when comparing genes overexpressed 255in resistant strains and genes underexpressed in resistant strains. In resistant strains, cellular 256macromolecule processes were enriched 1.2-fold (P = 0.035 for overexpressed, P = 0.036 for 257underexpressed genes). These macromolecule processes can include DNA replication and repair, RNA 258processes, and cell wall metabolism. As the GO nucleobase-containing compound metabolic process was 259also enriched across all resistant-associated datasets (1.2-fold with P = 0.020 in overexpressed genes, and 2601.3-fold with P = 4.4e-3 in underexpressed genes), this data suggests that resistant strains devote 261comparatively more resources to DNA/RNA maintenance than do sensitive strains. All ontologies 262overrepresented solely in underexpressed genes from resistant strains were associated with metabolic or 263biosynthetic processes, from alpha amino acid associated processes (1.7-fold, P = 0.011) to cellular 264nitrogen compounds (1.4-fold, P = 6.4e-6) and organic cyclic compounds (1.4-fold, P = 6.6e-5). These 265changes could imply that resistant strains conserve energy by decreasing expression of various metabolic 266processes. 267

Overrepresented ontologies in sensitive strains were generally not related to specific metabolic or

268biosynthetic processes. Many genes underexpressed in sensitive strains were associated with energy 269generation; anaerobic respiration (2.1-fold, P = 2.9e-3), aerobic respiration (2.4-fold, P = 5.3e-3), and 270electron transport chain (2.7-fold, P = 0.030) ontologies were overrepresented in this set of genes. 271Transport-related ontologies were over-represented in multiple sets of genes. For instance, ion transport 11 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 35

272was the only ontology overlap found between underexpressed genes in resistant strains and overexpressed 273genes in sensitive strains. Transport-associated ontologies (organic substance transport and single274organism transport) were overrepresented in resistant overexpression, sensitive overexpression, and 275sensitive underexpression gene sets. 276

Overall, this analysis reveals that gene expression changes promoting transport and altered energy

277metabolism predominate in sensitive strains, while metabolic and biosynthetic pathways modifications 278underlie the response in resistant strains. 279Enriched interactions in regulatory networks 280

Locating key regulatory players is another important aspect of decoding the non-genetic basis for

281tolerance and resistance phenotypes. To understand which regulators may be having more influence under 282any individual stress condition, we obtained sigma factor and transcription factor regulatory interaction 283files for E. coli K-12 from RegulonDB34. 284

For each stress condition, the number of differentially expressed genes known to interact with

285either a sigma factor or transcription factor was tallied to obtain the “actual” count of interactions. An 286expected probability of interactions for each regulatory factor was then obtained by dividing the number 287of genes a factor interacts with by the total number of gene interactions noted for all regulators. 288Significance of enrichment was assessed with a binomial test, as in Panther’s overrepresentation tool33. A 289Bonferroni correction was applied to all P-values reported here. We computed the fold enrichment in 290interactions by dividing the actual count by the expected count, then used hierarchical clustering to sort 291fold enrichment by the stress conditions and regulators and gain perspective on the sigma factors and 292transcription factors that are most likely to participate in a phenotype (Fig. 4). The complete datasets are 293provided in the Supporting Data File. 294

Sigma factors are proteins necessary to initiate transcription. The seven sigma factors in E. coli

295each enable RNA polymerase to bind with greater specificity to different promoter sequences, allowing 12 ACS Paragon Plus Environment

Page 13 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

296cells to broadly shift gene expression by altering sigma factor expression levels. The primary sigma factor 297RpoD is the most abundant sigma factor during normal exponential growth and interacts most strongly 298with the core RNA polymerase35. When cells are stressed, RpoD can be replaced by alternate sigma 299factors, including RpoS, RpoH, FliA, RpoE, RpoN, and FecI36. RpoS is known to be involved in the 300general stress response, starvation, and stationary phase. RpoH is the heat-shock sigma factor, while FliA 301controls hundreds of flagellar synthesis and motility-associated genes. RpoE is associated with membrane 302stress. RpoN controls nitrogen-related genes, and FecI is involved with iron transport. 303

In Fig. 4A, we include the fold enrichment in interactions by stress condition for seven sigma

304factors. Interestingly, when clustering the sigma factors, the motility factor FliA is placed farthest away 305from the other factors. The heatmap row corresponding to FliA shows that, for a large set of phenotypes, 306genes induced by FliA are overrepresented. In other words, more of the differentially expressed genes are 307regulated by FliA than would be expected, considering the number of genes FliA is known to regulate and 308the number of differentially expressed genes in each of these conditions. The phenotypes for which FliA 309was most highly enriched include glycerol (11.0-fold increase, P = 2.8e-16), gold (8.4-fold increase, P = 3107.7e-20), and geraniol (6.8-fold increase, P = 5.3e-20). Genes regulated by FliA were underrepresented 311most often in conditions such as ethanol (8.0-fold decrease, P = 4.1e-33) and heat (12.2-fold decrease, P 312= 7.6e-4). 313

In some stress conditions, we can compare the profiles from sensitive and resistant strains. For

314example, we have differential expression profiles from strains sensitive to ampicillin (labeled as 315ampicillin (S)) as well as resistant (labeled ampicillin). The sensitive ampicillin profile clusters with other 316sensitive profiles, including to the antibiotics kanamycin and rifampicin. These profiles all show 317enrichment of genes regulated by the heat shock factor RpoH (1.9-fold with P = 3.7e-7 for ampicillin, 3183.3-fold with P = 9.8e-37 for kanamycin, and 2.5-fold with P = 1.0e-9 for rifampicin). In contrast, the 319ampicillin-resistant profile is clustered most closely with profiles from cinnamaldehyde and geraniol 320exposure, which all have overrepresentation of genes regulated by FliA and RpoN. 13 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

321

Page 14 of 35

The most highly overrepresented sigma factor/phenotype pair was the ferric citrate transport

322sigma factor FecI in profiles from strains sensitive to the antibacterial compound 4,5-Dihydroxy-2323Cyclopenten-1-One (DHCP) (25.5-fold overrepresentation, P = 1.3e-5). The involvement of genes 324regulated by FecI was not noted in the original publication discussing these datasets37. This finding 325demonstrates the power of accumulating transcriptome datasets for re-evaluation in bringing to light new 326conclusions independent from interpretations in the available literature. 327

We performed similar inquiry to detect any overrepresented interactions with 203 different

328transcription factors. We show a limited set of the 35 transcription factors represented across the largest 329number of conditions in Fig. 4B. Here, we note that stresses belonging to the same category (e.g. 330solvent/biofuel, antibiotic, etc.) are often clustered together. For instance, similar pairs like organic acids 331succinic and itaconic acid, or the hormones epinephrine and norepinephrine, or the antibiotics ampicillin 332and kanamycin, are each clustered together. 333

Uncovering important transcription factors, and therefore promising bioengineering or therapeutic

334targets, can be achieved by focusing on particularly interesting regions. In the columns corresponding to 335acid profiles (either acid sensitive or acid tolerant), we note one area of highly overrepresented 336transcription factor interactions, which includes the transcription factors GadX, GadE-RcsB, GadW, 337AdiY, TorR, and FliZ. For acid gene expression profiles, genes controlled by these regulators are 338overrepresented between 7.3-fold and 21.2-fold, with P < 0.05 for all cases. Many of these transcription 339factors are expected; for instance, GadX and GadW are part of the same operon, known as the acid fitness 340island, which is regulated by a central activator of E. coli acid response, GadE-RcsB38,39. AdiY has also 341been noted to be induced under low pH40. However, the regulators TorR and FliZ have not been 342recognized as playing a role in response to acid treatment. Thus, this analysis allows for a deeper 343understanding of regulatory interactions associated with a complex stress condition. 344

The most overrepresented transcription factor/phenotype pair was MarR in strains exposed to

345romaine lettuce lystate, at 42.4-fold enrichment, P = 1.2e-4. MarR is a repressor for the mar operon. This 14 ACS Paragon Plus Environment

Page 15 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

346transcription factor was also found to be significantly enriched in E. coli treated with the arylamide 347foldamer PMX (27.8-fold, P = 0.038), tetracycline (19.8-fold, P = 1.6e-4), bleach (29.0-fold, P = 1.8e-5), 348or cis-1,2-dicholorethylene (18.0-fold, P = 0.017). 349Combining genome and transcriptome data to ascertain targets not yet mutationally exploited 350

A significant advantage of this database is the availability of both mutation-level and

351transcriptome-level data. With this information, it is possible to perform queries that were previously very 352labor-intensive, requiring manual search of available literature. We here provide the results of one such 353search: inquiring into genes that are commonly differentially expressed in many stress conditions that 354have not been explored as metabolic engineering targets (i.e. there are no mutations for these genes in the 355database). Table 3 lists the ten genes that were most often differentially expressed but not mutated in any 356instance in the database. Nine out of these ten genes are nonessential per the PEC database 357(https://shigen.nig.ac.jp/ecoli/pec/); rpmC, encoding a 50S ribosomal subunit protein, is the only essential 358gene. 359

All non-mutated but often differentially expressed genes shown in Table 3 are most likely

360participants in more universal stress responses, as they are observed to respond to a variety of dissimilar 361conditions. In some cases, we can speculate which phenotypes might be impacted by gene mutation by 362examining the stress conditions in which that gene often appears differentially expressed. For instance, 363the ribosome stability factor raiA was found to be differentially expressed upon treatment with four 364different metal ions cadmium, cobalt, nickel, and zinc. RaiA, also known as YfiA, prevents translation by 365blocking tRNA binding41. It was underexpressed in the presence of the metal ions and paraquat and 366overexpressed in the other conditions shown in Table 3, except for ethanol, in which it was 367underexpressed after exposure times greater than 700 hours and overexpressed after shorter times. pspC 368was differentially expressed in the presence of the disparate antibiotics ampicillin, fluoroquinolones, 369kanamycin, and norfloxacin. PspC, named phage shock protein for its role in bacteriophage response, is 370predicted to be part of a toxin-antitoxin pair with PspB42. The direction of differential expression change 15 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 35

371was inconsistent between conditions. It was overexpressed in kanamycin, fluoroquinolone, MMS, 372isobutanol, n-hexane, romaine lettuce and underexpressed in ampicillin, norfloxacin, UTI, and ethanol 373stress. The mechanism for involvement in stress response is unclear, but there is some evidence that PspC 374influences growth rate42. 375Estimating promiscuity in differential expression 376

When utilizing gene expression datasets to select targets for a specific bioengineering problem, it

377is possible to consider whether a promising gene candidate is responding exclusively to the specific stress 378condition being studied, or if the gene is likely to respond under unrelated conditions. Considering the 379promiscuity of gene expression is one way to distinguish genes that may be participating in a global, 380nonspecific stress response versus those that respond in a more targeted manner. Higher promiscuity 381implies that a gene was differentially expressed across many dissimilar conditions, whereas lower 382promiscuity implies that the gene was differentially expressed only under a few, similar conditions. The 383degree of promiscuity in gene expression can be estimated by searching across transcriptome datasets 384from unrelated stress conditions. Depending on the application, high promiscuity might be a trait that is 385either acceptable or undesirable. To estimate promiscuity, we counted the number of unique conditions in 386which each gene in the database was differentially expressed (a complete list is provided in the 387Supporting Data File). 388

As expected, we find that genes very commonly differentially expressed tend toward high

389promiscuity (Fig. 5A). Indeed, the number of times a gene is differentially expressed is positively 390correlated with the number of conditions in which it appears (Pearson correlation coefficient = 0.89, P < 3911e-100). A distribution of the number of conditions in which each gene was differentially expressed is 392shown in Fig. 5B. The median is six conditions and the maximum is 33 conditions. 393

The two genes differentially expressed across the largest number of distinct conditions are marA

394and ompC (Fig. 5C). marA was differentially expressed a total of 58 times in 33 conditions. It was

16 ACS Paragon Plus Environment

Page 17 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

395overexpressed in nearly all instances (93%). The ability of the mar regulon to broadly convey tolerance 396renders marA an interesting target for drug-resistance applications, but perhaps not as suitable of a target 397for long-term deployment of robust production strains, as it has been suggested that activation of the mar 398phenotype is costly43. ompC was differentially expressed 42 times in 31 conditions, and was 399underexpressed in a majority of cases (60%). ompC encodes an outer-membrane porin, recognized as 400contributing to stress response to a variety of environmental conditions44. Estimations of promiscuity can 401be used to avoid the selection of highly general stress-associated genes for applications that require a 402more targeted approach to engineer tolerance. 403

We located genes that were less promiscuous: differentially expressed in response to a limited set

404of stress categories. To identify genes responding most specifically to antibiotic stress, we considered 405only those genes differentially expressed in at least 10 datasets. We found yjhP to be 100% specific to 406antibiotic stress, being differentially expressed only in ampicillin, glyphosate, kanamycin, norfloxacin, 407and rifampicin (Fig. 5D). yjhP was overexpressed in 64% of cases, but the direction of expression was not 408specific to the stress; it was overexpressed in one ampicillin case but underexpressed in another. The 409function of this gene has not been validated. Sequence-based analysis classifies yjhP as a putative 410methyltransferase in the AdoMet_MTases superfamily, which use S-adenosyl L-methionine as a substrate 411and produce S-adenosyl-L-homocysteine45. While methylation and epigenetic regulation is well-studied in 412eukaryotes, methyltransferases have not yet been implicated in bacterial antibiotic tolerance. Recently, it 413was shown that a 5-methylcytosine methyltransferase does influence membrane stress response in Vibrio 414cholerae46, so it may be that yjhP plays a regulatory role via epigenetic modification. 415

When searching for genes that responded most specifically to solvent/biofuel stress, we limited

416the search to genes that were differentially expressed in at least 20 datasets. rpoZ and hldD were two 417genes that were highly specific to biofuel/solvent stress (Fig. 5E). Out of the 32 datasets in which rpoZ 418was differentially expressed, 94% were associated with either ethanol or isobutanol stress. This gene 419encodes the ω-subunit of RNA polymerase, and could influence global gene expression by modifying 17 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 35

420sigma 38 and sigma 70 factor binding to the core polymerase47. An rpoZ deletion mutant was found to 421have increased influence from sigma 3847. This agrees with our sigma factor overrepresentation results, as 422rpoZ was underexpressed in ethanol and isobutanol and genes regulated by sigma 38 were 423overrepresented in each of these conditions. hldD was differentially expressed in 25 instances, 93% of 424which were ethanol or isobutanol. This gene encodes an enzyme in the ADP-L-β-D-heptose biosynthesis 425pathway, which produces an important component of lipopolysaccharides48. hldD was underexpressed in 426ethanol and UTI, but overexpressed in isobutanol. 427

The least promiscuous genes would be those differentially expressed only in one stress condition.

428In Fig. 5E, we show four genes (asnU, asnV, asnT, and asnW) each differentially expressed in more than 42928 datasets, all of which correspond to ethanol treatment. These genes encode the four asparagine tRNAs, 430which all share the same 76-bp sequence in E. coli K12 MG1655. They were overwhelmingly 431overexpressed, with asnV and asnW were overexpressed in 100% of instances, while asnU was 432overexpressed in 97% of cases and asnT in 96% of cases. None of the asparagine tRNAs were mutated in 433any record, suggesting the potential to exploit these genes in novel ways for improving ethanol tolerance. 434While there is no literature demonstrating that asparagine usage directly improves ethanol tolerance in E. 435coli, researchers have noted that asparagine pretreatment potentiates the effects of alcohol consumption in 436rats49. Others studies in rats treated with amino acids have also impacted the effects of alcohol, leading to 437the hypothesis that amino acids might form a complex with acetaldehydes50. 438Conclusions 439

We here perform an analysis of E. coli tolerance and resistance that encompasses information

440from more than 170 gene expression datasets across 89 different stress conditions. Our results show that 441the accumulation of diverse transcriptome datasets can both validate previous discoveries, such as the 442importance of acid-related transcription factors in acid stress, or the identification of common stress443response regulators like marA. We were also able to gain a wealth of new information on E. coli 444resistance and tolerance. The most common differentially expressed genes were highlighted, along with 18 ACS Paragon Plus Environment

Page 19 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

445genes that were most commonly co-expressed, across a broad set of conditions. This provides a snapshot 446of a universal stress response, which could inform studies as wide-ranging as the evolution of antibiotic 447resistance, tolerance to desired biochemical products, and food safety. We scrutinize how resistant and 448sensitive profiles differ in the functional classes impacted, which could lead to new methods by which to 449engineer tolerance. We present promising gene and regulator targets that are likely to be players in a 450universal, general stress response (Fig. 2, Table 2, Table 3), as well as targets more likely to contribute to 451tolerance to specific stresses (in Fig. 4, Fig. 5D-E). 452

This work results in the first resource for E. coli that facilitates search and analysis of both

453genomic and transcriptomic-level information, along with metadata relevant to allow inquiry specific to 454myriad resistance phenotypes and experimental conditions. Data curation will be an ongoing process, and 455this resource will only grow in power as the number of included studies with transcriptome information is 456increased. Future efforts will continue to expand the types of biological data considered, including 457planned incorporation of proteomic and metabolomic data coupled with new tools to estimate the effect of 458protein mutations, and eventually design predictions for de novo engineering of resistant strains. Further 459analysis that focuses on any one of a conceivable set of conditions of interest will enable researchers to 460identify biological processes that promote tolerance and better understand the subtle discrepancies in non461genetic responses that confirm either a resistant or sensitive phenotype. We have made the data and 462software needed to leverage this resource publicly available, with the anticipation that the community will 463take advantage of this database for their own creative and unique research questions on resistance and 464tolerance. 465Supporting Information 466

The Supporting Text provides information and examples of how to load and access the database. The

467Supporting Data File contains additional information from which results presented in the main text are 468derived. This file contains datasets including:

19 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

469



470 471







480

Pairwise Counts: Number of instances in which each pair of genes was observed to be differentially expressed in the same strain



478 479

Promiscuity by condition: All genes and the conditions in which they were observed to be differentially expressed

476 477

Transcription Factor Overrepresentation: Fold enrichment in genes regulated by indicated transcription factor, along with adjusted P-values

474 475

Sigma Factor Overrepresentation: Fold enrichment in genes regulated by indicated sigma factor, along with adjusted P-values

472 473

Page 20 of 35

Pairwise Counts-Overexpression: Number of instances in which each pair of genes was observed to be differentially overexpressed in the same strain



Pairwise Counts-Underexpression: Number of instances in which each pair of genes was observed to be differentially underexpressed in the same strain

481Acknowledgements 482We are grateful to Prof. Hsueh-Fen Juan (National Taiwan University) and her group for hosting K.E.E. 483during the development of this project. This work was supported by the National Science Foundation 484(#1613746 to K.E.E), the Department of Energy Genome Sciences Project (#DE-SC008812 to R.T.G), 485and the National Science Foundation (MCB-1714564), DARPA Young Faculty Award (D17AP00024) 486and William Keck Foundation (#OCG6145B) to A.C. 487Author Contribution 488K.E.E., J.D.W., and A.C. wrote the manuscript. K.E.E. curated data related to gene expression and 489modified software for incorporation of gene expression data. J.D.W. wrote all the original software and 490curated data associated with mutation. D.T.N. assisted with compatibility between the mutation and gene491expression associated records. A.C. and R.T.G. supervised the project. 492Author Information

20 ACS Paragon Plus Environment

Page 21 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

493Corresponding author 494 *Email: [email protected] 495Present addresses 496‡(KEE) Theoretical Division, Los Alamos National Laboratory, PO Box 1663, Los Alamos, NM, USA 497§(JDW) Shell Biodomain, 3333 Texas 6, Houston, TX, USA 498Notes: The authors declare no competing financial interest. 499References 500(1) Sandoval, N. R., Kim, J. Y. H., Glebes, T. Y., Reeder, P. J., Aucoin, H. R., Warner, J. R., and Gill, R. 501T. (2012) Strategy for directing combinatorial genome engineering in Escherichia coli. Proc. Natl. Acad. 502Sci. 109, 10540–10545. 503(2) Warner, J. R., Reeder, P. J., Karimpour-Fard, A., Woodruff, L. B., and Gill, R. T. (2010) Rapid 504profiling of a microbial genome using mixtures of barcoded oligonucleotides. Nat. Biotechnol. 28, 856– 505862. 506(3) Nicolaou, S. A., Gaida, S. M., and Papoutsakis, E. T. (2010) A comparative view of metabolite and 507substrate stress and tolerance in microbial bioprocessing: From biofuels and chemicals, to biocatalysis 508and bioremediation. Metab. Eng. 12, 307–331. 509(4) Gibson, M. K., Forsberg, K. J., and Dantas, G. (2014) Improved annotation of antibiotic resistance 510determinants reveals microbial resistomes cluster by ecology. ISME J. 9, 1–10. 511(5) Gupta, S. K., Padmanabhan, B. R., Diene, S. M., Lopez-Rojas, R., Kempf, M., Landraud, L., and 512Rolain, J. M. (2014) ARG-annot, a new bioinformatic tool to discover antibiotic resistance genes in 513bacterial genomes. Antimicrob. Agents Chemother. 58, 212–220. 514(6) McArthur, A. G., Waglechner, N., Nizam, F., Yan, A., Azad, M. A., Baylay, A. J., Bhullar, K.,

21 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 35

515Canova, M. J., De Pascale, G., Ejim, L., Kalan, L., King, A. M., Koteva, K., Morar, M., Mulvey, M. R., 516O’Brien, J. S., Pawlowski, A. C., Piddock, L. J. V, Spanogiannopoulos, P., Sutherland, A. D., Tang, I., 517Taylor, P. L., Thaker, M., Wang, W., Yan, M., Yu, T., and Wright, G. D. (2013) The comprehensive 518antibiotic resistance database. Antimicrob. Agents Chemother. 57, 3348–3357. 519(7) Zankari, E., Hasman, H., Cosentino, S., Vestergaard, M., Rasmussen, S., Lund, O., Aarestrup, F. M., 520and Larsen, M. V. (2012) Identification of acquired antimicrobial resistance genes. J. Antimicrob. 521Chemother. 67, 2640–2644. 522(8) Rowe, W., Baker, K. S., Verner-Jeffreys, D., Baker-Austin, C., Ryan, J. J., Maskell, D., and Pearce, 523G. (2015) Search engine for antimicrobial resistance: A cloud compatible pipeline and web interface for 524rapidly detecting antimicrobial resistance genes directly from sequence data. PLoS One 10. 525(9) Jacoby, G. A. (2005) Mechanisms of resistance to quinolones. Clin. Infect. Dis. 41, S120-6. 526(10) Ezekiel, D. H., and Hutchins, J. E. (1968) Mutations affecting RNA Polymerase associated with 527Rifampicin Resistance in Escherichia coli. Nature 220, 276–277. 528(11) Leinonen, R., Sugawara, H., and Shumway, M. (2010) The sequence read archive. Nucleic Acids 529Res. gkp1019. 530(12) Edgar, R., Domrachev, M., and Lash, A. E. (2002) Gene Expression Omnibus: NCBI gene 531expression and hybridization array data repository. Nucleic Acids Res 30, 207–210. 532(13) Faith, J. J., Driscoll, M. E., Fusaro, V. A., Cosgrove, E. J., Hayete, B., Juhn, F. S., Schneider, S. J., 533and Gardner, T. S. (2008) Many Microbe Microarrays Database: Uniformly normalized Affymetrix 534compendia with structured experimental metadata. Nucleic Acids Res. 36, 866–870. 535(14) López-Maury, L., Marguerat, S., and Bähler, J. (2008) Tuning gene expression to changing 536environments: from rapid responses to evolutionary adaptation. Nat. Rev. Genet. 9, 583–593. 537(15) Timsah, Z., Ahmed, Z., Ivan, C., Berrout, J., Gagea, M., Zhou, Y., Pena, G. N. A., Hu, X., Vallien, 22 ACS Paragon Plus Environment

Page 23 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

538C., Kingsley, C. V, Lu, Y., and Hancock, J. F. (2015) Grb2 depletion under non-stimulated conditions 539inhibits PTEN , promotes Akt-induced tumor formation and contributes to poor prognosis in ovarian 540cancer. Oncogene 1–11. 541(16) Ben-Jacob, E., Lu, M., Schultz, D., and Onuchic, J. N. (2014) The physics of bacterial decision 542making. Front. Cell. Infect. Microbiol. 4, 154. 543(17) Solopova, A., van Gestel, J., Weissing, F. J., Bachmann, H., Teusink, B., Kok, J., and Kuipers, O. P. 544(2014) Bet-hedging during bacterial diauxic shift. Proc. Natl. Acad. Sci. U. S. A. 111, 7427–32. 545(18) Garcia-Bernardo, J., and Dunlop, M. J. (2015) Noise and low-level dynamics can coordinate 546multicomponent bet hedging mechanisms. Biophys. J. 108, 184–193. 547(19) Veening, J.-W., Smits, W. K., and Kuipers, O. P. (2008) Bistability, epigenetics, and bet-hedging in 548bacteria. Annu. Rev. Microbiol. 62, 193–210. 549(20) Wood, T. K., Knabel, S. J., and Kwan, B. W. (2013) Bacterial persister cell formation and dormancy. 550Appl. Environ. Microbiol. 79, 7116–7121. 551(21) Fernández, L., and Hancock, R. E. W. (2012) Adaptive and mutational resistance: role of porins and 552efflux pumps in drug resistance. Clin. Microbiol. Rev. 25, 661–81. 553(22) Otoupal, P. B., Erickson, K. E., Bordoy, A. E., and Chatterjee, A. (2016) CRISPR perturbation of 554gene expression alters bacterial fitness under stress and reveals underlying epistatic constraints. ACS 555Synth. Biol. 6, 94–107. 556(23) Lynch, M. D., Warnecke, T., and Gill, R. T. (2007) SCALEs: multiscale analysis of library 557enrichment. Nat. Methods 4, 87–93. 558(24) Freed, E. F., Winkler, J. D., Weiss, S. J., Garst, A. D., Mutalik, V. K., Arkin, A. P., and Gill, R. T. 559(2015) Genome-Wide Tuning of Protein Expression Levels to Rapidly Engineer Microbial Traits. ACS 560Synth. Biol. 4, 1244–1253. 23 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 35

561(25) Erickson, K. E., Otoupal, P. B., and Chatterjee, A. (2017) Transcriptome-level signatures in gene 562expression and gene expression variability during bacterial adaptive evolution. mSphere 2, 1–17. 563(26) Winkler, J. D., Halweg-Edwards, A. L., Erickson, K. E., Choudhury, A., Pines, G., and Gill, R. T. 564(2016) The Resistome: A Comprehensive Database of Escherichia coli Resistance Phenotypes. ACS 565Synth. Biol. 566(27) Krzywinski, M. et al. (2009) Circos: an Information Aesthetic for Comparative Genomics. Genome 567Res 19, 1639–1645. 568(28) Martinez, A., and Kolter, R. (1997) Protection of DNA during oxidative stress by the nonspecific 569DNA-binding protein Dps. J. Bacteriol. 179, 5188–5194. 570(29) Nair, S., and Finkel, S. E. (2004) Dps Protects Cells against Multiple Stresses during Stationary 571Phase Dps Protects Cells against Multiple Stresses during Stationary Phase. J Bacteriol 186, 4192–8. 572(30) Alekshun, M. N., and Levy, S. B. (1999) The mar regulon: multiple resistance to antibiotics and 573other toxic chemicals. Trends Microbiol. 7, 410–413. 574(31) Fong, S. S., Joyce, A. R., and Palsson, B. Ø. (2005) Parallel adaptive evolution cultures of 575Escherichia coli lead to convergent growth phenotypes with different gene expression states. Genome 576Res. 15, 1365–1372. 577(32) Price, M. N., Deutschbauer, A. M., Skerker, J. M., Wetmore, K. M., Ruths, T., Mar, J. S., Kuehl, J. 578V, Shao, W., and Arkin, A. P. (2013) Indirect and suboptimal control of gene expression is widespread in 579bacteria. Mol. Syst. Biol. 9, 660. 580(33) Mi, H., Muruganujan, A., Casagrande, J. T., and Thomas, P. D. (2013) Large-scale gene function 581analysis with the PANTHER classification system. Nat. Protoc. 8, 1551–1566. 582(34) Salgado, H., Peralta-Gil, M., Gama-Castro, S., Santos-Zavaleta, A., Muñiz-Rascado, L., García583Sotelo, J. S., Weiss, V., Solano-Lira, H., Martínez-Flores, I., and Medina-Rivera, A. (2013) RegulonDB 24 ACS Paragon Plus Environment

Page 25 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

584v8.0: omics data sets, evolutionary conservation, regulatory phrases, cross-validated gold standards and 585more. Nucleic Acids Res. 41, D203–D213. 586(35) Maeda, H., Fujita, N., and Ishihama, A. (2000) Competition among seven Escherichia coli sigma 587subunits: relative binding affinities to the core RNA polymerase. Nucleic Acids Res. 28, 3497–3503. 588(36) Österberg, S., del Peso-Santos, T., and Shingler, V. (2011) Regulation of alternative sigma factor 589use. Annu Rev Microbiol 65, 37–55. 590(37) Phadtare, S., Kato, I., and Inouye, M. (2002) DNA Microarray Analysis of the Expression Profile of 591Escherichia coli in Response to Treatment with 4 , 5-Dihydroxy-2-Cyclopenten-1-One. J. Bacteriol. 184, 5926725–6729. 593(38) Hommais, F., Krin, E., Coppee, J. Y., Lacroix, C., Yeramian, E., Danchin, A., and Bertin, P. (2004) 594GadE (YhiE): A novel activator involved in the response to acid environment in Escherichia coli. 595Microbiology 150, 61–72. 596(39) Ma, Z., Richard, H., Tucker, D. L., Conway, T., and Foster, J. W. (2002) Collaborative Regulation of 597Escherichia coli Glutamate-Dependent Acid Resistance by Two AraC-Like Regulators , GadX and GadW 598Collaborative Regulation of Escherichia coli Glutamate-Dependent Acid Resistance by Two AraC-Like 599Regulators , GadX and GadW (YhiW). J. Bacteriol. 184(24):70, 1123–1233. 600(40) Stim-Herndon, K. P., Flores, T. M., and Bennett, G. N. (1996) Molecular characterization of adiY, a 601regulatory gene which affects expression of the biodegradative acid-induced arginine decarboxylase gene 602(adiA) of Escherichia coli. Microbiology 1311–1320. 603(41) Vila-Sanjurjo, A., Schuwirth, B., Hau, C. W., and Cate, J. H. D. (2004) Structural basis for the 604control of translation initiation during stress. Nat. Struct. Mol. Biol. 11, 1054–1060. 605(42) Brown, J. M., and Shaw, K. J. (2003) A Novel Family of Escherichia coli Toxin-Antitoxin Gene 606Pairs. J. Bacteriol. 185, 6600–6608.

25 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 35

607(43) Garcia-Bernardo, J., and Dunlop, M. J. (2013) Tunable stochastic pulsing in the Escherichia coli 608multiple antibiotic resistance network from interlinked positive and negative feedback loops. PLoS 609Comput. Biol. 9, e1003229. 610(44) Pratt, L. A., Hsing, W., Gibson, K. E., and Silhavy, T. J. (1996) From acids to osmZ: Multiple 611factors influence synthesis of the OmpF and OmpC porins in Escherichia coli. Mol. Microbiol. 20, 911– 612917. 613(45) Marchler-Bauer, A., Derbyshire, M. K., Gonzales, N. R., Lu, S., Chitsaz, F., Geer, L. Y., Geer, R. C., 614He, J., Gwadz, M., Hurwitz, D. I., Lanczycki, C. J., Lu, F., Marchler, G. H., Song, J. S., Thanki, N., 615Wang, Z., Yamashita, R. A., Zhang, D., Zheng, C., and Bryant, S. H. (2015) CDD: NCBI’s conserved 616domain database. Nucleic Acids Res. 43, D222–D226. 617(46) Chao, M. C., Zhu, S., Kimura, S., Davis, B. M., Schadt, E. E., Fang, G., and Waldor, M. K. (2015) A 618Cytosine Methytransferase Modulates the Cell Envelope Stress Response in the Cholera Pathogen. PLoS 619Genet. 11, 1–24. 620(47) Geertz, M., Travers, A., and Mehandziska, S. (2011) Supercoiling in Coordinating Structural 621Coupling between RNA Polymerase Composition and DNA Supercoiling in Coordinating Transcription: 622a Global Role for the Omega Subunit? MBio 2, 1–11. 623(48) Kneidinger, B., Marolda, C., Graninger, M., Zamyatina, A., McArthur, F., Kosma, P., Valvano, M. 624A., and Messner, P. (2002) Biosynthesis pathway of ADP-L-glycero-beta-D-manno-heptose in 625Escherichia coli. J. Bacteriol. 184, 363–9. 626(49) Forney, R. B., Hughes, F. W., Richards, A. B., and Gates, P. W. (1963) Toxicity and depressant 627action of ethanol and hexobarbital after pretreatment with asparagine. Toxicol. Appl. Pharmacol. 5, 790– 628793. 629(50) Breglia, R. J., Ward, C. O., and Jarowski, C. I. (1973) Effect of selected amino acids on ethanol

26 ACS Paragon Plus Environment

Page 27 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

630toxicity in rats. J. Pharm. Sci. 62, 49–55. 631(51) Keseler, I. M., Collado-Vides, J., Santos-Zavaleta, A., Peralta-Gil, M., Gama-Castro, S., Muñiz632Rascado, L., Bonavides-Martinez, C., Paley, S., Krummenacker, M., and Altman, T. (2011) EcoCyc: a 633comprehensive database of Escherichia coli biology. Nucleic Acids Res. 39, D583–D590. 634

27 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 35

635 636Figure 1. Structure and contents of the database. (A) Hierarchy of classes in the database. A summary of 637attributes in each class are shown. (B) The proportion of gene expression changes, by year, that were measured with 638either microarray, RNAseq, or qPCR. (C) The percentage of gene expression changes associated with various stress 639categories. (D) The thirty-four E. coli strains with gene expression data, and the proportion of stress categories for 640which that gene expression data corresponds. Strains are grouped by their common derivative: either the laboratory 641strains B, K-12, or W, or a pathogenic strain. 642

28 ACS Paragon Plus Environment

Page 29 of 35

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Synthetic Biology

643 644Figure 2. Frequency and pairwise co-occurrence of differential gene expression for all gene expression 645datasets. Genes are plotted relative to their position on the E. coli MG1655 chromosome; genes not present in 646MG1655 are not included. The spike height in the outer radius is proportional to the frequency at which the 647corresponding gene was found to be differentially expressed. The inner radius links genes that were found to be 648differentially expressed in the same strain. The top 2080 pairwise linkages are displayed here. Plot was generated 649using Circos. 650

29 ACS Paragon Plus Environment

ACS Synthetic Biology

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 35

651 652Figure 3. Insight into dynamics of differentially expressed genes and functional classes impacted by stress. 653The number of instances in which the genes (A) dps, (B) manX, and (C) marA were differentially expressed, by the 654amount of exposure time used in the experiment. (D) Overrepresented biological process gene ontologies (GO) from 655select gene expression datasets. Overrepresented GO (P