Computational challenges in modeling gene regulatory events

ABSTRACT

Cellular transcriptional programs driven by genetic and epigenetic mechanisms could be better understood by integrating “omics” data and subsequently modeling the gene-regulatory events. Toward this end, computational biology should keep pace with evolving experimental procedures and data availability. This article gives an exemplified account of the current computational challenges in molecular biology.

Keywords:

• epigenetic mechanisms
• gene regulation
• systems biology
• transcription
• transcription factors

Introduction

Genetic and epigenetic mechanisms act together to organize the fine-tuned gene regulatory circuitries underlying cell-type identities and functions.^1,2 Misregulation of these pathways has been associated with diseases. Genetic variations are associated with a multitude of disorders, including developmental defects and cancers.³ Sequence-dependent targeting of the critical transcription factors at enhancers is required to generate a gene expression program underlying cell-type specification and function.^4,5 In line with this need, mutations in regulatory elements such as enhancers have been shown to disrupt their ability to bind host transcription factors, leading to misregulation of gene expression.⁶ With the advent of “epigenetics”, a number of mechanisms have been identified that regulate gene expression without influencing the genetic sequence. These mechanisms include the modification of histones, DNA methylation, and non-coding RNA, which modulate the chromatin structure and alter transcriptional outcomes.¹ The study of these epigenetic mechanisms at the genome-wide level in a high-throughput fashion is termed “Epigenomics”.⁶

The arrival of next-generation sequencing technology led to many critical discoveries in multiple branches of the life sciences by providing a comprehensive overview of the genomic, transcriptomic, and epigenomic maps of cell types of interest.⁷ These assays have facilitated the discovery of novel gene regulatory features.⁸ These methods can be implemented toward (1) obtaining gene expression profiles (RNA-seq), (2) interrogating chromatin and transcription factor landscapes (ChIP-seq), (3) reconstructing chromatin structures and long-range chromosomal interactions (e.g., ATAC-seq, DNAse-seq,. and HiC), and (4) uncovering the methylation and other modification states of DNA (e.g., MEDIP-seq and BiS-seq). In parallel, the evolution of proteomic tools has facilitated discoveries in the gene regulation field, including the identification of the differentially expressed proteins and protein–protein interactions involved in transcriptional regulation and the discovery of novel histone modifications and their modifiers, readers, and erasers.⁹ Computational tools for the analysis of both genomics and proteomics datasets are increasingly advancing and routinely updated.^9,10

These advances in the large data acquisition have also invited challenges that demand active systematic efforts to tackle them. In any cell at a given time point, a large number of transcription factors are functioning in concert with one another and the chromatin to guide the gene expression program underlying cell identity and functions.¹¹ Given the limitations of the current experimental measures, it is not possible to fully explore these dimensions. However, advanced bioinformatics and systems biology tools offer the opportunity to gain better insights into the complexity of gene regulatory events and to generate comprehensive quantitative models of these networks. Furthermore, there is a growing need to develop algorithms that integrate data sets from various “omics” studies, such as genomics and proteomics experiments, to improve the quantitative models of gene regulatory networks by integrating events, including protein–protein interactions and post-transcriptional (RNA) and translational modifications (protein). In the wake of the increasing data availability, there is a need to develop improved methods for analysis that take into account pre-existing information to allow deeper biological insights in the analysis.¹² Integration and making the data availability process smoother with a one-point entry access remain important challenges.¹³ This article outlines critical and immediate computational challenges and reports possible measures to overcome them with a special focus on gene regulation.

Challenges in molecular systems biology

Modeling of gene regulatory networks

Comprehensive and integrative modeling of gene regulatory networks is one of the key challenges in current molecular systems biology. Toward this end, transcription factor binding preferences and the effects of binding on chromatin and the gene regulatory network mediated by these factors need to be understood.

Systematic prediction of transcription factor binding sites

DNA binding factors can function both as activators and repressors of transcription, and their dynamics at various loci determine the ultimate transcriptional response. A set of transcription factors have also been described as so-called “pioneer factors” that can bind closed chromatin regulatory regions and induce subsequent activities to create accessible and active chromatin landscapes, thereby increasing gene expression.¹⁴ By contrast, settler factors function on genomic loci made accessible by pioneer factors to drive subsequent transcriptional changes.¹⁴

Most transcription factors bind in a highly sequence-specific manner, and the consensus sequence (i.e., the motif) of these factors is unique to them or their families.¹⁵ Interestingly, although motif sequences appear numerous times in the genome, only a small fraction of these motifs are actually utilized. Motif analysis can be complemented with the analysis of secondary DNA structures and sequential features in the binding sites and flanking regions to considerably increase accuracy of transcription factor binding prediction.^16,17 Thus, efforts should be made in these directions to improve the computational methods and allow predictions of transcription factor binding with high accuracy (Fig. 1A). Furthermore, the chromatin state of the target sites should be added as an additional feature in this analysis. This analysis would allow the cell-type specific binding patterns of transcription factors to be deciphered.¹⁸

