Introduction

Bromodomain-containing 7 (BRD7) is a tumor suppressor gene that was first cloned through cDNA representational difference analysis (RDA) and is down-regulated in nasopharyngeal carcinoma (NPC) [1–5]. Accumulating evidence indicates that BRD7 was also down-regulated in multiple types of human cancer, including breast [6], prostate [7], endometrial [8], colorectal [9], ovarian [10], and pancreas cancer [11], glioma [12], and osteosarcoma [13]. BRD7 can inhibit cell growth through multiple mechanisms, including cell cycle arrest and apoptosis [14, 15]. Recent work suggests that BRD7 acts as a transcriptional cofactor binding to BRCA-1 and p53 and is essential for the transcriptional activation of p53 target genes such as p21, MDM2, and TIGAR [16–18]. BRD7 is also a subunit of the SWI/SNF complex specific for polybromo BRG1-associated factor (PBAF) [19], modulates chromatin remodeling by binding to acetylated histone H3 [20, 21] to activate/inhibit some BRD7-dependent transcription [17, 22], and plays a critical role as a transcription regulator in several important tumor suppressor pathways [22–25]. However, given the divergent cellular roles of the BRD7 gene, the BRD7 target gene and its downstream gene regulating network have not yet been explored.

Chromatin immunoprecipitation (ChIP) followed by large-scale DNA analysis such as DNA-chip (ChIP-chip) or high-throughput sequencing (ChIP-seq) is a powerful way to identify the directly binding loci of a transcriptional factor throughout the genome [26, 27]. Digital gene expression (DGE) is a method that generates a digital output proportional to the number of transcripts per RNA by partially sequencing millions of randomly selected cDNA tags from relevant libraries. Differentially expressed genes can then be detected from variations in the counts of their cognate sequence tags. This method has the benefit of not requiring presynthesized oligonucleotide probes (as in microarrays), and allowing the direct enumeration of transcript molecules, which is directly comparable across different experiments [28]. In this study, we combine ChIP-seq and digital gene expression (DGE) profiling by RNA-sequencing upon the overexpression of BRD7 in human cells to identify BRD7-binding regions and its regulated genes. We identified four novel BRD7 specifically recognized motifs within BRD7-binding site peaks. We also identified a core set of 184 genes as BRD7 target genes that have BRD7-binding sites. We preliminarily constructed interaction networks of the differentially expressed genes regulated by BRD7. Our study provides an important resource for understanding the function of BRD7, especially its transcriptional activity.

Materials and methods

Cell lines, plasmids, and transfection

HEK293 and HeLa cells were obtained from the American Type Culture Collection (ATCC, Rockville, MD) and were maintained in Dulbecco’s modified Eagle’s medium (DMEM) supplemented with 10 % fetal bovine serum (FBS), 100 U/ml penicillin, and 100 μg/ml streptomycin at 37 °C constant temperature and 5 % CO2.

To perform the specificity ChIP assay, we constructed an over-expression BRD7 vector that has 3 tandem Flag-tags. First, the full-length open reading frame (ORF) of the human BRD7 gene was cloned using the following PCR primers: 5′-CGGCTAGCCAATGGGCAAGAAGCACACAAGA-3′ and 5′-GAAGGCCTCTCAACTTCCACCAGGTC-3′. Then, 3 tandem Flag-tags were synthesized and added to the 5′ end of the BRD7 ORF and cloned into the pIRESneo3 vector. For the control vector, only 3 tandem Flag-tags were added to the empty pIRESneo3 vector. Lipofectamine® 2000 reagent (Life Technologies, USA) was used for transfection cells.

qRT-PCR

