Repression by RB1 characterizes genes involved in the penultimate stage of erythroid development

Abstract

Retinoblastoma-1 (RB1), and the RB1-related proteins p107 and p130, are key regulators of the cell cycle. Although RB1 is required for normal erythroid development in vitro, it is largely dispensable for erythropoiesis in vivo. The modest phenotype caused by RB1 deficiency in mice raises questions about redundancy within the RB1 family, and the role of RB1 in erythroid differentiation. Here we show that RB1 is the major pocket protein that regulates terminal erythroid differentiation. Erythroid cells lacking all pocket proteins exhibit the same cell cycle defects as those deficient for RB1 alone. RB1 has broad repressive effects on gene transcription in erythroid cells. As a group, RB1-repressed genes are generally well expressed but downregulated at the final stage of erythroid development. Repression correlates with E2F binding, implicating E2Fs in the recruitment of RB1 to repressed genes. Merging differential and time-dependent changes in expression, we define a group of approximately 800 RB1-repressed genes. Bioinformatics analysis shows that this list is enriched for terms related to the cell cycle, but also for terms related to terminal differentiation. Some of these have not been previously linked to RB1. These results expand the range of processes potentially regulated by RB1, and suggest that a principal role of RB1 in development is coordinating the events required for terminal differentiation.

Keywords:

• cell cycle
• erythroid differentiation, pocket proteins
• Retinoblastoma-1 protein
• transcription

Introduction

RB1 is the first tumor suppressor gene discovered and founding member of the pocket protein family. RB1 and the related proteins, p107 and p130, are adaptors that interact with other proteins through a conserved protein fold. Of particular significance, pocket proteins bind E2F transcription factors and repress their activity (for a review see refs. 1,2). E2Fs are categorized as activators (E2F1, E2F2, and E2F3a) or repressors (E2F3b, E2F4, E2F5, E2F6, E2F7, and E2F8). E2F4 and E2F5 repress genes by recruiting pocket proteins, whereas E2F6-8 are pocket protein-independent repressors. Activator E2Fs can also function as repressors, by recruiting RB1.

E2Fs operate downstream of a complex regulatory network. Cyclins, cyclin-dependent kinases, and pocket proteins integrate extracellular signals with the cell cycle machinery. Elegant models based on in vitro studies have assigned specific activities to cyclins and Cdks; however, genetic experiments in mice have shown that there is redundancy in the network. Mice deficient for all D-type cyclins survive until midgestation but die of severe anemia and heart defects.³ Single D-type cyclin mice exhibit tissue specific defects; cyclins D2 and D3, but not D1, support erythropoiesis,⁴ and cyclin D3 regulates the number of terminal erythroid cell divisions.⁵ Mice with combined deficiencies of cyclins E1 and E2 survive at least until birth.⁶

D- and E-type cyclins form complexes with Cdk2, Cdk4, and Cdk6. Mice without Cdk4 and Cdk6 die midgestation of severe anemia.⁷ Embryos lacking Cdk2 and Cdk4 die of heart defects.⁸ Embryos lacking all 3 Cdks (Cdk2, Cdk4, and Cdk6) exhibit haematopoietic and heart defects but still survive to mid-gestation.⁹ By contrast, Cdk1 deficient embryos fail to reach the blastocyst stage; further, Cdk1 is required for S phase entry and RB1 phosphorylation in cells lacking other Cdks. Thus, Cdk1 has an integral role in cell cycle progression, whereas the interphase cyclin-Cdk complexes act in specific developmental contexts. Cyclin D-Cdk4/6 complexes are RB1-specific kinases.¹⁰ Given the role these complexes have in erythropoiesis, it follows that one or more of the pocket proteins is active in erythroid cells, and must be inhibited at some phase of erythroid development.

The role of RB1 in development has been assessed in mice. Germline RB1 deficiency causes developmental defects,^11-13 due to placental insufficiency.¹⁴ Erythropoiesis can be restored, although some developmental defects persist.^14-20 Erythropoiesis is partially rescued by concurrent E2F2 deficiency, which establishes E2F2 as a target of RB1-mediated repression.²¹ However, even as a consensus emerges regarding the phenotypic effects of RB1 deficiency in the erythroid lineage, questions remain. Among these is whether there is redundancy between the pocket proteins in erythroid cells, and whether the role of the pocket proteins is limited to cell cycle regulation. In the present study, we address these questions through the use of mouse genetics and expression profiling. Our results suggest that RB1-mediated repression of E2F target genes has a central role in erythroid development: regulating the transition of erythroid cells from precursors to a fully differentiated state.

Materials and Methods

Analysis of mice

Conditional Rb1^fl/fl mice,²² p107^−/− and p130^−/− mice,^23,24 and Epor-cre mice²⁵ were bred to obtain embryos with all 8 possible combinations of homozygous mutations. Timed mated pregnant females were anesthetized with inhaled isoflurane and euthanized by cervical dislocation. Embryos were harvested and fetal livers obtained by dissection. Animal studies were performed under protocol in accordance with the guidelines established by the Institutional Animal Care and Use Committees of St. Jude Children's Research Hospital and the New York Blood Center.

Genotyping

Embryos were genotyped by PCR (Qiagen Taq, Q-solution, dNTPs, and 10X coral load buffer), with the following primers and conditions.

Rb1: (94°C × 5 min) 1 cycle; (94°C × 30 sec, 58°C × 30 sec, 72°C × 40 sec) 35 cycles; (72°C × 8 min) 1 cycle.

5′-ggcgtgtgccatcaatg-3′; 5′-ctcaagagctcagactcatgg-3′