Figure 1. Computational modeling of cell type-specific transcription factor binding sites. Step 1: Construction of a database of classification models for every transcription factors by modelling features of bound and unbound motifs. For the flanking regions of these motifs, information is extracted from DNA structure (DNAshape), sequence features, and chromatin state dynamics to build the classifying model (i.e., SVM, RBF, or Regression). If the CV accuracy of the classification model is greater than 75% (arbitrary cutoff), the model is retained in the database. Step 2: Utilize the constructed model database to predict cell type-specific binding: for a set of genomic regions of interest, information is collected from DNAshape, sequence features, and chromatin states. Every classification model in the database is scanned. The factor for which the AUC is higher potentially binds these genomic regions. Thus, factors for which the AUC scores are higher are reported.

Classifying transcription factors as pioneers and settlers is very important because although both of these factors may bind similar target sites, they function via different mechanisms and for different purposes. Recently developed algorithms such as PIQ¹⁴ represent great advances in recognizing transcription factor binding sites and identifying pioneer-settler pairs by modeling the shape and the magnitude of chromatin accessibility at factor binding sites. These binding events can be validated using a ChIP-seq assay and used to calculate the accuracy of the model in an unbiased manner.¹⁹

Deciphering transcriptional responsiveness to transcription factor binding

An expected consequence of transcription factor binding is the transcriptional regulation of target genes. However, it has been extremely difficult to determine the actual transcriptional response of a transcription factor binding to its target site in the native environment. The reasons for such difficulty are multifold; first, analysts face the challenge of assigning a target gene with a transcription factor binding event. In today's genomics standards, in cases when the binding of a factor occurs at a promoter, it is assumed that the associated gene is regulated. However, if the binding of the factor is distal, then in most cases the nearest gene is considered to be regulated, which generates noise in the computational analysis.²⁰ In such cases, researchers would benefit from the availability of genome-wide maps of long-range chromosomal maps (e.g., Hi-C) and enhancer marks (e.g., H3K27ac) for the same cell type. A recently developed machine learning-implemented computational method (TargetFinder) predicts potential target genes of distal factor binding based on a range of features learned from the Hi-C data.²¹ A second issue is that the levels of binding do not correlate well with the actual transcriptional response even when the real target genes are deciphered.²² Emerging evidence suggests that both the DNA binding and transcriptional responses are a function of the DNA sequence and structure as well as the epigenetic state.²³ Therefore, the current big issue in quantitative genomics is to decipher the determinants of transcriptional factor targeting and transcriptional output.

Mathematical modeling of gene regulatory networks

With the identification of genome-wide targets of a given transcription factor and its influence on gene expression, one ultimately aims to generate a quantitative model of the gene regulatory networks underlying cell identity or functions. The reconstruction of gene regulatory networks from these minimal information datasets remains a challenging task. Methods for this systems biological characterization have evolved in recent years from traditional co-expression networks to dynamic modeling approaches.²⁴

Reconstruction of networks based on gene expression profiles

Co-expression networks are reconstructed on the assumption that genes that are temporally co-expressed are functionally related and are regulating the transcriptional state of genes that respond sequentially thereafter. Popular tools such as WGCNA use this principle to reconstruct networks.²⁵ Conversely, ARACNE also reconstructs gene regulatory network based on gene expression profiles but removes indirect nodes in the network using data processing inequality (DPI) to obtain results that are close to real networks.²⁶

Although this model is helpful for the prediction of gene functions, it does not provide resolution of singular gene regulatory events.²⁷

Combination of binding and transcriptional readout factor-based methods

Popular tools such as GeneProf and EpiRegNet^28,29 are based on studies in which the binding of a transcription factor/epigenetic regulator and/or the transcriptional response and/or the co-expression in a time course data set can be integrated to elucidate sequential networks. Because these methods employ co-expression or information from a ChIP-sequencing or motif analysis, a combinatorial response of multiple active transcription factors cannot be determined.

Mathematical modeling integrating biological knowledge

The knowledge-based modeling protocol is the most robust approach to reconstruct realistic gene regulatory events. Biological knowledge can be integrated into a mathematical model that will provide a measure of accuracy to the reconstruction of the network.^30,31 For example, if the binding sites of a transcription factor (X) target (Y) are deciphered by either ChIP-seq or prediction and the target exhibits significantly changed expression during the progression of a particular biological process, then the network node X–Y can be reconstructed. The summation of all nodes for all active transcription factors can be modeled using a probabilistic approach, such as Bayesian modeling. A similar approach is implemented in MARA (Motif Activity Response Analysis); however, this approach is limited to investigating promoter regions and offers no possibility of integrating ChIP-seq and RNA-seq-derived information.³²

Therefore, given the current lack of standard methods to reconstruct gene regulatory networks, it is essential to develop methods that integrate information from ChIP sequencing and RNA sequencing experiments in combination with motif energies to reconstruct gene regulatory events.

Integration of multiple level datasets

Most biological processes involve a combinatorial contribution of transcriptional regulation, alternative splicing, post-translational protein modification, and complex protein–protein interactions. Thus, to deeply understand a biological process it is necessary to integrate data sets spanning these readouts, including transcriptomics, epigenomics, and proteomics.³³