Total RNA was extracted and purified using TRIzol (Invitrogen) according to the manufacturer’s instructions. The purified RNA was reverse transcribed using AMV reverse transcriptase (Promega, USA). Quantitative real-time PCR was performed using SYBR® Premix Ex Taq™ (TaKaRa, Dalian, China). The primers used for gene expression analysis were as follows: BRD7, 5′-TGGAGATGTCATTGCCTGAAGA-3′ and 5′-CCCTGGTGGCTCTACTTCTG-3′; BIRC2, 5′-TGGAGATAGGGTAGCCTGCTT-3′ and 5′-GGAAAATGCCTCCGGTGTTC-3′; BIRC3, 5′-TGATGAAAAGCGCCAACACG -3′ and 5′-AGAAACCAGCACGAGCAAGA-3′; TXN2, 5′-CTGCCTGTACCCGGAAGTGA-3′ and 5′-CTGAGCCATCTCCCTGCAATG-′3; NOTCH1, 5′-CCACTGTGAGAACAACACGC-3′ and 5′-CACAAGGTTCTGGCAGTTGG-3′; and GAPDH, 5′-TCTAGACGGCAGGTCAGGTCCACC-3′ and 5′-CCACCCATGGCAAATTCCATGGCA-3′. GAPDH was chosen as the reference gene for the normalization of all gene expression results. The average of three independent analyses for each gene was calculated. The fold changes were calculated through relative quantification \(( 2^{{ - \varDelta \varDelta C_{\text{t}} }} )\) [29]. All reactions were run in triplicate and repeated in three independent experiments.

Western blotting analysis

Western blotting was performed as described previously [30–32]. Cell lysates were separated by electrophoresis on 10 % sodium dodecyl sulfate (SDS) polyacrylamide gels, and the separated proteins were transferred to a polyvinylidene fluoride (PVDF) membrane (Millipore, Billerica, MA). The immunoblots were incubated sequentially with the primary antibodies mouse anti-FLAG M2 monoclonal antibody (Sigma-Aldrich, St. Louis, USA) and anti-GAPDH (Abcam corporation, Cambridge, UK), followed by horseradish peroxidase-coupled secondary antibodies. Signals were generated by chemiluminescence with ECL substrate reagent (Pierce, Rockford, USA).

Immunofluorescence

Cells were fixed in 4 % paraformaldehyde for 20 min, permeabilized with 0.5 % Triton X-100 for 3 min, and blocked with phosphate-buffered saline (PBS) containing 7 % fetal bovine serum for 30 min. Cells were incubated with mouse anti-FLAG M2 monoclonal antibody (Sigma-Aldrich, St. Louis, USA) at 4 °C for 90 min. The cells were further reacted with the FITC-labeled secondary antibody goat anti-mouse IgG (1:200, Gibco) for 90 min at room temperature, followed by nuclear staining using Hoechst 33258 for 5 min. The sections were washed 3 times with PBS after each incubation. The AX-80 analytical microscope system (Olympus, Tokyo, Japan) was used for observation [30, 31].

MTT assay

The cells were seeded at 1.0 × 103 cells per well into 96-well plates. The cells were incubated with MTT [3-(4, 5-dimethylthiazolyl-2)-2,5-diphenyltetrazolium bromide] reagent (Sigma-Aldrich, St. Louis, USA) at 0.5 mg/ml. The medium was removed, and formazan crystals were solubilized in isopropanol. Absorbance was measured at 570 nm. The assays were performed in triplicate.

Flow cytometric analysis

To measure the cell cycle distributions, the cells were fixed in 70 % (v/v) ethanol at 4 °C overnight. Then, the fixed cells were digested with 100 mg/ml of RNase A in PBS and stained for 30 min with 50 mg/ml of propidium iodide (Sigma-Aldrich, St. Louis, USA) in the dark [33].

To measure apoptosis, the cells were washed with cold PBS and resuspended in 500 µl of annexin V binding buffer (PE annexin V apoptosis detection kit; BD Biosciences, San Diego, CA) containing phycoerythrin (PE)-labeled annexin V and 7-amino-actinomycin (7-AAD) and were incubated for 10 min at room temperature in an area shielded from light.

