乳腺癌中p53调控增强子的特征与功能分析

张 茵,王文著,谭政堂,李昶蓥,岳俊杰,郭志云*

(1.西南交通大学生命科学与工程学院,中国四川 成都 610031;2.中国人民解放军军事科学院军事医学研究院生物工程研究所,中国北京 100071)

增强子是一类具有组织特异性的顺式调控元件,通过富集多种转录因子调控基因的表达,并且不受距离和方向的限制。先前研究表明,大多数活性增强子在转录因子介导下可转录出RNA,即增强子RNA(enhancer RNA,eRNA)[1]。已有研究证实eRNA的表达量和增强子的活性相关,增强子的表达失调会导致包括乳腺癌在内的多种癌症的发生[2]。

p53(tumor protein 53)作为一种肿瘤抑制转录因子参与调控绝大多数肿瘤的发生发展[3]。研究表明,大部分p53在染色体上的结合位点除了具有保守的p53模体序列外,还具有增强子的信号特征,并且这些p53结合的区域具有增强子的活性[4]。随后相关研究进一步证实,在阿霉素诱导DNA损伤的情况下,p53的表达量显著升高,而且大部分p53都结合在增强子区域并起着调节增强子活性的能力[5]。

然而,p53如何介导增强子调控乳腺癌的发生与发展以及p53调控增强子的特征目前尚不清楚。为此,本文分析了乳腺癌细胞MCF-7中p53的染色质免疫沉淀测序(chromatin immunoprecipitation sequencing,ChIP-seq)数据,发现结合p53的增强子(Enhp53)与未结合p53的增强子(Enhno-p53)在表达量、组蛋白修饰和转录因子结合方面存在显著差异。同时,我们分析了MCF-7细胞中p53活性升高前后的RNA-seq数据,在近端调控和远端调控两个层面上共识别了148对Enhp53-mRNAs。这些mRNA的功能富集分析结果表明,受Enhp53调控的mRNA与肿瘤显著相关,并显著影响了乳腺癌病人的总生存时间。综上所述,p53通过介导增强子调控参与乳腺癌通路,从而影响肿瘤的发生与发展,这一结果为进一步探讨转录因子与增强子的调控提供了理论依据与方法基础。

1 材料与方法

1.1 材料