Wild type = 235 bp. Unrecombined floxed allele = 283 bp.

p107: (94°C × 5 min) 1 cycle; (94°C × 10 sec, 60°C × 1 min, 68°C × 2 min) 32 cycles; (72°C × 10 min) 1 cycle.

5′-tgtcctgagcatgaacagac-3′; 5′-tcgctggcagtctgagtcag-3′; 5′-acgagactagtgagacgtgc-3′

Wild type = 280 bp. Mutated allele = 330 bp.

p130: (94°C × 5 min) 1 cycle; (94°C × 10 sec, 60°C × 1 min, 68°C × 2 min) 32 cycles; (72°C × 10 min) 1 cycle.

5′-tacatagtttccttcagcgg-3′; 5′-gaagaacgagatcagcagc-3′; 5′-acggatgtcagtgtcacg-3′

Wild type = 230 bp. Mutated allele = 330 bp.

Epor-cre: (94°C × 5 min) 1 cycle; (94°C × 10 sec, 60°C × 1 min, 68°C × 2 min) 32 cycles; (72°C × 10 min) 1 cycle.

5′-gtgtggctgccccttctgcca-3′; 5′-ggcagcctgggcaccttcac-3′; 5′-caggaattcaagctcaacctca-3′

Wild type = 431 bp. Mutated allele = 679 bp.

Embryonic blood smears

Embryos were bled by puncturing the carotid artery with a 2 microliter glass microcapillary tube and the blood expelled to prepare a peripheral blood smear. Blood smears were stained with Wright-Giemsa, and the erythrocyte morphology examined by light microscopy. Blood smears were viewed with a Nikon Eclipse E600 microscope and Nikon Plan Apo 40× objective, NA 0.95. Images were acquired with a Nikon DXM1200 digital camera and ACT-1 software.

Protein assays and antibodies

For Western blots, proteins were separated with NuPAGE Bis-Tris gels (Life Technologies) and transferred to 0.45 μm PVDF membrane. We used the following primary antibodies: RB1 (554136), Ter119 (553670), N-terminal p130 (610262) from BD Transduction Labs; p107 and C-terminal p130 (C-20) from Santa Cruz; β-actin (A5441) from Sigma. Secondary antibodies: Donkey anti-rabbit IgG HRP (NA934V) and Sheep anti-mouse IgG HRP (NA931V) from GE Healthcare Life Sciences.

Cell culture