The specimens were analyzed by fluorescence-activated cell sorting (FACS) using a FACSCalibur apparatus (BD Biosciences), acquiring 10,000 events.

TUNEL assay

The TUNEL assay was also used to observe apoptosis. The cells were fixed in 4.0 % paraformaldehyde for 20 min, permeabilized with 0.5 % Triton X-100 for 3 min, and then subjected to the TUNEL assay by the Dead End™ Fluorometric TUNEL system (Promega, Madison, WI) according to the manufacturer’s specifications. The cells were then incubated simultaneously with Hoechst 33258 for nuclear staining for 15 min at 37 °C. The sections were washed 3 times with PBS after each incubation. The AX-80 analytical microscope system (Olympus, Tokyo, Japan) was used for observation.

Chromatin immunoprecipitation (ChIP)

The ChIP assay was performed using an EZ-ChIP kit (Upstate Biotechnology, USA). Briefly, the cells were fixed with 1 % paraformaldehyde in flasks at room temperature for 10 min and subsequently washed with ice-cold PBS. Then, the cells were lysed with 50 mM Tris HCl, at a pH of 7.4, 150 mM NaCl, 1 mM EDTA, and 1 % Triton X-100. The lysate was incubated with anti-FLAG M2 monoclonal antibody (Sigma-Aldrich, St. Louis, USA) at 4 °C overnight. The DNA was purified with a QIAquick DNA purification kit (Qiagen, Duesseldorf, Germany) according to the manufacturer’s instructions.

ChIP-seq analysis

