The transcriptional network that controls growth arrest and differentiation in a human myeloid leukemia cell line

Abstract

Using deep sequencing (deepCAGE), the FANTOM4 study measured the genome-wide dynamics of transcription-start-site usage in the human monocytic cell line THP-1 throughout a time course of growth arrest and differentiation. Modeling the expression dynamics in terms of predicted cis-regulatory sites, we identified the key transcription regulators, their time-dependent activities and target genes. Systematic siRNA knockdown of 52 transcription factors confirmed the roles of individual factors in the regulatory network. Our results indicate that cellular states are constrained by complex networks involving both positive and negative regulatory interactions among substantial numbers of transcription factors and that no single transcription factor is both necessary and sufficient to drive the differentiation process.

Figure 1

Motif Activity Response Analysis (MARA). (a) CAGE tags are mapped to the human genome and their expression is normalized; vertical lines represent TSS positions, and their height is proportional to the normalized expression. (b) Mapped tags are clustered into promoters on the basis of their relative expression, and neighboring promoters are joined into promoter regions. (c) A window of −300 to +100 flanking each promoter region is extracted, multiply aligned and the MotEvo algorithm is used to predict binding sites for known motifs. (d–f) Observed expression of all promoters (d) and predicted site-counts (e) are used to infer motif activities (f). (g) The statistical significance of the regulatory edge from motif to promoter is calculated based on correlation of the promoter expression and motif activity profiles.

Figure 2

Statistical significance and consistency across replicates of the inferred motif activity profiles. Each dot corresponds to a motif. The significance of each motif in explaining the observed expression variation is quantified by the z value of its activity profile (horizontal axis, see Methods). The consistency of the inferred activity profile of each motif is quantified by the fraction of the variance (FOV) in the activity profile across all six replicates (three biological replicates for both CAGE and Illumina), which is reproduced in each replicate (vertical axis, see Methods).

Figure 3

Inferred time-dependent activities of the key regulatory motifs. (a) The time-dependent activity profile of the E2F1-5 regulatory motif as inferred from CAGE (left) and microarray (right) data. The three biological replicates are shown in red, blue and green. (b) The 30 most significant motifs with consistent activity profiles across all replicates (CAGE and microarray) were clustered into nine sets of motifs with similar dynamics. Each panel shows the activity of the members of the cluster (colored curves), the names of motifs contributing and the cluster average activity profile (black).

Figure 4

Predicted core regulatory network of the 30 core motifs. An edge X→Y is drawn whenever the promoter of at least one of the TFs associated with motif Y has a predicted regulatory edge for motif X (z value ≥1.5) and the edge has independent experimental support. The color of each node reflects its cluster membership and the size of the node reflects the significance of the motif. Edges confirmed in the literature, by ChIP or by siRNA are shown in red, blue and green, respectively. In cases where there are multiple lines of support only one evidence type is shown. Supplementary Table 5 shows all predicted edges and their experimental support. GO terms significantly enriched among target genes are shown as white nodes with black edges. FOS/JUN (FOS,B,L1_JUNB,D), CREB (ATF5_CREB3), GABPA (ELK1,4_GABPA,B2).

Figure 5

Validation of predicted target promoter sets using siRNA knockdowns. (a) Difference in the average log expression ratio upon knockdown between predicted target promoters and predicted nontargets (vertical axis) as a function of the z-value cut-off on target prediction (horizontal axis, more stringent cut-offs are on the right) for knockdown of the TF genes MYB (red), SNAI3 (orange), RUNX1 (green) and EGR1 (light blue). (b) As in a but now for knockdown of SPI1 followed by 1 h without treatment (light blue), 24 h without treatment (dark blue), 1 h of PMA treatment (orange) and 24 h of PMA treatment (red). All straight lines are linear regression fits. (c) Pearson correlation coefficients between the average log expression ratio difference of targets and nontargets and the cut-off on target predictions (horizontal axis). Red bars indicate correlation coefficients larger than 0.75 in absolute value; green bars, absolute values between 0.5 and 0.75; and blue bars, less than 0.5. (d) Significance (z value) of the difference in log expression ratio between predicted targets and nontargets (cut-off z = 1.5) for all 28 TFs associated with a motif, measured as a z value (number of standard errors). Red bars correspond to significant changes, that is, greater than two standard errors; green bars, changes between 1 and 2 standard errors; and blue bars, changes less than 1 standard error. siRNA knockdowns were carried out in biological triplicate and knockdown was assessed by qRT-PCR (Supplementary Table 7).

Figure 6

Most significant motif activity changes (as measured by z value, red bars) for four TF gene knockdowns that induce motif activity changes that have a differentiative overlap with the PMA time course of more than 50%. The corresponding motif activity changes observed in the PMA time course are shown as gray bars.