Fetal liver cells (FLC) were disaggregated in 1 ml of phosphate-buffered saline by passage through a 21 gauge needle. FLC were cultured at a final concentration of 5 × 10⁵ cells/ml in complete medium (30% fetal bovine serum, 1% deionized bovine serum albumin, 0.001% monothioglycerol, 2 mM glutamine, and penicillin-streptomycin in Iscove's modified Dulbecco medium). Cells were maintained in a humidified incubator at 37°C, 5% CO2 and harvested by pipeting at the indicated times. Erythroid colonies were cultured in methylcellulose (M3234 for CFU-E; M3534 for BFU-E) from Stem Cell Technologies. M3534 contains rmSCF, rmIL3, rhIL6, rh-insulin, and human transferrin (iron saturated). M3234 contains no cytokines. We added rhEPO to both to a final concentration of 4 units/ml. CFU-E and BFU-E were scored at 2 days and 7 days, respectively.

Flow cytometry

For cell cycle analysis, we centrifuged 1 × 10⁶ FLC at 400g for 5 min and aspirated the supernatant. The cells were resuspended in 1 ml of PI solution (0.05 mg/ml propidium, 0.1% sodium citrate, 0.1% Triton X-100) and stored at 4°C in the dark. Samples were treated with 4 μg/ml RNase for 30 min at room temperature. Samples were filtered through 40 micron mesh just prior to analysis. Samples were run on a FACSCalibur analyzer (BD Biosciences) and the data processed with ModFit software (Verity Software House).

For Brdu incorporation studies, 10 μM Brdu was added to FLC in culture. For time course experiments, Brdu was added 45 min before each time point. For pulse chase, Brdu was added for 45 min initially, then washed out, and the cells harvested at hourly intervals thereafter. Cells were stained with the FITC Brdu flow kit (BD Transduction Labs) in accordance with the manufacturer's instructions. Samples were run on an LSR analyzer (BD Biosciences) and the data processed with FACSDiva software (BD Biosciences).

For routine flow cytometry, 2 × 10⁶ FLC were blocked with γ-globulin and stained with CD71-PE and Ter119-APC (BD PharMingen). Samples were run on an LSR analyzer (BD Biosciences) and the data processed with FACSDiva software. Sorting was performed on FACSVantage and FACS Aria cell sorters (BD Biosciences).

Microarray procedure

FLC were obtained from E13.5-14.5 embryos. 3-4 × 10⁶ FLC were washed once in phosphate-buffered saline and total RNA extracted with 1 ml of RNA-Bee (Tel-Test) and 0.1 ml of chloroform, in accordance with the manufacturer's instructions. RNA quality was verified with an Agilent 2100 Bioanalyzer. Total RNA (1-10 μg) was processed in the Hartwell Center microarray core facility following the Affymetrix eukaryote one-cycle target labeling protocol (701025, Rev. 6). Biotin-labeled cRNA (15 μg) was hybridized overnight at 45°C to Mouse Genome 430 2.0 GeneChip arrays which correspond to more than 39,000 transcripts. After staining and washing, the arrays were scanned, and expression values were summarized using the MAS5 algorithm as implemented in the GCOS v1.4 software (Affymetrix). Signals were normalized for each array by scaling to a 2% trimmed mean of 500. Detection calls (present, absent, and marginal) were determined using the default parameters.

Microarray data analysis

In the first experiment, we performed expression profiling on 34 fetal liver samples. All Rb1 genotypes were floxed/floxed. There were 12 samples without Epor-cre (RB1 replete) and 22 samples with Epor-cre (RB1 deficient). To maximize the number of samples available for the RB1 comparison, we used all 34 samples, without regard to the genotypes of the other pocket proteins. Because RB1 had a greater effect, for p107 and p130 we excluded RB1-deficient samples, leaving 12 samples for the comparison. The second experiment was a time course (0, 24, and 48 hours). In this experiment, there were 12 RB1-replete samples and 15 RB1-deficient samples. All raw data files (cel) are submitted to the GEO database (GSE67285). These experiments comply with MIAME (minimal information about a microarray experiment) standards.

Butterfly plots were generated on the Gene Set Enrichment Analysis server.²⁶ The data was permuted 1,000 times with randomized column headings to establish 1% and 5% significance thresholds. Hierarchical clustering was performed on the GenePattern server.²⁷ Next, the data were ranked employing the Signal-to-Noise (S2N) metric. The genes were divided into 16 bins (1,000 genes/bin), based on rank order, with the genes most derepressed in the absence of RB1 in bin 1, and those most repressed in bin 16. Finally, the data were resorted according to their original hierarchical clustering position. The S2N metric is the difference in the mean of the 2 groups in a comparison, divided by the sum of the standard deviations. K-means clustering of human expression data^28,29 was performed on the GenePattern server. Microarray data and external data sets were merged with Microsoft Access and Excel, 2010. False discovery rates (FDR) were determined by comparing gene sets corresponding to a group of interest (e.g., E2F1 bound genes) against the rank order of the expression data, based on the K-S statistic (GSEA server). Functional annotation was performed on the DAVID server, employing the default settings.³⁰

Results

p107 and p130 do not compensate for RB1 in erythroid cells

A plausible explanation for the modest effect of RB1 deficiency on erythropoiesis in vivo is that the other pocket proteins supply redundant activity in erythroid cells. Pocket protein redundancy exists during development,^23,31-33 and in retinoblastoma^34,35 and other cancers.³⁶ To address this, we generated mice lacking pocket proteins in the erythroid lineage. We bred mice with a conditional Rb1 allele with mice that express Cre in an erythroid-restricted pattern.^22,25 RB1 accumulates in erythroid cells in its hypophosphorylated state during terminal differentiation (Figure S1A).³⁷ Rb1^fl/fl;Epor-cre FLC did not express RB1 (Figure S1B); further, they lacked RB1 from the CD71+Ter119- stage onward (Figure S1C). Rb1-floxed mice were bred with p107 and p130 mice.^23,24 FLC from p130^−/− embryos expressed a faster migrating protein that reacts with carboxyl-terminal p130 antibody (Figure S1B). The mutant p130 allele used in our studies has a similar effect on fibroblast immortalization as another null allele;^38,39 however, we cannot exclude the possibility that this protein has retained some p130 function.

Because combined p107 and p130 deficiency is associated with neonatal lethality,²⁴ we limited our studies to midgestational (E13.5-14.5) embryos. Genetic crosses produced 250 embryos representing all possible combinations of RB, p107, and p130 deficiency (Table 1). All combinations, including the triple knockout (TKO), were obtained in the expected proportions. Differences in the number of FLC and erythroid colony-forming units were not significant (data not shown). Thus, regarding erythroid develop-ment, p107 and p130 do not provide redundant activity.

Table 1. Midgestational Rb1 family erythroid TKO embryos are viable

Genotype	Predicted	Actual
Rb1(fl/fl);pl07(+/ˆ);p130(+/ˆ);Epor(+/+)	50	46
Rb1(fl/fl);pl07(+/ˆ);p130(+/ˆ);Epor(+/Cre)	50	55
Rb1(fl/fl);pl07(−/−);p130(+/ˆ);Epor(+/+)	47	41
Rb1(fl/fl);pl07(−/−);p130(+/ˆ);Epor(+/Cre)	47	50
Rb1(fl/fl);pl07(+/ˆ);p130(−/−);Epor(+/+)	18	15
Rb1(fl/fl);pl07(+/ˆ);p130(−/−);Epor(+/Cre)	19	22
Rb1(fl/fl);pl07(−/−);p130(−/−);Epor(+/+)	9	13
Rb1(fl/fl);pl07(−/−);p130(−/−);Epor(+/Cre) TKO	9	8
Total		250

The symbol ˆ indicates plus or minus.

Cell cycle withdrawal is a feature of terminally differentiating cells. Isolated FLC showed an increase in S phase cells in the absence of RB1, which was unchanged in the absence of all 3 pocket proteins (Figure S2A, TKO). FLC in S phase are poised to exit the cell cycle, so to monitor this process we cultured the cells for 48 hours. RB1-replete cells exited the cell cycle efficiently, regardless of the status of the other pocket proteins, whereas RB1-deficient cells arrested in S phase (Figure S2B). The TKO defect was equivalent to that of cells lacking RB1 alone. Similar results were obtained with FLC sorted for maturity (Figure S2C). Thus, p107 and p130 do not contribute to cell cycle withdrawal during terminal erythroid differentiation.

We confirmed the existence of an enucleation defect in the absence of RB1,¹⁹ and noted that it affects the definitive lineage based on morphology and is no worse in the absence of RB1 and p107 (Figure S3). Given the increase in cells with S phase DNA content, we examined the proliferation of RB1-deficient FLC in culture. RB1-deficient cells incorporated more Brdu initially than controls (54.4 ± 3.8% versus 48.0 ± 1.6%), which reversed itself by the end of the time course (18.9 ± 3.3% vs. 22.7 ± 2.5%). At the same time, RB1-deficient cells with >2N DNA content accumulated in the Brdu negative gate (Figure S4, arrow). These results suggest that deregulation of the G1-S transition causes an initial increase in Brdu incorporation, which is followed by the activation of S phase checkpoints, and decreased Brdu incorporation at later time points. Consistent with this interpretation, Brdu pulse-chase showed impaired S phase progression in the absence of RB1 (Figure S5).

RB1 deficiency is associated with large-scale changes in gene expression in erythroid cells

To extend existing concepts of RB1 function, we performed expression profiling of all genotypes. We did not observe major differences in erythroid subpopulations by surface marker expression (data not shown); therefore, we used unfractionated cells. We compared 34 fetal livers, employing the Signal-to-Noise (S2N) metric.²⁶ Significance levels were set by random permutation of columnar data (GSEA server), and the results presented in butterfly plots.²⁶ About 300 genes were upregulated in RB1-deficient fetal livers above a 1% threshold (Fig. 1A, barbed arrow). Expression of these genes increased a mean 1.4 ± 0.3 fold in the absence of RB1. By contrast, no genes were downregulated below the 1% threshold. Genes in p107- and p130-deficient fetal livers did not deviate significantly above or below the 1% thresholds. Furthermore, no effect was observed in TKO fetal livers, beyond that attributable to RB1 deficiency. We conclude that RB1 deficiency significantly upregulates a large number of genes, and that p107 and p130 do not have significant genome-scale effects on expression. However, loss of p107 did decrease expression of Serinc3, which was independent of RB1.

Figure 1. Genome-wide effects of RB1 deficiency. (A) Butterfly plots comparing microarray results from pocket protein deficient FLC to controls employing the Signal-to-noise (S2N) metric. The blue line corresponds to downregulated genes (the 300 most negative S2N values), the red line to upregulated genes (the 300 most positive S2N values), in rank order. The green, tan, and black lines represent 50%, 5%, and 1% thresholds, respectively, generated by permutation of randomized columnar data. The barbed arrow points a shoulder of several hundred genes expressed above the 1% threshold in the absence of RB1. A second arrow points to the position of the gene most repressed in the absence of p107, Serinc3. (B) Butterfly plots comparing microarray results from pocket protein deficient FLC to controls, in culture for 0, 24, and 48 hours. All time points were combined in calculating the S2N metric. The labels and scale are the same as in panel A.

RB1 repressed genes are well expressed but downregulated late in development

RB1 activity increases during erythroid differentiation; therefore, to examine the relationship between RB1-mediated repression and differentiation, we performed expression profiling on FLC in culture (0, 24, and 48 hours), with an additional 27 samples. Combining the time points, we confirmed genome-scale changes in expression, except, perhaps reflecting greater RB1 activity at later times, more genes were upregulated in the absence of RB1 (Fig. 1B, barbed arrow). Once again, p107 had little effect, except on Serinc3. Thus, our subsequent analysis was limited to RB1.

To characterize RB1 repressed genes we organized the microarray data by 2 independent criteria and merged the results graphically. First, the genes (32,170 probe sets reduced to 16,047 unique genes) were organized by hierarchical clustering, based on the expression pattern of the control samples at 0, 24, and 48 hours (5, 4, and 3 lanes of data, respectively), employing the Broad Institute algorithm (Gene Pattern).²⁷ Six clusters were assigned colors (gray for no pattern) with lighter shades representing regions of transition between clusters. Second, the data were divided into 16 bins (1,000 genes per bin), based on rank order, employing the S2N metric, with the genes most upregulated by loss of RB1 in bin 1 and those most downregulated in bin 16. Third, the genes in each bin were resorted according to their original hierarchical clustering position. Consequently, the genes are clustered by their native expression pattern along the vertical axis and by their response to the loss of RB1 along the horizontal axis (Fig. 2 and Supplemental File S1). This two-dimensional representation reveals several distinctive gene clusters. Of these, the genes in cluster C are most enriched in bin 1 (upregulated following loss of RB1), followed by those in clusters A and D.

Figure 2. Clustering analysis. (Top) Heat map of microarray results comparing RB1-deficient FLC in culture to controls. The results show the expression patterns of genes in control cells cultured for 0 hours (5 lanes), 24 hours (4 lanes), and 48 hours (3 lanes). Z-scores are gene normalized over the time course, using only the control samples. The data is divided into 16 bins. The genes are binned according to their regulation by RB1, employing the S2N metric (combining all time points), with those most derepressed in the absence of RB1 in bin 1 and those most repressed in bin 16. There are 1,000 unique genes per bin, except for bin 16, which has 1,047 genes. See supplemental file S1 for greater detail. (Bottom) Color coding and labeling of the major expression pattern clusters.

Next, we examined the expression patterns of the clusters in detail. The mean Z-scores of the genes in cluster C, and to a lesser extent clusters A and D, are downregulated between 24 and 48 hours (Fig. 3A and S6). By contrast, the genes in cluster B are downregulated in the initial 24 hours, and those in cluster F are upregulated. Bin 1 genes (clusters A, C, and D) were derepressed in the absence of RB1, although they retained a biphasic expression pattern. The expression ratio of this group of genes (sans gene normalization) at 48 hours compared to 24 hours was 0.61 ± 0.15 in the presence of RB1 and 0.87 ± 0.20 in the absence of RB1; thus, approximately two-thirds of the decrease in expression at 48 hours is due to RB1. Cluster E, a small group enriched in ribosomal genes, was not considered further. The median expression of the genes in clusters A, C, D, and F was similar, whereas that of cluster B was one-third as much. Present and absent calls provide an independent assessment of expression. Most of the genes in bin 1 were flagged present in all 27 microarrays (Fig. 3B). Consistent with their lower median expression, the genes in cluster B exhibited a mix of present and absent calls. The latter, however, occurred mostly at 24 and 48 hours, suggesting that the genes in this cluster were undergoing RB1-independent silencing. Based on these analyses, we conclude that RB1 downregulates hundreds of well-expressed genes late in erythroid development.

Figure 3. RB1-repressed genes are well expressed but downregulated late in development. (A) Charts show the mean Z-scores of genes in each of the major clusters in the control samples at 0 hours, 24 hours, and 48 hours. The chart at the lower right shows the Z-scores of the genes in bin 1 (clusters A, C, and D), in the control and RB1 deleted (KO) samples, at the same time points. Beneath the charts, the table shows the number of genes and relative median expression of all the genes in each cluster. (B) Bar graphs show the number of genes in each bin and cluster that satisfy the indicated condition. Present and absent calls are determined by Affymetrix software. The color schemes are consistent throughout. The bar graph at the lower right shows the percentage of genes in each cluster with 3 absent calls, at each time point.

Comparisons with external data sets

To assess the correlation between RB1-mediated repression, transcription factor binding, and histone modifications, we compared our results to several external ChIP-Seq and ChIP-chip data sets (Fig. 4). Alignment of our results with the top 2,000 E2F1 peaks in a breast cancer cell line,⁴⁰ showed that E2F1-bound genes are significantly enriched in bin 1. For E2F4-bound genes in primary murine tissues,⁴¹ or genes bound by DREAM complex components (of which E2F4 is one) in quiescent human glioblastoma cells,⁴² the degree of enrichment was even greater. KLF1 regulates cell cycle related and lineage marker genes in erythroid cells.^43,44 Consistent with this dual role, KLF1-bound genes were modestly enriched in bin 1; however, there was no such enrichment for genes that exhibited stage-specific KLF1 binding (data not shown). GATA1-bound genes,⁴⁵ by contrast were not significantly enriched.

Figure 4. Repression by RB1 correlates with E2F binding. Comparison to external ChIP-chip, ChIP-Seq, and alternative splicing data sets. The bar graphs show the number of genes in each bin and cluster bound by the protein indicated at the top. H3K27me3 (down) represents genes exhibiting downregulation of this histone mark at the Ter119- to Ter119+ transition. Color schemes are consistent with previous figures.

KDM5A is an RB1-interacting protein and histone lysine demethylase,^46-48 which inhibits differentiation.⁴⁹ KDM5A has been linked to the regulation of mitochondrial and cell cycle genes.⁵⁰ In our data set, genes belonging to these categories co-clustered by temporal expression pattern but differed in their response to RB1. Genes bound by KDM5A at the transcription start site⁵⁰ were not enriched in bin 1, suggesting that the response to RB1 may be dictated in part by this factor. A panel of histone marks has been examined for changes during erythroid development at the Ter119- to Ter119+ transition.⁵¹ Of these, only genes exhibiting downregulation of H3K27me3 showed evidence of enrichment (Fig. 4 and S7). However, RB1 recruitment may elicit histone modifications after the Ter119 transition, coincident with repression.

We also compared our results to 2 expression data sets. First, terminal erythroid differentiation is characterized by genome-scale changes in mRNA splicing,⁵² and we observed modest enrichment for alternatively-spliced genes in bin 1. Second, we compared our results to those from human CD34-derived erythroid cells, cultured and sorted by developmental stage.^28,29 We sorted the latter into 6 groups by k-means clustering (GenePattern). The genes in k-means cluster 1 were highly significantly enriched in bin 1 (Fig. 5). The expression of these genes increases progressively until the polychromatophilic erythroblast stage, then decreases. Of the 1201 human genes exhibiting this expression pattern, the orthologs of about one-third are RB1-repressed in mice. Whether RB1 is responsible for this expression pattern in human cells remains to be determined.

Figure 5. Murine orthologs of human genes downregulated in orthochromatophilic erythroblasts are repressed by RB1. The line charts show human RNA sequencing data gene normalized by Z-score and organized by k-means clustering into 6 clusters. Eight sorted human populations are shown: CD34+, BFU-E, CFU-E, proerythroblast (PRO), early basophilic erythroblast (E-BASO), late basophilic erythroblast (L-BASO), polychromatophilic erythroblast (POLY), and orthochromatophilic erythroblast (ORTHO). The bar graphs show the distribution of the murine orthologs of the genes in each cluster by bin. The color scheme is consistent with other figures.

Bioinformatics predictions of RB1 regulated processes

To gain insight into the role of RB1 in development, we assembled a list of RB1 repressed genes. To increase the stringency of our analysis we merged differential and time-dependent changes in gene expression; thus, we took the intersection of the genes in clusters A, C, or D, and bins 1 or 2. This approach yielded 801 RB1-repressed genes (group ACD (1-2)) (Supplemental File S2). Functional annotation clustering was performed on the DAVID bioinformatics server (Table 2).³⁰ RB1 repressed genes were enriched for terms involved in the cell cycle, DNA replication, and DNA repair. However, the terms could also be related to processes involved in terminal differentiation, such as chromatin condensation, mRNA splicing, nucleocytoplasmic transport, and modification-dependent protein catabolism. By comparison, the genes in cluster B were enriched for terms related to the extracellular matrix and cell adhesion, and those in cluster F were enriched for terms related to vesicular trafficking and apoptosis.

Table 2. DAVID functional annotation clustering

Pattern	Cluster	Cat	ES	GO:ID	Term	LH/PH	Benjamini
ACD (1-2)	1	BP	38.5	0007049	cell cycle	146/600	1.9E-57
	2	CC	34.2	0005694	chromosome	103/353	2.8E-59
801 genes	3	BP	23.7	0006259	DNA metabolic process	95/398	2.9E-35
	4	CC	14.0	0005654	nucleoplasm	71/587	7.7E-16
	5	BP	10.8	0006397	mRNA processing	47/258	8.8E-12
	6	CC	10.2	0015630	microtubule cytoskeleton	61/435	2.1E-16
	8	CC	4.5	0005643	nuclear pore	12/58	1.4E-04
	9	BP	4.2	0051276	chromosome organization	67/378	7.6E-17
	10	MF	4.1	0003887	DNA-directed DNA pol activity	10/30	6.5E-04
	12	CC	4.0	0005657	replication fork	11/24	1.5E-07
	13	BP	3.9	0050000	chromosome localization	6/8	7.6E-04
	14	BP	3.6	0019941	mod-dependent protein catabolic process	48/483	8.0E-04
ACD (4-16)	1	CC	28.6	0005730	nucleolus	104/304	9.5E-32
	2	BP	21.5	0034660	ncRNA metabolic process	75/198	5.5E-23
1599 genes	4	CC	12.9	0005743	mitochondrial inner membrane	69/286	5.8E-12
	5	BP	10.1	0006396	RNA processing	115/1099	2.2E-07
	6	BP	9.0	0006399	tRNA metabolic process	36/109	1.8E-08
	7	CC	8.4	0005761	mitochondrial ribosome	23/39	5.8E-12
	8	MF	6.8	0000166	nucleotide binding	265/2088	2.9E-06
	9	CC	6.2	0005839	proteasome core complex	11/19	3.3E-05
	10	MF	6.2	0003743	translation initiation factor activity	21/59	1.5E-05
	11	BP	5.7	0022900	electron transport chain	33/101	1.6E-07
B (1-16)	1	CC	9.4	0005789	endoplasmic reticulum membrane	51/145	1.9E-08
	2	CC	7.3	0005604	basement membrane	28/69	8.1E-06
1991 genes	3	MF	6.2	0030247	polysaccharide binding	38/122	3.5E-05
	4	BP	5.9	0016054	organic acid catabolic process	28/79	1.3E-04
	5	BP	5.7	0050817	coagulation	25/68	2.7E-04
	6	BP	5.7	0001568	blood vessel development	61/243	5.7E-05
	7	BP	5.3	0007155	cell adhesion	107/524	4.6E-05
	8	BP	4.5	0016125	sterol metabolic process	28/76	6.0E-05
	9	BP	4.2	0006487	protein amino acid N-linked glycosylation	14/24	2.6E-04
	10	MF	3.7	0004867	serine-type endopeptidase inhibitor activity	28/90	7.7E-04
	11	CC	3.7	0034358	plasma lipoprotein particle	13/27	3.3E-03
F'FF” (1-16)	1	CC	9.3	0005773	vacuole	63/198	1.5E-07
	2	BP	8.6	0006952	defense response	114/376	3.3E-10
2477 genes	3	BP	6.5	0007264	small GTPase med signal transduction	76/252	2.3E-06
	4	CC	6.4	0031982	vesicle	116/507	7.2E-06
	5	BP	6.2	0006915	apoptosis	110/455	1.1E-04
	6	BP	5.6	0050900	leukocyte migration	24/42	2.0E-06
	8	BP	5.2	0030029	actin filament-based process	55/175	4.3E-05
	10	BP	3.9	0001775	cell activation	77/232	1.8E-08
	11	MF	3.7	0017048	rho GTPase binding	13/24	9.2E-03

Pattern indicates hierarchical cluster pattern (included bins in parentheses); Cluster, DAVID cluster number; Cat, DAVID category; ES, enrichment score; GO:ID, gene ontology identification number; LH, list hits; PH, population hits; Benjamini, statistical significance corrected for multiple hypothesis testing.

Included in the mRNA processing term are members of the Ser/Arg-rich splicing factor family; collectively, these genes are repressed by RB1 (Fig. 6). These proteins regulate alternative splicing,⁵³ a process which is active during erythroid differentiation,⁵² and is implicated in hematological disease.^54,55 Other genes shared the biphasic expression pattern, but were unaffected or even activated by RB1 (ACD (4-16)). These were enriched for terms related to noncoding- and tRNA processing, translation, mitochondrial ribosomes, and the electron transport chain (Table 2). Within this group, NADH:ubiquinone oxidoreductase genes were downregulated in the absence of RB1 (Figure S8).

Figure 6. Ser/Arg-rich splicing factors are repressed by RB1 during terminal erythroid differentiation. Heat map of Ser/Arg-rich splicing factor expression in control and RB1-deficient FLC in culture at 0 hours, 24 hours, and 48 hours. Each row corresponds to a unique probe set. Z-scores are gene normalized over the time course (using control and RB1 deleted samples).

Discussion

In mice, germline mutation of Rb1 is associated with cell death in the nervous system, anemia, and embryonic lethality,^11-13 secondary to placental insufficiency.¹⁴ RB1 deficiency limited to the embryo proper causes unscheduled DNA replication in the brain, lens defects, and pituitary adenomas, but does not generally affect development.^15,16,56, RB1-deficient neonates have defects in skeletal muscle development,⁵⁷ due at least in part to a mitochondrial defect.⁵⁸ In the retina, RB1 is required for rod photoreceptor development;⁵⁹ retinoblastoma arises from cells that develop in the absence of RB1.⁶⁰ RB1 is also required for Purkinje cell¹⁵ and lens fiber development.⁶¹ In sum, RB1 deficiency causes specific developmental defects; however, RB1 is not generally required for the implementation of lineage-specific genetic programs.

The limited developmental defects caused by RB1 deficiency may be attributable to the existence of RB1-related proteins. The RB1 family is differentially expressed during development,⁶² and in several tissues p107 and p130 supply redundant function. These include skin,³¹ heart,³² muscle,⁶³ and embryonic central nervous system, lens, and blood vessels.³³ p107 and p130 can prevent immortalization,^38,39 and supply tumor suppressor activity in the retina,^34,35 and other tissues.³⁶ In contrast, p107 and p130 are not redundant with RB1 during terminal erythroid differentiation.

Regarding the lack of redundancy, it is instructive to consider how the functional role of the pocket proteins has evolved. Drosophila dE2F2-RBF complexes maintain genes in a repressed state during S phase; dE2F1 reverses this effect at cell cycle- but not developmentally-regulated genes.^64,65 This repressive effect is mediated through the dREAM complex (Drosophila, RBF, E2F, Myb-interacting proteins).⁶⁶ Human DREAM has a similar composition, with the role of RBF filled by p130.⁴² In contrast to the Drosophila complex, human DREAM is principally associated with cell cycle regulated genes. In mammals, the developmental role of dREAM could be met by RB1, the principal repressor of activator E2Fs, and most rapidly evolving member of the pocket protein family.⁶⁷ Regarding the unique requirement for RB1 in some tissues but not others, the answer may lie in the divergent activities and tissue-specific expression patterns of the E2F family.

The requirement for RB1 in erythroid cells could be related to increased E2F2 expression in late stage erythroblasts.²¹ E2F2 expression is preceded by E2F4; paradoxically, since E2F4 is considered a repressor, E2F4-deficient embryos are transiently anemic due to impaired FLC proliferation.^68-70 Once E2F2 is expressed, it becomes the principal RB1-associated E2F.²¹ E2F2 deficiency causes hypoproliferation of erythroid cells and erythropoietic abnormalities.^21,71 By contrast, E2F2 deficiency causes hyperproliferation of T cells and autoimmunity.⁷² Whether E2F2 activates or represses cell cycle regulators depends on the absence or presence of RB1, respectively.⁷³ Consistent with this interpretation, genetic experiments in mice indicate that the erythroid defects caused by RB1 deficiency are diminished by the loss of E2F2.²¹ This relationship may explain the lack of redundancy with p107 and p130, because they do not interact with E2F2 and therefore are unable to suppress its activity.

Depending on the spectrum of E2F proteins present, the other pocket proteins may play a role earlier in hematopoiesis. Indeed, deficiency of all 3 pocket proteins causes a myeloproliferative disorder.⁷⁴ It is interesting that all of the pocket proteins are able to establish quiescence, but only RB1 is able to establish senescence, in human fibroblasts.⁷⁵ Perhaps all the pocket proteins are capable of maintaining haematopoietic progenitor quiescence, whereas RB1 is uniquely required during terminal erythroid differentiation, which like senescence is characterized by changes in chromatin structure, and is normally an irreversible transition.

If RB1 supplies the required pocket protein activity during terminal erythroid differentiation, as our studies suggest, then why is RB1 deficiency associated with a mild erythroid defect in mice? One possibility is that differentiation does not require repression of E2F activity; however, E2F targets are still partially suppressed in RB1-deficient cells, suggesting that there is compensation. One mechanism of RB1-independent E2F repression is competition with repressor E2Fs for binding to E2F sites. E2F8 is an RB1-independent repressor; mice with combined deficiency of RB1 and E2F8 have hemolytic anemia,⁷⁶ which is caused at least in part by unrestrained E2F2 activity.⁷⁷ In Drosophila, S phase entry triggers ubiquitin-dependent degradation of dE2F1.⁷⁸ In mammals, phosphorylation of E2F-bound DP by cyclin A-Cdk2 prevents E2F-DP binding to DNA.⁷⁹ Thus, redundant mechanisms exist that mitigate the effect of RB1 deficiency on E2F target gene expression.

Loss of p107 and p130 did not have a broad effect on gene transcription. However, p107 deficiency did affect expression of the gene Serinc3. Serinc3 is a member of a family of proteins that incorporate serine into membranes for phosphatidylserine and sphingolipid synthesis.⁸⁰ Notably, RB1 had no effect on Serinc3 expression. Although p107 selectively represses certain genes,⁸¹ this is an example of p107-specific gene activation, which has not been previously described. Interestingly, the Serinc3 locus (Tde1) is a target of insertional mutagenesis in Arf-null mice.⁸² Whether there is crosstalk between p107 and the p53 pathway, as this observation implies, is unknown.

By contrast, RB1 deficiency has a major effect on erythroid gene transcription, causing derepression of hundreds of genes. At the single gene level, the magnitude of the effect is modest but statistically significant. Compounded over hundreds of genes, there are likely to be cellular consequences; this point is illustrated by the cell cycle defects present in RB1-deficient erythroid cells. From a network perspective, it is logical to think of RB1 as a horizontally-applied brake on E2F-driven biological processes; along this line, an accurate compilation of RB1-repressed genes should improve our understanding of the genes and processes involved in terminal differentiation.

RB1 recruits chromatin modifying enzymes to genes through its interaction with E2Fs.⁸³ In the absence of RB1, E2Fs are effectively converted from repressors to activators, derepressing E2F targets. In this regard, E2Fs and RB1 may act as a switch during development, first activating then repressing E2F targets.⁷³ During terminal erythroid differentiation, Cdk4/6-associated RB1 kinase activity decreases, and hypophosphorylated RB1 increases,³⁷ suggesting the switch is linked to a change in the phosphorylation state of RB1. Comparison of our results with those from human erythroid cells ^28,29 further suggests that an RB1-regulated switch may be a conserved aspect of development.

Repression by RB1 correlates with the presence of bound E2F1, E2F4, and DREAM complex components. In Drosophila, dE2F-dDP complexes are required for RBF recruitment.⁸⁴ In human senescent fibroblasts, RB1 recruitment correlates with the presence of E2F binding sites.⁷⁵ In murine FLC, DP2 knockdown interferes with downregulation of E2F2 target genes.⁸⁵ Together, these findings support the idea that E2Fs mediate RB1 recruitment and gene repression in erythroid cells.

KLF1 appears to have a dedicated role regulating the expression of RB1-repressed genes, which is not limited to direct targets, since E2F2 itself is regulated by KLF1.^86,87 GATA1-bound genes were not significantly enriched for RB1-repressed genes, but GATA1 interacts with RB1 and positively regulates its activity;⁸⁸ thus, depending on the relative abundance of E2Fs and RB1, KLF1 and GATA1 could have reinforcing or opposing effects on the expression of E2F target genes.

By merging differential and temporal gene expression patterns, we identified 801 RB1-repressed genes in erythroid cells. Approximately one-third of the genes we identified are bound by E2Fs in other cell types.^40-42 Cell cycle related genes and confirmed E2F targets are enriched in the list;⁸⁹ however, so are genes and processes not known to be regulated by E2F and RB1. Because they share differential and temporal gene expression patterns, we anticipate that many of the novel genes will also prove to be direct targets of RB1-mediated repression in erythroid cells, but this will require validation. In the interim, we expect that this gene list will be a useful resource for studies of the mechanisms of terminal erythroid differentiation.

The large number of RB1-repressed genes is consistent with a developmental role for RB1, possibly one that involves coordinating terminal differentiation with cell cycle withdrawal. A recent survey found that 2,250 genes are alternative spliced during terminal erythroid differentiation.⁵² Splicing events occurred in late stage erythroblasts, concurrent with extensive cellular remodeling. RB1 repressed Ser/Arg-rich splicing factors, which are key regulators of alternative splicing.⁵³ Ser/Arg-rich splicing factors are themselves inactivated by a separate mechanism, alternative splicing-nonsense mediated decay, which further argues that their downregulation is important for terminal differentiation.⁵²

RB1 repressed genes involved in transcription and chromatin structure, such as GFI1B, EZH2, and SUZ12. The PRC2 subunits EZH1 and EZH2 undergo an expression switch during terminal erythroid differentiation, which is regulated by GATA2 and GATA1.⁹⁰ Repression by RB1 adds another dimension to their regulation. Ubiquitin-dependent catabolism is another process represented on the list; KCTD6 is one of several Cullin 3 adaptors identified and is implicated in ankyrin protein turnover and hereditary spherocytosis.⁹¹

Genes involved in mitochondrial biogenesis and respiration are downregulated in the absence of RB1.²⁰ We showed this effect persists throughout terminal differentiation, and extended it to mitochondrial ribosomes, noncoding RNA metabolism, and translation. Regarding the role of RB1 in oxidative phosphorylation, the finding that mitochondrial biogenesis is diminished in RB1-deficient muscle cells is concordant with our results.⁵⁸ Conversely, findings that oxidative phosphorylation is repressed by RB1,⁹² and degradation-defective cyclin E increases mitochondrial biogenesis are not.⁹³ Potentially, RB1 activates these genes by reversing KDM5A-mediated repression.⁵⁰ Additional experiments are needed, since the metabolic derangements associated with RB1 deficiency are related to those present in cancer cells.⁹⁴

RB1 regulates the expression of genes involved in the cell cycle, chromatin organization, alternative splicing, transcription, and ubiquitination. These processes set the stage for the final phase of erythroid development, which includes cellular remodeling, enucleation, and the expression of erythroid-specific protein isotypes. We anticipate that processes like these will be regulated by pocket proteins in a wide range of tissue types, but that the specific genes may differ. For some combinations of RB1 targets and tissues, down-regulation may be essential for terminal differentiation to proceed to completion. Identification of these genes will yield insights into development, and possibly into the causes of cancer associated with mutations in the RB1 pathway.

Disclosure of Potential Conflicts of Interest

No potential conflicts of interest were disclosed.

Supplemental Material

Supplemental data for this article can be accessed on the publisher's website.

Acknowledgments

The authors acknowledge Drs. Anton Berns, Tyler Jacks, Ursula Klingmüller, and Michael Dyer for the conditional Rb1, p107 and p130 null, and Epor-cre mouse strains. The authors thank Drs. Mohandas Narla, Xiuli An, and David Bodine for sharing data. We thank the Laboratory Animal Research Service of the New York Blood Center. We also thank the Animal Resource Center, and the Flow Cytometry and Cell Sorting Shared Resource of St. Jude Children's Research Hospital.

Funding

This work was supported by grants from the American Society of Hematology and the F.M. Kirby Foundation (PAN), the New York Blood Center, and the American, Lebanese, and Syrian Associated Charities (ALSAC).