Systematic Discovery of Archaeal Transcription Factor Functions in Regulatory Networks through Quantitative Phenotyping Analysis

Abstract

Gene regulatory networks (GRNs) are critical for dynamic transcriptional responses to environmental stress. However, the mechanisms by which GRN regulation adjusts physiology to enable stress survival remain unclear. Here we investigate the functions of transcription factors (TFs) within the global GRN of the stress-tolerant archaeal microorganism Halobacterium salinarum. We measured growth phenotypes of a panel of TF deletion mutants in high temporal resolution under heat shock, oxidative stress, and low-salinity conditions. To quantitate the noncanonical functional forms of the growth trajectories observed for these mutants, we developed a novel modeling framework based on Gaussian process regression and functional analysis of variance (FANOVA). We employ unique statistical tests to determine the significance of differential growth relative to the growth of the control strain. This analysis recapitulated known TF functions, revealed novel functions, and identified surprising secondary functions for characterized TFs. Strikingly, we observed that the majority of the TFs studied were required for growth under multiple stress conditions, pinpointing regulatory connections between the conditions tested. Correlations between quantitative phenotype trajectories of mutants are predictive of TF-TF connections within the GRN. These phenotypes are strongly concordant with predictions from statistical GRN models inferred from gene expression data alone. With genome-wide and targeted data sets, we provide detailed functional validation of novel TFs required for extreme oxidative stress and heat shock survival. Together, results presented in this study suggest that many TFs function under multiple conditions, thereby revealing high interconnectivity within the GRN and identifying the specific TFs required for communication between networks responding to disparate stressors. IMPORTANCE To ensure survival in the face of stress, microorganisms employ inducible damage repair pathways regulated by extensive and complex gene networks. Many archaea, microorganisms of the third domain of life, persist under extremes of temperature, salinity, and pH and under other conditions. In order to understand the cause-effect relationships between the dynamic function of the stress network and ultimate physiological consequences, this study characterized the physiological role of nearly one-third of all regulatory proteins known as transcription factors (TFs) in an archaeal organism. Using a unique quantitative phenotyping approach, we discovered functions for many novel TFs and revealed important secondary functions for known TFs. Surprisingly, many TFs are required for resisting multiple stressors, suggesting cross-regulation of stress responses. Through extensive validation experiments, we map the physiological roles of these novel TFs in stress response back to their position in the regulatory network wiring. This study advances understanding of the mechanisms underlying how microorganisms resist extreme stress. Given the generality of the methods employed, we expect that this study will enable future studies on how regulatory networks adjust cellular physiology in a diversity of organisms.

Keywords: Archaea; functional ANOVA; phenomics; transcription factors.

FIG 1

TF candidate selection pipeline. Genes encoding proteins with a putative DNA binding domain were annotated using sequence databases PFAM (37) and PSI-BLAST (38), structural predictions (26), and protein functions from COG (39). These annotations were stored in the Systems Biology Experiment and Analysis Management System (SBEAMS) database (36), resulting in 130 putative TFs. Transcriptome analysis across 1,495 experimental conditions (43) and GRN network inference models (7, 15) were then used to generate predictions regarding TF functions. Details of GRN predictions and criteria for selection of the final collection of 27 TFs are given in Table S1 in the supplemental material. DBD, DNA binding domain.

FIG 2

FANOVA modeling and statistical ranking of TF knockout mutant growth phenotypes across five environmental conditions, standard growth conditions, paraquat stress, peroxide, low salt, and heat shock. (A) Mutants with the largest difference in growth compared to the Δura3 control strain under each condition. Raw data from growth trajectories for individual cultures under standard conditions (thin gray lines) were fit and compared with growth under stress conditions (thin blue lines) using FANOVA (see Materials and Methods). Solid lines indicate the mean of the fit to all replicate trajectories, and shaded regions are the 95% confidence interval of the fit. (B) Functional difference (OD_Δ) of the TF knockout strain relative to the isogenic Δura3 parent strain. The mutants with the top three scoring phenotypes according to $||{\text{OD}}_{\Delta }||$ are shown. (C) Summary statistical metric ranking ($||{\text{OD}}_{\Delta }||$) for those mutants with strongly different growth trajectories compared to that of the Δura3 strain ($||{\text{OD}}_{\Delta }||$≥ 0.337 across all conditions).

