RNA sequencing of MCF-7 breast cancer cells identifies novel estrogen-responsive genes with functional estrogen receptor-binding sites in the vicinity of their transcription start sites

Abstract

Estrogen receptor α (ERα) is a key transcription factor in breast cancer, which plays an essential role in the pathophysiology of the disease by regulating the expression of various target genes. In the present study, we performed deep RNA sequencing (RNA-seq) as an unbiased high-throughput technique for comprehensive transcriptome analysis in ERα-positive human breast cancer MCF-7 cells, to facilitate the elucidation of ERα regulatory gene networks. From the 17,336 mapped RefSeq genes from the sequenced fragments of the cell samples treated with estrogen time dependently, substantial numbers of sequence reads were observed in 3,386 genes (>100 tags per million reads per sample at any of the six time points studied). ERα occupancy within and in the proximal regions of the genes (<10-kb upstream and downstream regions) was significantly enriched in the subgroup of the 3,386 genes compared to the whole 17,336 RefSeq genes. Of the 3,386 genes, we focused on 29 genes, which included ERα occupancy adjacent to their transcription start sites and whose expression was estrogen dependently altered by >3-fold. Knockdown studies using siRNAs specific to the 29 genes validated that prototypic ERα targets V-myc myelocytomatosis viral oncogene homolog and cyclin D1 promote both proliferation and migration of MCF-7 cells and further identified novel candidate ERα targets EIF3A and tumor protein D52-like 1, which will also facilitate the proliferation or migration of MCF-7 cells. Taken together, the present findings provide a valuable dataset that will elucidate ERα regulatory mechanisms in breast cancer biology, based on the integrative analysis of RNA-seq combined with the genome-wide information for ERα occupancy.

Fig. 1

Alteration in the gene expression level by 17β-estradiol (E₂). a Annotation of RNA fragments at each time point after E₂ treatment. RNA-seq was performed using RNAs prepared from MCF-7 cells that were treated with E₂ for 0, 2, 4, 8, 12, and 24 h. Of all mapped fragments, 6.73 × 10⁷ (73.9 %) fragments were mapped onto exons, 11.5 × 10⁷ (12.6 %) were mapped onto introns, and 1.22 × 10⁷ (13.4 %) were mapped onto the regions where no RefSeq genes were detected (Intergene). b Venn diagram of E₂-regulated RefSeq genes (fold change, >2). Genes whose expression levels are low (less than 1 tpm at any time point) were excluded and, of the 3,439 genes, 776 (22.6 %) were up-regulated and 92 (2.7 %) were down-regulated by E₂. Only one gene was both up- and down-regulated. c GREB1, a known estrogen-responsive gene, was found to be up-regulated by E₂ treatment. Black bars represent reported estrogen receptor-binding sites (ERBSs) [13], which are located within 10 kb from GREB1. d PGR, a known estrogen-responsive gene, was also found to be up-regulated by E₂ treatment. Black bars represent reported ERBSs located within 10 kb from PGR

Fig. 2

Estrogen-regulated genes validated by qRT-PCR. Of the 869 genes that were detected as candidate estrogen-responsive genes in RNA-seq analysis, 29 were located within 10 kb from reported ERBSs and validated for their estrogen responsiveness by qRT-PCR. Numbers represent tpm at each time point, and colors represent relative expression levels compared to those at 0 h in RNA-seq

Fig. 3

Locational relationship between estrogen-responsive genes and ERBS. a An ERBS was located in the upstream of MLL5. Lower panel shows the result of RNA-seq, indicating that MLL5 is up-regulated by E₂. b An ERBS was located in the second exon of the FAM84B. Lower panel shows the result of RNA-seq, indicating that FAM84B is up-regulated by E₂. c An ERBS was located in the first intron of insulin-like growth factor binding protein 4 (IGFBP4). Lower panel shows the result of RNA-seq, indicating that IGFBP4 is up-regulated by E₂. d ERBSs were located in the upstream, fourth intron, and downstream of cyclin D1 (CCND1). Lower panel shows the result of RNA-seq, indicating that CCND1 is up-regulated by E₂

Fig. 4

Knockdown analysis of novel estrogen-responsive genes by siRNA. a Knockdown effect of each siRNA against targeting genes. MCF-7 cells were transfected with 5 nM siRNA that target the E₂-reponsive candidate genes and siControl with Lipofectamine 2000 for 48 h. Gene expression levels with siRNA were detected by qRT-PCR and normalized with siControl. b Assessment of cell growth using MTS assay. Cell viability was determined using MTS assay. MCF-7 cells were transfected with 5 nM siRNA specific for the indicated genes and control siRNA for 4 days. Results are represented as mean ± SD of three experiments. Fold change of absorbance at 490 nm of each siRNA sample was normalized with that of siControl sample. c Trans-well migration assay. MCF-7 cells transfected with each siRNA were incubated for 24 h and plated to uncoated filters for 48 h. The migrating cells were fixed, stained, and counted for three fields of the filter. d Assay for ER-mediated transcription. ERE-luciferase reporter assay. ERE-luciferase reporter plasmid and 5 nM siRNA were transfected in MCF-7 cells with or without 100 nM E₂ for 24 h. Fold change of luciferase activity in response to E₂ was calculated for each siRNA. Student’s t test was performed between each siRNA sample and siControl sample (*P < 0.05, **P < 0.01, ***P < 0.005)

Fig. 5

EIF3A is a novel estrogen-responsive gene that affects the growth and migration of MCF-7 cells. a Genome view shows RNA fragments mapped on EIF3A. b Quantification of mapped fragments onto EIF3A indicated the E₂ responsiveness. c qRT-PCR confirmed the E₂ responsiveness of EIF3A. d MTS assay indicated that the growth of MCF-7 cells was decreased by the knockdown of EIF3A. e Trans-well assay indicated that the migration of MCF-7 cells was decreased by the knockdown of EIF3A. Student’s t test was performed between each siRNA sample and siControl sample (**P < 0.01)

Fig. 6

TPD52L1 is also a novel estrogen-responsive gene that affects the migration of MCF-7 cells. a Genome view shows RNA fragments mapped on TPD52L1. b Quantification of mapped fragments onto TPD52L1 indicated the E₂ responsiveness. c qRT-PCR confirmed the E₂ responsiveness of TPD52L1. d MTS assay indicated that the growth of MCF-7 cells was not affected by the knockdown of TPD52L1. e Trans-well assay indicated that the migration of MCF-7 cells was decreased by the knockdown of TPD52L1. Student’s t test was performed between each siRNA sample and siControl sample (***P < 0.005)