All standard protocols for Illumina sequence preparation, sequencing, and quality control were performed by BGI (Shenzhen, China). In brief, DNA fragments recovered from a conventional ChIP procedure were quantified and the integrity was verified, followed by end repair, adaptor ligation, and size selection to construct libraries of the BRD7 overexpression and control groups. The DNA libraries were validated and sequenced using the Illumina HiSeq 2000 sequencer, and approximately 30 million sequencing reads were obtained per sample. The ChIP-seq data have been deposited in the Gene Expression Omnibus (GEO) database under accession number GSE65981. After being aligned to the human reference genome Hg19, sequencing tags were processed using MACS (BGI, China) to perform peak scanning, and the peaks were mapped using the UCSC Genome Browser (http://genome.ucsc.edu/). MEME (Motif-based sequence analysis tools, http://tools.genouest.org/tools/meme/) was used for de novo motif exploring, and MAST (Motif Alignment & Search Tool, http://tools.genouest.org/tools/meme/cgi-bin/mast.cgi) was used for motif alignment to exclude duplicate motifs.

ChIP-PCR validation

The ChIP-seq data were validated by ChIP-PCR. Briefly, ChIP was performed in separate samples of BRD7-overexpressing or control cells, including HEK293 and HeLa, and then the BRD7 target genes identified by ChIP-seq were validated by regular PCR. The validated genes and their PCR primers were as follows: BIRC2, 5′-CACAGTGGTTGCTGTCAAGG-3′ and 5′-CTGGCCCTTTAACCCTGTCT-3′; BIRC3, 5′-ATGACATTTTAGGGACATGGTGTT-3′ and 5′-ATGTGTTGACCCCGTGTTCT-3′; TXN2, 5′-TACCTACGTGGAGGCCACTC-3′ and 5′-ACCCACGCTTCAGGTTCTAC-′3; NOTCH1, 5′-GGAATCATATTCCCTCGTGCG-3′ and 5′-CAACGACGGGTCAGCAAA-′3; and GAPDH, 5′-TCTAGACGGCAGGTCAGGTCCACC-3′ and 5′-CCACCCATGGCAAATTCCATGGCA-3′. GAPDH was chosen as a control. The PCR results were analyzed by 2 % agarose gel electrophoresis.

Digital gene expression (DGE) profiling

The RNA-sequencing-based DGE was performed by the BGI Company. Briefly, mRNA was purified from the BRD7-overexpressing and control cells and fragmented for first-strand cDNA synthesis with a random primer by reverse transcriptase. After synthesis, second-strand cDNA synthesis was performed using DNA polymerase I and RNase H and end repair, and the cDNA fragments were ligated with adapters. These products were purified and amplified by PCR to create the final cDNA library. The cDNA libraries were then sequenced on an Illumina HiSeq 2000 system. Clean reads were aligned to the Hg19 human reference genome sequence and were annotated with specific genes. The normalized expression of every gene was measured using tags per million reads (TPM). The DGE data have been deposited in the Gene Expression Omnibus (GEO) database (Accession number: GSE53656). Genes with an absolute expression ratio > 2 between the BRD7-overexpressing and control cells and a false discovery rate (FDR) < 0.001 were determined to be significantly differential expression genes.

Bioinformatic analysis

To understand the biological significance of genes identified by ChIP-sequencing and DGE profiling, an interaction network analysis [34] was performed with ingenuity® pathway analysis (IPA, http://www.ingenuity.com/) [35, 36]. The database for annotation, visualization, and integrated discovery (DAVID 6.7, http://david.abcc.ncifcrf.gov/) [37–39] was used for functional clustering and pathway analysis. To visualize the integrated network of the ChIP-sequencing and DGE data, Cytoscape was used to map the interaction network based on the public gene/protein interaction network database STRING.

Statistical analysis

The data were analyzed by two-way repeated-measures analysis of variance (ANOVA). P values < 0.05 were considered statistically significant. All results were expressed as the mean ± standard deviation (SD). Curve fitting and analyses were performed using GraphPad Prism (GraphPad Software, San Diego, CA).

Results

BRD7 inhibits cell proliferation and induces apoptosis in HEK293 and HeLa cells

Prior to ChIP-sequencing, we generated two cell lines (HEK293 and HeLa) that consistently expressed the BRD7 gene (HEK293-BRD7 and HeLa-BRD7) and two corresponding control cells transfected by control vector (HEK293-Vector and HeLa-Vector) as the experimental model. The over-expression of BRD7 was validated by qRT-PCR (Fig. 1a) and western blotting (Fig. 1b). Immunofluorescence also showed that ectopically expressed BRD7 protein was mainly localized in the nuclei of HEK293 and HeLa cells (Fig. 1c).

Fig. 1
figure 1

Stable over-expression of the BRD7 gene in HEK293 and HeLa cells. a Real-time PCR validated the successful over-expression of BRD7 in HEK293 and HeLa cells. Vector indicates cells transfected with empty vector, and BRD7 indicates cells transfected with BRD7 over-expression plasmid that have 3 tandem Flag-tags fused to the 5′ terminal of the BRD7 ORF sequence, **p < 0.01. b Western blotting validated the overexpression of BRD7 at the protein level in HEK293 and HeLa cells. To evaluate the over-expression plasmid transcript and to translate the flag-BRD7 fusion protein, flag antibody was used to detect ectopically expressed BRD7. c Immunofluorescence staining confirmed that ectopically expressed BRD7 was mainly located in the nucleus of HEK293 and HeLa cells. Images were acquired at 200×. Scale bar 50 μm

Full size image

The MTT assays showed that BRD7 inhibits cell growth and proliferation (Fig. 2a). Flow cytometric evaluation displayed that stable BRD7-overexpressing cells had an increased proportion of G1 and G2 phases and a decreased proportion of the S phase relative to control cells (Fig. 2b). Both flow cytometry (Fig. 2c) and TUNEL staining (Fig. 2d) demonstrated that the stable overexpression of BRD7 increased apoptotic cell populations compared to the controls. These results confirmed that BRD7 acts as a tumor suppressor gene and can inhibit cell growth and induce apoptosis in either the human embryonic kidney cell line HEK293 or the human cancer cell line HeLa. We noticed that overexpression of BRD7 in HEK293 arrested cell cycle progression significantly at both G1 and G2/M phases, while only G2/M arrest is significant in Hela cells (though proportion of G1 cells was also increased); we thought that BRD7 may have slightly different functions in different cell types and environments.

Fig. 2
figure 2

Overexpression of BRD7 inhibits cell proliferation and induces apoptosis. a Cell viability assays measured by MTT indicated that BRD7 obviously inhibited the cell growth in either HEK293 or HeLa cells. b The cell cycle distribution examined by flow cytometry showed that BRD7 decreased the proportion of cells in S phase and increased the proportion of cells in G1 and G2 phases. c Apoptosis assays indicated that BRD7 induced cell apoptosis. d TUNEL staining also confirmed that BRD7 induced cell apoptosis. Images were acquired at 100×. Scale bar 50 μm. *p < 0.05; **p < 0.01;***p < 0.001

Full size image

Motif analysis of the BRD7-binding peaks of ChIP-seq

Because the stable BRD7 expression level was higher in HEK293 cells than in HeLa cells, we chose HEK293 cells that stably transfected BRD7 to perform ChIP-seq and DGE analysis, and HeLa cells were used for validation. ChIP-seq yielded approximately 30 million individual sequencing reads in each run (raw data). Only reads that uniquely mapped onto the human genome were used for the subsequent analysis (Table 1). We mapped a total of 22,546,083 and 21,961,299 sequence tags uniquely to the human genome for BRD7-overexpressing and control cells, respectively. A total of 156 BRD7-binding site peaks were identified in the human genome (supplemental Table S1).

Table 1 Summary of ChIP-seq clean data mapping to the human reference genome
Full size table

Using multiple EM for motif elicitation (MEME), a de novo motif discovery program [40], we analyzed the consensus motif present in BRD7-binding peaks. We identified 4 novel motifs (Fig. 3a). Among these 4 motifs, motif 1 was found within all binding events (E value 1.5E−46), and the other 3 motifs, namely motifs 2, 3, and 4, were identified within 24.36 % (E value 1.9E−28), 57.05 % (E value 8E−18), and 26.28 % (E value 2.1E−13) of BRD7-binding sites, respectively.

Fig. 3
figure 3

The BRD7-binding peaks were co-localized with histone modification sites. a Four novel BRD7-binding motifs identified within BRD7-enriched ChIP-Seq peaks. b Representative BRD7-binding peaks and their coordinate histone modification sites located in BIRC2, BIRC3, TXN2, and NOTCH1 gene loci. H3K4me1, H3K4me3, and H3K27ac modification data from the Encyclopedia of DNA Elements (ENCODE) project were integrated in the UCSC Genome browser, and BRD7-binding peaks from our ChIP-seq data were also mapped to the UCSC Genome browser. c These BRD7-binding sites were confirmed by ChIP-PCR in an independent ChIP assay, and GAPDH served as the negative control

Full size image

The BRD7-binding peaks were co-localized with histone modification sites

Using the UCSC Genome browser, we visualized 156 BRD7-binding peaks; these peaks were located in or adjacent to 20 kb of 184 human genes (supplemental Table S2). BRD7-binding peaks were co-localized with histone modification sites, such as the acetylation of histone 3 at lysine 27 (H3K27) and the methylation of histone 3 at lysine 4 (H3K4). Representative BRD7-binding peaks and their coordinate histone modification sites located in BIRC2, BIRC3, TXN2, and NOTCH1 gene loci are shown in Fig. 3b. These BRD7-binding sites were also confirmed by ChIP-PCR (Fig. 3c).

Construction of a regulating network of BRD7 target genes

Because 156 BRD7-binding peaks identified by ChIP-seq were located in or adjacent to 20 kb of 184 human genes, these 184 genes may be directly targeted by BRD7 and may be transcribed by BRD7. We annotated and analyzed the biological function and the interactions within the network of these genes by IPA. A network annotated as “Cell Death and Survival” was constructed (Fig. 4). This network involved 22 genes playing important roles in cell cycle and apoptosis regulation, such as members of the inhibitor of apoptosis protein (IAP) family, BIRC2 and BIRC3, which participate in inhibiting apoptosis through caspases, ubiquitin and Akt signaling pathways, and NOTCH1, a transmembrane receptor, which may play multiple roles during development and carcinogenesis by controlling cell fate decisions through its downstream network.

Fig. 4
figure 4

The interactive network of BRD7 target genes. A total of 184 BRD7 target genes identified by ChIP-seq were annotated and analyzed according to their biological function by ingenuity pathway analysis (IA), and a network annotated as “Cell Death and Survival” was constructed. Genes marked in red are BRD7 target genes, and others are related genes that have no BRD7-binding peaks. (Color figure online)

Full size image

DGE profile and functional annotation of BRD7-regulated genes

To explore BRD7 downstream genes, including directly and indirectly regulated genes, gene expression levels in stably transfected BRD7 and control HEK293 cells were profiled by RNA-sequencing-based DGE. There were 6,646,235 and 6,621,819 clean reads that uniquely mapped to specific human genes in BRD7 stably over-expressing and control cells, respectively. Gene expression profiling was measured and normalized by tags per million reads (TPM), and 1648 genes were significantly differentially expressed, at 2.0-fold or greater (FDR < 0.001); in the BRD7-overexpressing cells compared to the vector control, 560 genes were up-regulated and 1088 were down-regulated (supplemental Table S3).

To identify the biological pathways regulated by BRD7, we classified these dysregulated genes based on enriched Gene Ontology (GO) terms such as the cellular component, molecular function, and biological processes, and we mapped these genes onto the KEGG canonical pathways using the DAVID program. The apoptosis pathway has the greatest number of genes, followed by RNA metabolism, cell cycle, cofactor metabolism, and chemotaxis (Fig. 5a). Figure 5b illustrates dysregulated genes in the apoptosis pathway. Interestingly, 104 non-coding RNA sequences were dysregulated by BRD7, 33 were up-regulated (such as LINC00472 [41], GNAS-AS1 [42], OIP5-AS1, etc.), and 71 were down-regulated (such as MALAT1 [43–46], CDKN2B-AS1 [47], LINC00310, LINC00341, etc.) (Supplemental Table S3).

Fig. 5
figure 5

Signaling pathway analysis of differentially expressed genes regulated by BRD7. a A total of 1648 genes were significantly differentially expressed after BRD7 overexpression. We classified these dysregulated genes based on enriched Gene Ontology (GO) terms and mapped these genes onto the KEGG canonical pathways using the DAVID software. The top five categories and pathways are presented. b Dysregulated genes in the canonical apoptosis signaling pathway are shown. Red represents induction, and green represents repression. (Color figure online)

Full size image

Integrative analysis of the ChIP-seq and DGE data

To identify direct and functionally relevant BRD7 target genes, we analyzed the ChIP-seq dataset against transcriptional profiling data. A total of 184 BRD7 target genes identified by ChIP-seq, as well as 1648 BRD7-regulated genes identified by DGE profiling, were entered into the interaction database STRING to find their functional interacting relationships. Next, a regulatory network of BRD7 target genes and regulated genes was constructed and illustrated by Cytoscape (Supplemental Figure S1); in this network, BRD7 was surrounded by its target genes, while other BRD7-regulated genes were laid in the outer cycle. Because some BRD7 target genes were not significantly differentially expressed (less than 2.0-fold followed by BRD7 over-expression), we filtered out these genes and their downstream interacting genes from Supplemental Figure S1 and focused on those genes that are differentially expressed (2.0-fold or greater) and involved in the “cell death and survival” network (Fig. 4) and their downstream interacting genes in the regulating network illustrated in Fig. 6. A number of differentially expressed genes were involved in the regulation of the cell cycle and apoptosis, such as BID, CASP6, CASP3, UBE2L3, and AIFM1, which interact with potential BRD7 target genes such as BIRC2, BIRC3, TXN2, and NOTCH1.

Fig. 6
figure 6

Regulatory network of BRD7 target genes and their potential downstream genes identified by ChIP-seq and DGE profiling. A network of BRD7 target genes and their downstream genes, containing BRD7 target genes that were differentially expressed (2.0-fold or great) and involved in the “cell death and survival” network and their downstream interacting genes, was constructed. Red represents induction, and green represents repression. (Color figure online)

Full size image

Validating expression of BRD7 downstream target genes

To confirm the expression patterns measured by DGE, four differentially expressed genes, namely BIRC2, BIRC3, TXN2, and NOTCH1, were chosen for validation using qRT-PCR in BRD7-overexpressing or control HEK293 and HeLa cells. The results showed that BIRC2, BIRC3, and NOTCH1 were significantly decreased by BRD7, while TXN2 was induced by BRD7 (Fig. 7), and this finding is consistent with our DGE data.

Fig. 7
figure 7

Expression of the BRD7 target genes BIRC2, BIRC3, NOTCH1, and TXN2 was validated by real-time PCR in HEK293 and HeLa cells

Full size image

Discussion

Accumulating evidence demonstrates that inactivation of the BRD7 gene contributes to the development of many human cancers by transcriptional dysregulation [1, 2, 6–13]. BRD7 regulates signaling pathways that involve cell growth, apoptosis, cell cycle, and mobility [14–21]. As a component of SWI/SNF chromatin-remodeling complexes, BRD7 interacts with other transcriptional factors, such as BRD2 [15], BRCA1 [22], and p53 [16], and transcriptionally regulates their target genes and downstream signaling pathways. However, given the divergent cellular roles of the BRD7 gene, the BRD7 target gene and its downstream gene regulating network have not yet been explored. Therefore, uncovering the function of BRD7, especially its downstream genes and regulating network, will shed light on the development of cancers and may provide a novel strategy for cancer therapy from the view of transcription dysregulation.

ChIP analysis is one of the most common approaches to studying the binding sites of transcription factors (TFs) and the mechanisms of TF functions. Owing to the tremendous progress in next-generation sequencing technology, ChIP-seq offers higher resolution, less noise, and greater coverage than its array-based predecessor ChIP-chip. With the decreasing cost of sequencing, ChIP-seq has become an indispensable tool for studying gene regulation and epigenetic mechanisms [48, 49]. DGE profiling with massive parallel sequencing achieved high sensitivity and reproducibility for transcriptome profiling. Although it lacks the ability of detecting alternative splicing events compared to RNA-SEQ, it is much more affordable and clearly out-performed microarrays in detecting lower abundant transcripts [50]. While not suffering from some of the disadvantages of hybridization-based detection (ChIP-chip or microarray), sufficient sequencing depth and an appropriate peak calling or tag cognizing algorithm are essential for ensuring robustness of conclusions derived from ChIP-seq or DGE data. Sequencing-based experiments generate large quantities of data, and effective computational analysis is also crucial for uncovering biological mechanisms.

In this study, we used ChIP-Seq to identify BRD7-enriched regions in human cells stably transfected with the BRD7 gene and found 156 BRD7-binding peaks located within or adjacent to 184 genes. These peaks were enriched with 4 novel motifs, and, to our knowledge, these are the first reported BRD7-specific binding motifs. BRD7 is a member of the bromodomain-containing protein family. Bromodomain has specific binding affinity for acetylated lysines on the N-terminal tails of histones [51]. Bromodomain proteins may facilitate the accession of transcription factors to chromatin by modulating chromatin remodeling or the acetylation of histones [52]. We mapped 156 BRD7-binding peaks to histone modification sites of the human genome, based on data derived from the ENCODE project, and found that most BRD7-binding peaks were indeed co-localized with histone modification sites. These results further confirm the BRD7 function on histone modification. Although there were a relatively small number of BRD7 direct binding signatures, it is consistent with the function of BRD7 on cell cycle and apoptosis previously reported [2, 14, 17]. Based on IPA analysis, many of these BRD7 target genes were involved in the cell death and survival gene regulating network.

We also undertook a comprehensive transcriptome analysis using DGE profiling by RNA-sequencing to identify differentially expressed genes regulated by BRD7. A total of 1649 dysregulated genes were identified in BRD7-overexpressing cells, including 560 up-regulated genes and 1088 down-regulated genes. We integrated the ChIP data with the gene expression data and found that 41 BRD7-binding genes are dysregulated by BRD7 by at least twofold. To further examine the molecular mechanisms underlying the BRD7 gene, a gene regulatory network was constructed for the dysregulated genes and BRD7 target genes based on their regulatory relationships. Some hub nodes in the network, such as BIRC2, BIRC3, NOTCH1, and TXN, and their downstream signaling pathways are related to cell growth and apoptosis. Through ChIP-PCR and qRT-PCR, we further confirmed our ChIP-Seq and DEG results of the BIRC2, BIRC3, NOTCH1, and TXN2 genes, in both HEK293 and HeLa cells, and identified these genes as direct target genes of BRD7.

BIRC2 and BIRC3, members of the inhibitor of apoptosis proteins (IAPs) group, inhibit apoptosis by binding to tumor necrosis factor receptor-associated factors and reducing the activity of their downstream caspases [53, 54], or activating the ubiquitination of some apoptosis-related factors [55]. They are also involved in the regulation of G2/M cell cycle process [56, 57]. TXN2 encodes a redox protein located in mitochondria [58]. TXN2 protein activates caspase-3 by binding to the regulatory region of pro-caspase-3, which advances cellular apoptosis in human lymphocytes and liver cancer cells [59]. In HeLa cells, TXN2 functions mainly by preventing the generation of reactive oxygen species (ROS) in mitochondria to reduce apoptosis [60]. NOTCH1 encodes a transmembrane receptor. Recent studies have demonstrated that overexpression of NOTCH1 can advance the occurrence and development of breast cancer [61], lung cancer [62], glioma [63], and acute lymphoblastic leukemia [64]. Our results suggest that these BRD7 target genes, namely BIRC2, BIRC3, NOTCH1, and TXN2, may play important roles in BRD7-induced apoptosis and cell cycle signaling and that thoroughly investigating their downstream genes and regulating networks will contribute to a more complete understanding of the function of BRD7, especially in carcinogenesis.

Interestingly, there were 104 non-coding RNA sequences dysregulated by BRD7, most of them being long non-coding RNA sequences (lncRNA). LncRNAs compose a large set of RNA molecules that exceed 200 nt in length, completely lack or possess limited protein-coding capacity, and represent a substantial portion of the transcriptome [65–72]. LncRNAs widely regulate gene expression at the epigenetic, transcriptional, and post-transcriptional levels. Substantial evidence indicates that the aberrant expression or dysfunctional activities of lncRNAs are correlated with tumor initiation and progression. For example, we first reported that LINC00472 and GNAS-AS1 were up-regulated by BRD7; these results are consistent with their tumor suppressor function. Previous studies have reported that they were down-regulated in breast cancer [41] and colorectal cancer [42]. Additionally, MALAT1 [43–46] and CDKN2B-AS1 (ANRIL) [47] were down-regulated by BRD7 and act as oncogenes in multiple types of cancers.

In summary, our unbiased genome-scale screening for BRD7 target genes combining ChIP-seq and DGE analyses provides an important resource for elucidating the biological function of BRD7. This study is the first to report the genome-wide chromatin occupancy of BRD7 and the construction of a BRD7-regulated gene network. We identified 4 novel BRD7-binding motifs and identified BIRC2, BIRC3, NOTCH1, and TXN2 as BRD7 target genes. We also identified lncRNAs that are regulated by BRD7, and these results imply that lncRNAs are a novel component of the BRD7 regulatory network and are worth further investigation.