FIG 3

Phenotype network analysis reveals three major classes of TF mutants and extensive cross-regulation of stress responses by TFs. Node and edge attributes are given in the key at the bottom of the figure. OD_Δ numbers in the edge color legend refer to the maximum (fast growth) or minimum (slow growth) value of the mean of the posterior prediction from the FANOVA model across the entire growth time course. $||{\text{OD}}_{\Delta }||$numbers for edge thickness refer to the median value of the distribution from $||{\text{OD}}_{\Delta }||$ boxplots (see Fig. 2C, Fig. S3, and Materials and Methods).

FIG 4

TFs that regulate each other have similar OD_Δ phenotype trajectories. (A and B) Heat maps depict hierarchical clustering of OD_Δ trajectories for the 27 TF knockout mutants under paraquat (PQ) (A) and peroxide (B) conditions. TFs transcriptionally regulated by RosR are known (50) and indicated by red text. Colors in the dendrogram represent different clusters. The color scale indicates the mean of the posterior OD_Δ distribution for each mutant across the growth time course (x axis). (C) Each of the seven mutants regulated by RosR are statistically enriched for strongly correlated phenotypes ($\widehat{\rho }$≥ 0.4) with the ΔrosR mutant under PQ and standard growth conditions relative to other conditions (red text). Colored squares represent significant correlations (P ≤ 0.001), and white squares represent nonsignificant correlations. See the color scale to the right of the figure for correlation values.

FIG 5

Phenotype validation: a novel heat shock function for the copper-responsive regulator CopR. (A) Slow growth under heat shock conditions in the ΔcopR mutant (left graph) is complemented by expression of the copR gene in trans (right graph). (B) Cytoscape gene regulatory network depicting the significant overlap between genes regulated by CopR in response to copper overload (CopR node at the top) and in wild-type cells in response to heat shock (heat node at the bottom) (55) and copper (copper node on the left) (22). Node colors (representing expression levels in wild-type cells under heat shock) and edge line types are indicated in the keys. Gray boxes behind groups of nodes represent arCOG functional categories (61). P values of enrichment were calculated by the hypergeometric distribution test. Genes of unknown function are not shown here for clarity but are given in Fig. S5. (C) Quantitative RT-PCR gene expression of cctA and nirJ genes in the knockout strain compared to the Δura3 control strain. The levels of expression before heat shock (0′) and 30 min after induction of heat shock (30′) are shown. Expression is normalized relative to a control gene whose expression does not change during heat shock (see Materials and Methods). Error bars represent standard errors of the means (SEM) of three biological replicate cultures, each with three technical replicate trials.

FIG 6

Phenotype validation: a novel oxidative stress function for the cold shock family protein CspD1. (A) Slow growth under oxidative stress conditions in the ΔcspD1 mutant (left graph) is complemented by expression of the cspD1 gene in trans (right graph). (B) Line plot depicting the expression of the cspD1 gene (left axis) during fluctuations in oxygen concentrations (gray line, right axis). cspD1 expression is compared to that of the gene encoding known oxidative stress regulator RosR (left axis). (C) Line plot depicting the expression of genes requiring CspD1 for appropriate dynamic expression in response to oxygen. Each line represents the mean expression value of the groups of genes indicated in the legend. Expression data and annotations for individual genes are given in Table S4. WT, wild type. (D) Expression profiles for 106 of the 132 CspD1-dependent genes differentially expressed under PQ conditions in wild-type cells (15). Thin lines represent the expression of individual genes. Thick lines show the mean of induced or repressed genes. (E) Overlap between EGRIN predictions for CspD1 target gene regulatory influences (7) and differentially expressed genes (DEG) in the ΔcspD1 background. (F) Similar to panel E except that predictions under oxidative stress conditions are shown (15).