Modeling the dynamics of a 21‐gene, 50‐edge gene regulatory network controlling the transcriptional response to cold shock in Saccharomyces cerevisiae using GRNmap
Authors: Johnson K. Grace , Williams Natalie E. , Fitzpatrick Ben G. , Dahlquist Kam D. 
Published: June 2021

Gene expression is regulated by proteins called transcription factors which can either repress or activate a gene's transcriptional output. A gene regulatory network (GRN) consists of a set of transcription factors that regulate the level of expression of genes encoding other transcription factors. The dynamics of a GRN show how gene expression in the network changes over time. Previously in the lab, a MATLAB software package called GRNmap was developed that uses ordinary differential equations to model the dynamics of medium‐scale GRNs from budding yeast, Saccharomyces cerevisiae. The program estimates production rates, expression thresholds, and regulatory weights for each transcription factor in the network based on DNA microarray data, and then performs a forward simulation of the dynamics of the network. DNA microarray data for 6189 yeast genes was obtained from the Dahlquist lab where they subjected yeast to cold shock at 13°C and measured gene expression at three time points (after 15, 30, and 60 minutes of cold shock). We performed LOESS normalization using the limma package in R and used a modified ANOVA to determine which genes had a log2 fold change significantly different than zero at any of the timepoints studied. Using GRNmap, we estimated the parameters of a GRN derived from the YEASTRACT database, consisting of 21 nodes (transcription factors) and 50 edges (regulatory relationships). To answer our fundamental question, whether the network accurately models what actually occurs during cold shock in the yeast cell, we analyzed the results of the modeling as follows. For each gene we evaluated how well the simulated data generated by the model fit the expression profile measured by the microarrays. Factors contributing to the goodness of fit included whether the transcripton factor itself exhibited significant dynamics (in the ANOVA test) and whether the regulators of that transcription actor also showed significant dynamics. This analysis showed that, while a few genes were modeled well, others were not, potentially because they are regulated by other transcription factors that were not present in the network modeled. To investigate this further, we compared the results of the modeling of the YEASTRACT‐derived network to 10 random networks that had the same nodes and and the same number of edges, but had randomized connections between nodes. We found that the random networks had larger values for the least squares error in the estimation and a different structure and degree‐distribution than the YEASTRACT‐derived network, suggesting that the YEASTRACT‐derived network modeled the transcriptional dynamics better. To improve the performance of the model and account for missing regulators, we are now evaluating a family of YEASTRACT‐derived networks by paring down an initial network of 35 nodes systematically down to 15 nodes by removing the next least significant transcription factor one by one from the network. From this analysis we expect to gain insight into the gene regulatory network that controls the cold shock response in yeast. Our working code is available on the GRNmap page (http://kdahlquist.github.io/GRNmap/), and visualization of the network is available on GRNsight (http://dondi.github.io/GRNsight/).

Support or Funding Information

This work was partially supported by NSF award 0921038 (B.G.F., K.D.D.), a Kadner‐Pitts Research Grant (K.D.D., K.G.J., N.E.W.), and a Loyola Marymount University Honors Summer Research Fellowship (K.G.J., N.E.W.).