Modeling protein target search in human chromosomes

Several processes in the cell, such as gene regulation, start when key proteins recognize and bind to short DNA sequences. However, as these sequences can be hundreds of million times shorter than the genome, they are hard to find by simple diffusion: diffusion-limited association rates may underestimate in vitro measurements up to several orders of magnitude. Moreover, the rates increase if the DNA is coiled rather than straight. Here we model how this works in vivo in mammalian cells. We use chromatin-chromatin contact data from Hi-C experiments to map the protein target-search onto a network problem. The nodes represent DNA segments and the weight of the links are proportional to measured contact probabilities. We then put forward a diffusion-reaction equation for the density of searching protein that allows us to calculate the association rates across the genome analytically. For segments where the rates are high, we find that they are enriched with active gene starts and have high RNA expression levels. This paper suggests that the DNA's 3D conformation is important for protein search times in vivo and offers a method to interpret protein-binding profiles in eukaryotes that cannot be explained by the DNA sequence itself.

Figure 1

Modeling protein search on DNA as search on a weighted network. (a) Hi-C map: the colors represent the number contacts (in ${log}_{10}$) between 10 kilo-basepair (kb) DNA fragments in a part of human chromosome 21. Dark pixels indicate many contacts (dynamic range: $\sim {10}^0–{10}^5$). (b) Schematic representation of the model. Key parameters: jump rates between nodes $i$ and $j$ (${\omega }_{ij}$), unbinding rate to the bulk ($k_{\mathrm{off}}$), and rebinding rate ($k_{\mathrm{on}}$). Red circles represent searching proteins. (c) Coarsed network representation of the Hi-C map in (a). Each node represents a 160 kb fragment. The link weights $v_{ij}$ are proportional to the number of Hi-C contacts. We assume that ${\omega }_{ij}\propto v_{ij}$. Node numbering refers to positions along the DNA.

Figure 2

Predicted genome-wide association rates ${\widehat{K}}_a$ for chromosomes 1-21. ${\widehat{K}}_a$ vary by several orders of magnitude but the 95% confidence interval is a few percent of the mean $\langle {\widehat{K}}_a\rangle =1$ (${\widehat{K}}_a=1.0\pm 0.0027$). We calculated ${\widehat{K}}_a$ from Eq. (6) (${\widehat{k}}_{\mathrm{off}}\ll 1)$ with ${\widehat{k}}_{\mathrm{off}}=0.002, {\widehat{k}}_{\mathrm{on}}=0.001$, and ${\rho }_0=0.5$. ${\widehat{k}}_{\mathrm{off}}\gg 1$ show similar behavior, see Appendix pp5

Figure 3

Nodes with many gene starts have higher predicted association rates than nodes with few gene starts. We define the gene starts as the Transcription Start Sites (TSSs). The curves represent predicted association rates ${\widehat{K}}_a$ to nodes with active TSSs (gray), inactive TSSs (green), and any TSS type (pink). The active TSSs have higher ${\widehat{K}}_a$ than the genome-wide average (dashed), whereas nodes with inactive TSSs (green) are below (except one data point). The symbols represent the average ${\widehat{K}}_a$ (Eq. (6), ${\widehat{k}}_{\mathrm{off}}\gg 1$) and the coloured areas show the 95% confidence interval. Parameters (dimensionless, see Appendix pp3): ${\widehat{k}}_{\mathrm{off}}=0.002, {\widehat{k}}_{\mathrm{on}}=0.001$, and ${\rho }_0=0.5$. We omitted data points with less than 7 TSSs per node. The ${\widehat{k}}_{\mathrm{off}}\gg 1$ case has the same trend, see Appendix pp5.

Figure 4

Highly transcribed nodes are are easy to find. (a) Genome-wide distribution of ${\widehat{K}}_a$ for all nodes divided into 20 groups based on their RNA expression levels. White circles: the median; horizontal bars: the mean; the dashed line; genome-wide average $\langle {\widehat{K}}_a\rangle =1$. (b) Predicted ${\widehat{K}}_a$ as function of RNA expression level for nodes with at least one active TSS (grey) and no active TSSs (purple). Same grouping procedure as in (a). Nodes with active TSSs tend to be above the genome-wide average (18 points above $\langle {\widehat{K}}_a\rangle$), while most nodes with no active TSSs are below (6 points above $\langle {\widehat{K}}_a\rangle$). The shaded areas show the 95% confidence interval. Parameters: ${\widehat{k}}_{\mathrm{off}}=0.002, {\widehat{k}}_{\mathrm{on}}=0.001$, and ${\rho }_0=0.5$. We used data for cell line GM2878 but find similar results for K562, see Appendix pp7.