MCF-7细胞中的活性增强子来自于HACER数据库(http://bioinfo.vanderbilt.edu/AE/HACER/)[6],MCF-7细胞中p53的ChIP-seq数据来自于SRA(Sequence Read Archive)数据库(SRR287800)[7]。9种组蛋白修饰(包括 H3K4me1、H3K4me2、H3K4me3、H3K9ac、H3K9me2、H3K9me3、H3K27ac、H3K27me3、H4K20me1)和76种转录因子的ChIP-seq数据来自于ENCODE(The Encyclopedia of DNA Elements)项目。MCF-7细胞在Nutlin处理前后的差异表达RNA来自于 Léveillé等[8]的研究。MCF-7细胞在空间上的染色质相互作用数据来自于4DGenome数据库(https://4dgenome.research.chop.edu/)[9]。

1.2 MCF-7细胞中p53的ChIP-seq数据分析

使用FastQC[10]对测序数据进行质量检测,并用Cutadapt[11]切除引物和质量分数低于10的片段。使用Bowtie2[12]将Nutlin处理过的p53 ChIP-seq数据匹配到hg19的人类基因组上,并用SAM-tools[13]过滤未匹配的测序片段,随后使用MACS2[14]的callpeak进行peak calling,最后通过deepTools[15]得到p53的bigWig文件。

1.3 Enhp53的识别以及组蛋白修饰和转录因子结合分析

将MCF-7细胞中来自HACER数据库的增强子与p53在染色体上的结合位点进行位点匹配,得到结合p53的增强子(Enhp53),即该增强子受p53调控。通过ENCODE项目获得MCF-7细胞中组蛋白修饰以及转录因子的ChIP-seq数据,进而使用bwtool[16]分析Enhp53和Enhno-p53中点上下游1 kb内的组蛋白修饰和转录因子结合情况,并将Enhp53和Enhno-p53的转录因子结合信号的均值进行比较,比值大于1.2的转录因子被认定为具有协同p53调控增强子的功能。

1.4 Enhp53调控的差异表达mRNA的识别

本文通过两种方法识别Enhp53调控的差异表达mRNA。首先,通过4DGenome数据库得到MCF-7细胞在空间上的染色质相互作用数据,若增强子和mRNA分别位于相互作用的两个片段内,则认为该增强子调控mRNA。其次,若mRNA位于一个增强子的上下游100 kb内,则认为mRNA受该增强子调控[17]。

1.5 受Enhp53调控的mRNA的功能富集分析以及生存分析

利用R包clusterProfiler[18]对mRNA进行GO(Gene Ontology)和KEGG(Kyoto Encyclopedia of Genes and Genomes)富集分析。从TCGA(The Cancer Genome Atlas)数据库获取乳腺癌的RNA-seq数据以及临床数据,绘制Enhp53调控的mRNA的Kaplan-Meier生存曲线[19],P＜0.05即认为该mRNA对乳腺癌患者的总生存时间具有显著影响。

2 结果

2.1 MCF-7细胞中Enhp53的识别

为了探讨MCF-7细胞中受p53调控的增强子的情况,我们首先从HACER数据库中得到MCF-7细胞中8 714个增强子,随后对MCF-7细胞中的p53 ChIP-seq数据进行分析,得到7 784个p53的结合位点,通过对增强子和p53的结合位点进行位点匹配,最终共识别出459个Enhp53。通过分析发现,我们识别的增强子普遍具有活性增强子的信号特征,例如:增强子chr7:579 426～579 712具有高水平的H3K27ac和H3K4me1信号以及低水平的H3K4me3信号,这与文献[20]的报道相符,并且在增强子区域富集p53和脱氧核糖核酸酶(deoxyribonuclease,DNase)(图 1)。

图1 增强子chr7:579 426～579 712的基因组信息示例图黄色标记增强子的位置,绿色标记增强子上的各种信号特征,蓝色代表组蛋白修饰在基因组上的信号分布,红色代表结合在DNA上的转录因子的信号分布。Fig.1 Diagram of genomic information of enhancer chr7:579 426～579 712Yellow marks the position of the enhancer,green marks the various signal features on the enhancer,blue represents the distribution of histone modification signals,and red represents the signal distribution of transcription factors bound to DNA.

2.2 Enhp53具有更强的转录活性

先前研究表明,p53可结合在增强子上调控肿瘤的发生与发展[3]。为了研究p53的结合对增强子功能的影响,我们以是否结合p53为标准将增强子分为两类:Enhp53和Enhno-p53(图2A)。首先,我们对两类增强子的表达量进行比较,结果表明Enhp53的表达量显著高于Enhno-p53(图2B)。研究报道,活性增强子往往具有多种组蛋白修饰[20],为了研究p53在增强子上的结合是否会引起增强子组蛋白修饰信号的变化,我们对Enhp53和Enhno-p53进行了组蛋白特征分析。结果表明,在9种组蛋白修饰中,5种在增强子上出现了明显的双峰,分别 是 H3K4me1、H3K4me2、H3K4me3、H3K9ac 和H3K27ac,这5种组蛋白修饰信号在Enhp53上显著高于 Enhno-p53(图 2C～G)。

2.3 Enhp53结合转录因子调控增强子活性

为了进一步探明增强子结合转录因子的情况以及p53协同哪些转录因子共同调控增强子,我们比较了两类增强子上转录因子的结合情况。结果表明,这76种转录因子均在两类增强子上有结合信号。通过进一步分析转录因子在两类增强子上的信号比值,共识别了36个协同p53调控增强子的转录因子(图3A),其中GATA3(GATA binding protein 3)、FOXA1(forkhead box A1)以及DPF2(double PHD fingers 2)是增强子上受p53影响最显著的3个转录因子(图3B～D)。此外,参与增强子和靶基因空间成环的关键转录因子CTCF(CCCTC-binding factor)也在Enhp53上具有更强的信号(图 3E)。

图2 两类增强子的p53结合信号、表达量以及组蛋白修饰信号的差异(A)两类增强子结合p53的差异;(B)两类增强子的表达量差异,P值由t检验计算得出;(C～G)两类增强子的多种组蛋白修饰信号差异。Fig.2 Differences between the two types of enhancers in p53-binding signals,expression levels and histone modification signals(A)The difference between the two types of enhancers in p53-binding signals;(B)The difference in expression of the two types of enhancers.The P value is calculated by the t test;(C～G)The differences in multiple histone modification signals on two types of enhancers.

2.4 Enhp53调控的差异表达mRNA的识别以及功能分析

为了探讨Enhp53如何参与mRNA的差异表达,我们分析了MCF-7细胞在Nutlin处理前后的RNA-seq数据,得到1 817个差异表达的mRNA。其中,998个mRNA显著上调,819个mRNA显著下调。为了尽可能全面地识别Enhp53调控的mRNA,我们通过远端与近端的增强子-mRNA调控模式识别Enhp53调控的mRNA,结合4DGenome数据库中染色质相互作用数据和增强子上下游100 kb两种方法共识别出148组Enhp53-mRNAs,涉及差异表达的mRNA共120个。在这些mRNA中,多个mRNA的表达量与乳腺癌患者的总生存时间显著相关,例如 FOS、FOSL1(FOS like 1)、ARC(activity-regulated cytoskeleton-associated protein)、BTG2(BTG anti-proliferation factor 2)等。其中,在Nutlin处理后升高2.14倍的FOS受到两个增强子调控,分别是chr14:75 649 674～75 650 038以及chr14:75 721 756～75 722 023。通过TCGA中乳腺癌患者的临床数据以及RNA-seq数据绘制FOS的Kaplan-Meier生存曲线,结果表明,FOS的表达对患者生存时间有显著影响,FOS高表达的患者明显具有更长的总生存时间(图4A)。此外,我们对Enhp53调控的mRNA进行了GO和KEGG富集分析,结果表明这些基因在DNA损伤、细胞凋亡以及p53介导的信号转导过程中发挥着重要作用,并显著富集在多种癌症发生和p53相关的通路中(图 4B,C)。

图3 两类增强子的转录因子结合信号差异(A)两类增强子的转录因子结合信号的比值;(B～E)两类增强子的多种转录因子结合信号差异。Fig.3 Differences in transcription factor binding signals between two types of enhancers(A)The ratio of transcription factor binding signals of the two types of enhancers;(B～E)Differences in multiple transcription factor binding signals of two types of enhancers.

3 讨论

图4 FOS的生存曲线以及受Enhp53调控的mRNA的GO和KEGG分析(A)FOS的Kaplan-Meier生存曲线;(B)受Enhp53调控的mRNA的GO分析;(C)受Enhp53调控的mRNA的KEGG分析。Fig.4 The survival curve of FOS and GO and KEGG analyses of mRNAs regulated by Enhp53(A)Kaplan-Meier survival curve for FOS;(B)GO analysis of mRNAs regulated by Enhp53;(C)KEGG analysis of mRNAs regulated by Enhp53。

本研究发现,Enhp53和Enhno-p53在表达量、组蛋白修饰以及转录因子结合等方面均存在显著差异。Enhp53的表达量显著高于Enhno-p53,暗示p53的结合进一步促进增强子的活性,从而促进下游基因的表达。H3K4me1、H3K4me2、H3K4me3、H3K9ac和H3K27ac均是已被证实的位于活性增强子上的组蛋白修饰[21],本文研究显示，这5种组蛋白修饰信号在Enhp53上显著高于Enhno-p53。上述结果表明,p53结合增强子后可能通过改变这些组蛋白修饰信号从而起到调节增强子活性的作用。在对比分析中,Enhp53的H3K4me3信号值略高于H3K4me1,造成这一结果的原因可能是ENCODE项目中不同组蛋白修饰的ChIP-seq数据来自不同的实验室和样本,存在一定的实验误差,因此直接对H3K4me3和H3K4me1的信号值进行比较并不准确。此外,HACER数据库中的增强子是通过CAGE(cap analysis of gene expression)实验得来的,CAGE识别增强子依靠双向转录的RNA,并不依赖组蛋白修饰信号,因此可能包含着未知基因的启动子。已有研究表明,在多个组织和细胞中普遍存在的增强子往往具有更高的H3K4me3信号[22],因此可以针对MCF-7细胞中增强子的特异性进行进一步的研究。尽管如此,p53的结合对增强子H3K27ac和H3K4me1的影响显著大于H3K4me3(图 2)。

之前的研究表明,增强子通过结合转录因子与靶基因的启动子相互作用,从而调控基因的表达[23]。文中对两类增强子上多种转录因子的结合信号进行了比较,发现GATA3、FOXA1以及DPF2是增强子上受p53影响最显著的3个转录因子(图3)。已有研究表明,一部分p53预先结合在染色体不可及区域,在DNA损伤导致p53升高的情况下,这些p53结合位点的染色体可及性从不可及变成可及[5],而GATA和FOXA家族被证实可以促进染色质的开放[24],因此GATA3、FOXA1等转录因子很有可能参与到p53结合位点的可及性变化过程中。DPF2在细胞凋亡过程中发挥着重要作用[25],而p53在不同的刺激下会调节细胞周期和细胞凋亡[26];通过查询BioGRID数据库(https://thebi ogrid.org/),我们发现DPF2与p53结合,因此DPF2和p53可能在细胞凋亡过程中发挥协同作用。此外,CTCF在Enhp53上具有更强的信号,这一结果暗示p53可能在介导增强子与靶基因空间成环方面起到一定的作用。综上所述,p53通过与多种转录因子直接或间接相互作用调控增强子活性与染色体的可及性,并且可能参与CTCF介导的染色体成环。

p53与增强子、靶基因可形成复杂的调控网络。研究发现,p53除了结合增强子外,也可能同时结合在靶基因的启动子区域[5]。为了研究结合在启动子区域的p53对Enhp53的影响,我们将Enhp53调控的靶基因的启动子与p53 ChIP-seq数据进行了位点匹配,结果显示25%的受Enhp53调控的靶基因启动子区有p53结合。我们推测这部分基因可能与Enhp53、p53形成前馈环路,也可能在三维结构上发生增强子与靶基因启动子成环从而拉近两者距离使p53同时调控增强子与靶基因启动子,当然这些结论仍需要进一步的实验验证。

另外,本文通过近端和远端两种调控模式识别了148对Enhp53-mRNAs,其中差异表达的FOS被两个增强子调控,且FOS与乳腺癌病人的总生存时间显著相关(图4A)。已有研究表明,在肿瘤细胞中过表达FOS可增强细胞凋亡[27],这暗示FOS可与p53联合作用,通过调节相关基因的表达来发挥抗肿瘤作用。以上研究表明,p53通过介导增强子调控参与乳腺癌通路,从而影响肿瘤的发生发展,该结果为进一步探讨转录因子与增强子调控提供了理论依据与方法基础。