Anticipating response function in gene regulatory networks

Abstract

The origin of an ordered genetic response of a complex and noisy biological cell is intimately related to the detailed mechanism of protein-DNA interactions present in a wide variety of gene regulatory (GR) systems. However, the quantitative prediction of genetic response and the correlation between the mechanism and the response curve is poorly understood. Here, we report in silico binding studies of GR systems to show that the transcription factor (TF) binds to multiple DNA sites with high cooperativity spreads from specific binding sites into adjacent non-specific DNA and bends the DNA. Our analysis is not limited only to the isolated model system but also can be applied to a system containing multiple interacting genes. The controlling role of TF oligomerization, TF-ligand interactions, and DNA looping for gene expression has been also characterized. The predictions are validated against detailed grand canonical Monte Carlo simulations and published data for the lac operon system. Overall, our study reveals that the expression of target genes can be quantitatively controlled by modulating TF-ligand interactions and the bending energy of DNA.

Keywords: GCMC simulation; gene regulation; protein–DNA networks; response function; statistical mechanics.

Figure 1.

(a) Schematic view of the complex structure of nucleosomes/nucleoid. (b) The zoomed view of a section of nucleosomes/nucleoid. Different types of protein–DNA interactions control the populations of a specific configuration. We denote RNAP–DNA, TF–DNA, TF–RNAP and TF–ligand by ε_RD, ε_TD, ε_TR, ε_TL. The nearest-neighbour interactions for TF and RNAP on the DNA lattice are w_TT and w_RR. w_l represents the long-distance interactions among RNAP and TF molecules. The colour-coding schemes for RNAP, TF and ligand are green, orange and purple. The formation of the DNA loop is also shown here. (c) Protein–DNA interactions are modelled as a lattice of possible binding sites that TFs can occupy. Two possible configurations out of an enormous number of microstates are shown for representation. (d) A correspondence between the bead and 3-D structures of protein and DNA is shown.

Figure 2.

Average occupation number of RNAP and RNAP–TF complex, $\overline{o}$ for non-interacting and interacting systems as a function of RNAP activity. Solid lines and the symbol circles represent theoretical and simulation results. The calculations are done for three different values of interaction energies as shown by three different line colours, red, blue and green. Schematic lattice configurations of short segments of protein–DNA complexes are shown in the respective panels. Panel (a) simple binding between RNAP and DNA. Panel (b), when RNAP binds specifically and non-specifically to DNA. The value of binding interaction for RNAP with non-specific DNA binding sites, ${ϵ}_{\mathrm{RD}}^O=-2k_{B}T$ is considered. Three values of ε_RD are considered for both analyses. Panel (c), simultaneous binding of RNAP and TFs on DNA. The binding interactions of RNAPs and TFs with DNA, ε_RD = −4.68k_BT, ε_TD = −8.88k_BT and λ_T = 1.5 × 10⁻⁶ are considered. Three values of ε_TR are considered. (d) The nearest-neighbour interactions (w_RR) are introduced when both sites are occupied by RNAP. Three values of w_RR are considered for the calculation.

Figure 3.

Average occupation number ($\overline{o}$) of RNAP, TF–RNAP and TF–RNAP–L complexes for non-interacting systems as function of λ_T and λ_L. Solid lines and the symbol circles represent theoretical and simulation results. Schematic lattice configurations of short segments of protein–DNA complexes are shown in the respective panels. (a) Recruitment of RNAP by TF, (b) stimulation of RNAP due to the binding of TF, (c) cooperative stimulation by the dimeric TFs, (d) activation due to the binding of a ligand to TFs. Following binding interaction parameter are used: ε_RD = −4.68k_BT, ε_TD = −8.88k_BT, λ_R = 10⁻⁴, λ_T = 1.5 × 10⁻⁶.

Figure 4.

Various complexes formed in our MCS. The RNAP, TF and ligand are represented by green, orange and magenta beads. The arc in the figure corresponds with the DNA loop. Complexes formed on a linear DNA segment are shown in the first row. In the subsequent rows, we show the possible complexes formed due to DNA loop formation.

Figure 5.

Fraction of different species formed on linear DNA segments for both of the stiff (solid lines) and flexible (dashed lines) DNAs are shown as a function of ligand activity, λ_L. The following values of energy of interactions for various complex formations are considered, ε_RD = −4.68k_BT, ε_TD = −8.88k_BT, ε_TR = −3k_BT, ε_TT = −2k_BT, w_l = −0.2k_BT, w_RR = −5k_BT w_TT = −1k_BT, k_l = 3.8(k_BT/d²). The calculations are done at fixed values of activities for RNAP and TF, those are given by λ_R = 10⁻⁴, λ_T = 1.5 × 10⁻⁶, respectively.

Figure 6.

The population of higher-order oligomeric species formed on the flexible DNA due to the looping. The following values of energy of interactions for various complex formations are considered, ε_RD = −4.68k_BT, ε_TD = −8.88k_BT, ε_TR = −3k_BT, ε_TT = −2k_BT, w_l = −0.2k_BT, w_RR = −5k_BT, w_TT = −1k_BT, k_l = 3.8(k_BT/d²). The calculations are done at fixed values of activities for RNAP and TF, those are λ_R = 10⁻⁴, λ_T = 1.5 × 10⁻⁶, respectively. The average number of loops per configuration is shown at the extreme bottom right panel of the figure.

Figure 7.

Quantification of functional responses in GR systems. (a) NFC for the activation as a function of λ_T. (b) NFC for repression as a function of λ_T. (c) Control activation with a ligand, i.e. DR^A as a function of λ_L. (d) Control of repression with a ligand, i.e. DR^R as a function of λ_L.

Figure 8.

Protein–DNA interaction network for lac operon system at thermodynamic equilibrium. Different binding regions of the promoter of the lac operon gene are depicted by different shades of colour in the lattice. O₁, O₂, O₃ are three operators, where LacR can bind. CRP and RNAP bind to their respective binding regions. Different shapes are used to represent the different biomolecules in the figure. The letters inside the shapes such as R, L, C, I and c correspond with the RNAP, LacR, CRP, IPTG and cAMP molecules in the figure.

Figure 9.

Fold change (FC) and dose–response (DR) as a function of biomolecular activities, λ_L, λ_C, λ_I and λ_c. Solid lines represent theoretical results, and the results obtained from 20 independent GCMC simulation runs are shown by shaded colours. (a) Repression as a function of λ_L. The calculations are done at three values of λ_I. (b) Activation as a function of λ_C. Here the calculation is also done at three constant values of λ_c. (c) Activation as a function of λ_I. The calculations are done for both stiff and flexible DNAs. (d) Activation as a function of λ_c. Here, we also consider two cases such as with and without loop cases. The results obtained from GCMC simulations for the interacting system, where we consider nearest-neighbour interactions among RNAP molecules (w_RR = −5k_BT), are shown by solid green and blue colours in each of the panels.