Figure 5

Association rate increases with the fraction of transcribed DNA. This is similar to the RNA expression level which also increases with the association rate. Parameters: ${\widehat{k}}_{\mathrm{on}}=0.001$ and ${\widehat{k}}_{\mathrm{off}}=0.02$.

Figure 6

Comparison of the analytical approximations to the exact, time-integrated association rate Eq. (D1), using Hi-C data from chromosome 21 at 160 kb resolution. (a) The blue area represents the fast target-search regime where ${\widehat{k}}_{\mathrm{off}}\ll 1$. The red area shows the opposite regime, ${\widehat{k}}_{\mathrm{off}}\gg 1$, where the system is close to its steady-state and target finding is slow. We put the target in the middle of the system $a=N/2$. In the two lower panels, we show the association rates' profile to all possible targets on chromosome 21 when ${\widehat{k}}_{\mathrm{off}}\ll 1$ (b) and ${\widehat{k}}_{\mathrm{off}}\gg 1$ (c). Note the green and red circled dots in (b) and (c), respectively. These correspond to the encircled parameter values in panel (a).

Figure 7

Same as Fig. 3 in the manuscript but for large ${\widehat{k}}_{\mathrm{off}}$ instead of small. We recover the same trend: easy-to-find regions tend to have many gene starts.

Figure 8

The association rate vs the nodes' strength $V_a$ (sum of all link weights). (a) All data points follow the universal law ${\widehat{K}}_a={\widehat{k}}_{\mathrm{on}}+{\gamma }_{a_1}{V}_a$ [Eq. (C4)]. The density ${\rho }_0=0.5$ is kept fixed in all four cases as we increase $k_{\mathrm{on}}$. (b) Association rate during steady state (${\widehat{k}}_{\mathrm{off}}\gg 1$) calculated from Eq. (C8). The behavior is not universal as it depends on the number of nodes $N$ for each chromosome contact network. The dashed line represent the analytical formulas. As $N$, we used the mean chromosome size.

Figure 9

Same as Fig. 4 but for cell type K562 instead of GM12878. The figure shows that highly transcribed nodes—in terms of high levels of RNA expression—are are easy to find (low association rates) even if there are active genes there or not. See the caption of Fig. 4 in the manuscript for a more details on what the curves represent.

Figure 10

Genome-wide association rates ${\widehat{K}}_a$ at 40 kb resolution. (a) $K_a$ for slow unbinding rates predicted by Eq. (C7). This complements Fig. 2 showing $K_a$ for fast unbinding rates. As in Fig. 2, most values deviate only from the average by a few percent ($0.7531\pm 0.0018$, 95% confidence interval). We used these parameters: ${\widehat{k}}_{\mathrm{off}}=2, {\widehat{k}}_{\mathrm{on}}=0.001$ and ${\rho }_0=0.0005$. (b) $K_a$ predicted by Eq. (C7) (slow unbinding rate) after we removed 25% of the weakest contacts. The genome-wide average is $\langle {\widehat{K}}_a\rangle =1\pm 0.0028$ (95% confidence interval). We used same parameters as in Fig. 2: ${\widehat{k}}_{\mathrm{off}}=2, {\widehat{k}}_{\mathrm{on}}=0.001$ and ${\rho }_0=0.0005$.

Figure 11

(Left) A boxplot of the association rates for 40kb bins with significant DNase-seq signal (Right) Association rates for all 40-kb bins across the genome. The $y$ axis is cut at ${\widehat{K}}_a=4,$ which removes 13 outliers out of 66 452 data points.

Figure 12

(Top) Boxplot of the gc content for 40kb bins with a single transcription start site (TSS). The bins are separated in two groups—active and inactive—based on the transcriptional activity surrounding each TSS. As a reference, we show the gc content for all Hi-C bins (“genome-wide”) (Bottom) The range of association rates for the same three groups as to the left. The $y$ axis is cut at ${\widehat{K}}_a=3$ which removes about 10 outliers.