Integration of datasets arising from similar biological contexts

A range of data sets exist for cell types that have been extensively studied, such as ES cells. For these cell types, ChIP-seq for all available transcription factors and histone modifications as well as DNA methylation maps could be integrated to generate a comprehensive overview of the genome-wide and locus-wide transcription factor, epigenetic and transcriptional landscapes. These datasets have proven very useful for obtaining a better understanding of how a transcription factor communicates with different flavors of regulatory factors and chromatin feature combinations at distinct loci.³⁴ Furthermore, genetic data sets, such as DNA structures, the sequence conservation index and sequence patterns, can be integrated. All of this information could be made available in a cell type-specific binary file, which could be used for subsequent genomics studies (Fig. 2).

Figure 2. Integration of multilevel datasets. For every regulatory region (i.e., promoter and enhancer) in any cell type, data from different experiments, such as protein–protein interaction, chromatin state and transcription state in a temporal (t)/treatment (T)-specific manner, should be recorded in files. These files should be integrated in every genomic analysis to provide regulatory events enriched at the sites critical for a given analysis.

Challenges for the integration of multilevel datasets

Many attempts have been made to comprehensively collect and integrate multiple data types for a particular cell type and process. For example, the TCGA (The Cancer Genome Atlas project) is aimed toward the integration of the quantification of mRNA, miRNA, protein, DNA methylation, copy number alterations, and somatic chromosomal aberrations in cancer cells.³⁵ For this integrative analysis, it is necessary to implement careful statistical analysis due to the challenges in normalizing quantitative data from different experiments. Recently, penalized regression models were used to integrate “omics” data from three sources (genetic variations, DNA methylation, and gene expression) to reveal genes that contributed to bladder cancer.³⁶ Attempts have also been made to implement network-based approaches that integrate protein–DNA and protein–protein interactions.³⁷ The recently published tool “Network Analyst” also uses a network model to integrate differential gene expression data with protein–protein interaction maps.³⁸ However, because the techniques and their readouts are not directly comparable due to their different sensitivities and varied dynamic ranges, it has been difficult to generate a comprehensive computational model for the regulation of gene expression.³⁹ Therefore, it is necessary to develop robust statistical procedures capable of integrating cell type-specific omics data and to construct a database of these models to facilitate integrative research.

Optimizing current tools and coping with data availability

An additional computational challenge in large-scale genomics studies is the need to routinely update the procedures used to analyze the genomic data. For instance, peak calling is a critical step in the analysis of ChIP-seq data sets in which there is a need for standardized methods that ideally should take into account several parameters, such as differential peak profiles, abundance of the factor/mark on chromatin, and antibody affinity. Methods to sensitively incorporate these inputs will markedly improve the information gain from genomics datasets. Additionally, better alignment and analysis tools need to be developed to analyze the repetitive DNA sequences that constitute a large part of the genome.⁴⁰ Furthermore, different laboratories use different tools to analyze similar data sets, which makes it harder to perform comparative genomics. In summary, the genomics community faces critical challenges in light of the absence of strict standards for experimental and computational procedures, which can cause severe problems in the proper interpretation of large-scale data and the reproducibility of the findings.⁴¹

Recently, single-cell genomics has become an extremely powerful tool to identify distinct cell populations based on distinct expression and/or chromatin profiles. In addition to single-cell RNA-seq, it is possible to perform ATAC-seq to interrogate chromatin accessibility and Hi-C to analyze higher order chromatin structures at the single-cell level.⁴² These data are generally very noisy and require a large number of replicates to reach conclusions. Thus, integrating “omics” data from single cells remains a challenge for varied reasons, including heterogeneity and complicated experimental and computational methods. A number of laboratories are investing in the development of tools to improve data analysis and interpretation from these datasets.⁴³

An important discussion in recent years has been how to handle the accumulating large-scale data sets more efficiently, and efforts have been made to create one point entry for all of these databases. Many databases, such as CISTROME, ENCODE, and NCBI epigenomics, have been made available. Furthermore, consortia such as STATegra are actively engaged in developing methods for the integration of “omics” data.⁴⁴ However, many of the databases are not updated on a real-time basis, which pushes researchers to manually search for the required datasets. These databases should also encourage the integrated analysis of existing data sets with an easy workflow.

Conclusion

With the growth of experimental procedures for the high-throughput mapping of molecular events as well as increasing data availability and computational resources, it is necessary for computational biologists to put systematic efforts into quality checks, the development of robust but standard analytical tools, the integration of all existing datasets and the development of predictive models that make the best use of the existing genomics data resources. The current challenge in computational biology is to decipher the cell type-specific transcription factor, chromatin and transcriptional landscapes and to integrate these data to evolve quantitative models of gene regulatory networks. At the same time, it is extremely important to routinely update computational tools to deal with the rapidly evolving experimental methods and data availability. An important aspect is to develop standard experimental and computational methods for similar data sets worldwide and to provide an easy platform for integrative analysis.

Disclosure of potential conflicts of interest

No potential conflicts of interest were disclosed.