聯(lián)合應(yīng)用RNA-seq和ATAC-seq尋找FOXQ1轉(zhuǎn)錄因子的下游靶基因

伊夢(mèng)杰，張錦錦，郭 強(qiáng)，唐 慧*

(1．昆明理工大學(xué) 醫(yī)學(xué)院， 昆明 650504；2.云南省第一人民醫(yī)院，昆明理工大學(xué)附屬醫(yī)院，云南省消化內(nèi)鏡臨床醫(yī)學(xué)中心,云南省臨床病毒學(xué)重點(diǎn)實(shí)驗(yàn)室，昆明市腫瘤分子與免疫防治重點(diǎn)實(shí)驗(yàn)室，昆明 650032；3.中國(guó)科學(xué)院昆明動(dòng)物研究所，昆明 650223)

結(jié)直腸癌(Colorectal cancer, CRC)是一種由多致病因素導(dǎo)致且預(yù)后很差的惡性消化道腫瘤，是世界第三大常見(jiàn)腫瘤，其死亡率在所有癌癥中居于第四位(位于肺癌、肝癌、胃癌之后)[1]。叉頭框Q1(Forkhead box Q1,FOXQ1)是叉頭框(Forkhead box,FOX)基因家族的成員之一，基因定位于6p23-25，編碼含403個(gè)氨基酸的FOXQl蛋白，作為一類核轉(zhuǎn)錄因子，可以穩(wěn)定結(jié)合到靶基因啟動(dòng)子區(qū)域的GC盒等核心元件，控制下游基因轉(zhuǎn)錄活性從而發(fā)揮生物效應(yīng)[2]。2010 年Kaneda等[3]研究發(fā)現(xiàn)FOXQ1在CRC中異常高表達(dá)，近年來(lái)大量研究亦證實(shí)FOXQ1在卵巢癌、乳腺癌、膀胱移形細(xì)胞癌、胃癌、肝癌、非小細(xì)胞肺癌和神經(jīng)腦膠質(zhì)瘤等多種腫瘤中異常表達(dá)[4-6]，與多種腫瘤的發(fā)生、發(fā)展密切相關(guān)，在許多腫瘤中具有明確的促腫瘤生長(zhǎng)轉(zhuǎn)移的功能。因此，鑒定新的FOXQ1下游靶基因，豐富由FOXQ1參與并介導(dǎo)的信號(hào)通路信息，有望為腫瘤的靶向治療提供新的靶點(diǎn)。

ATAC-seq是一種檢測(cè)染色質(zhì)易接近性的測(cè)序方法，可以檢測(cè)樣本間染色質(zhì)易接近性的變化情況。ATAC-seq只需要很少的細(xì)胞(50 000個(gè))就能檢測(cè)基因組中所有活躍的調(diào)控序列[7]，DNA探針(作為轉(zhuǎn)座子發(fā)揮作用)通過(guò)酶促反應(yīng)(轉(zhuǎn)座酶Tn5)被整合到基因組的開(kāi)放區(qū)域，然后通過(guò)測(cè)序來(lái)鑒定這些區(qū)域[8]。本實(shí)驗(yàn)聯(lián)合應(yīng)用RNA-seq與ATAC-seq尋找由FOXQ1敲低引起的染色質(zhì)易接近性改變所導(dǎo)致的表達(dá)改變基因，即FOXQ1轉(zhuǎn)錄因子調(diào)控的潛在下游基因。

1 材料和方法

1.1 細(xì)胞、儀器和試劑

本實(shí)驗(yàn)所用細(xì)胞系均購(gòu)于中科院上海細(xì)胞庫(kù)，DMEM高糖培養(yǎng)液、RPMI 1640培養(yǎng)液購(gòu)于美國(guó)Corning公司，胎牛血清購(gòu)于美國(guó)Gibco公司。Puromycin購(gòu)于北京索萊寶科技有限公司，F(xiàn)OXQ1一抗購(gòu)于abcam公司，HRP-Rb-anti-goat二抗購(gòu)于Cell Signaling Technology公司，β-actin一抗、HRP-goat-anti-mouse二抗均購(gòu)于Proteintech公司。RNAzol RT RNA Isolation Reagent購(gòu)于美國(guó)MRC公司，引物由寶生物工程有限公司合成。核酸定量?jī)x為Thermo公司產(chǎn)品，LightCycler 480實(shí)時(shí)熒光定量?jī)x為Roche公司產(chǎn)品，WB垂直基礎(chǔ)電泳儀為伯樂(lè)公司產(chǎn)品。

1.2 方法

1.2.1 實(shí)驗(yàn)分組

根據(jù)課題組前期工作并結(jié)合文獻(xiàn)報(bào)道，選取在CRC細(xì)胞系中FOXQ1表達(dá)量較高的SW480和DLD1細(xì)胞系進(jìn)行FOXQ1基因的敲低實(shí)驗(yàn)，實(shí)驗(yàn)分為兩組共4株細(xì)胞。FOXQ1基因敲低組：SW480-shFOXQ1和DLD1-shFOXQ1；對(duì)照組：SW480-shControl和DLD1-shControl。

1.2.2 FOXQ1基因低表達(dá)穩(wěn)定細(xì)胞系的構(gòu)建

1.2.2.1 慢病毒表達(dá)質(zhì)粒的構(gòu)建、篩選與擴(kuò)增

根據(jù)FOXQ1序列，得到siRNA靶序列，設(shè)計(jì)合成3對(duì)shRNA干擾序列(見(jiàn)表1)，命名為shRNA-A、shRNA-B、shRNA-C，并設(shè)計(jì)一條對(duì)照序列shRNA-control；應(yīng)用Addgene的pSPAX2慢病毒包裝系統(tǒng)，PLKO.1-puro-shFOXQ1質(zhì)粒包裝構(gòu)建獲得lenti-shFOXQ1慢病毒，每組重組質(zhì)粒經(jīng)測(cè)序鑒定無(wú)誤后，用去內(nèi)毒素大提試劑盒提取質(zhì)粒，測(cè)濃度；酶切鑒定無(wú)誤后，在HEK293T細(xì)胞中擴(kuò)增，之后用LB培養(yǎng)基篩選、擴(kuò)增[9]。

表1 shRNA引物序列Table 1 shRNA primer sequences

1.2.2.2 慢病毒的感染及陽(yáng)性細(xì)胞的篩選

將待感染的SW480和DLD1細(xì)胞鋪在24孔板中培養(yǎng)，生長(zhǎng)密度達(dá)70%～80%時(shí)將病毒濃縮液加至細(xì)胞中，培養(yǎng)12～15 h后更換為完全培養(yǎng)基。待細(xì)胞長(zhǎng)滿時(shí)，將24孔板中的細(xì)胞傳代至六孔板中進(jìn)行初步擴(kuò)大培養(yǎng)，并加入含有1 μg/ml嘌呤霉素的培養(yǎng)基篩選陽(yáng)性細(xì)胞。篩選獲得的陽(yáng)性細(xì)胞分別應(yīng)用qRT-PCR和WB在mRNA和蛋白水平驗(yàn)證FOXQ1基因的敲低效率[10]。

1.2.3 RNA-seq

1.2.3.1 收集細(xì)胞

將穩(wěn)定轉(zhuǎn)染的細(xì)胞傳代培養(yǎng)穩(wěn)定生長(zhǎng)后，分別收集生長(zhǎng)狀態(tài)良好的DLD1、SW480基因敲低組與對(duì)照組細(xì)胞進(jìn)行后續(xù)RNA-seq實(shí)驗(yàn)。

1.2.3.2 RNA的提取和檢測(cè)

提取基因敲低組與對(duì)照組細(xì)胞總RNA，用Nanodrop 2000檢測(cè)所提總RNA的濃度和純度后用瓊脂糖凝膠電泳檢測(cè)RNA的完整性，最后利用Agilent 2100測(cè)定RIN值。

1.2.3.3 RNA的富集和測(cè)序

使用結(jié)合有poly-T寡核苷酸的磁珠從總RNA中分離出含有poly-A的mRNA，加入片段化緩沖劑將其打成片段，將片段化的mRNA逆轉(zhuǎn)錄成cDNA并純化；經(jīng)末端修復(fù)、poly-A添加、測(cè)序接頭連接及AMPure XP beads 篩選后得出大小合適的片段，進(jìn)行 PCR 擴(kuò)增，建立測(cè)序文庫(kù)并進(jìn)行文庫(kù)質(zhì)控[11]。

1.2.3.4 上機(jī)測(cè)試。

對(duì)構(gòu)建合格的測(cè)序文庫(kù)進(jìn)行雙末端(Paired-end)測(cè)序，本測(cè)序?qū)嶒?yàn)組和對(duì)照組均設(shè)置了三次生物學(xué)重復(fù)，測(cè)序工作由上海嘉因生物科技有限公司完成。

1.2.3.5 RNA-seq的數(shù)據(jù)分析

RNA-seq數(shù)據(jù)下機(jī)后首先進(jìn)行空載去除、接頭去除等數(shù)據(jù)預(yù)處理后產(chǎn)出原始數(shù)據(jù)，緊接著利用FastQC軟件進(jìn)行數(shù)據(jù)質(zhì)量控制，利用生物信息學(xué)軟件STAR、HTSeq和DESeq2對(duì)測(cè)序結(jié)果進(jìn)行參考序列比對(duì)、表達(dá)量統(tǒng)計(jì)、差異基因篩選等分析[12-14]，并通過(guò)R軟件繪制基因聚類分析圖、火山圖等。分別建立DLD1和SW480細(xì)胞FOXQ1敲低前后的轉(zhuǎn)錄譜。

1.2.4 ATAC-seq

1.2.4.1 收集細(xì)胞

分別收集DLD1、SW480基因敲低組和對(duì)照組細(xì)胞，計(jì)數(shù)50 000個(gè)細(xì)胞，離心去上清后依次用預(yù)冷的PBS、lysis buffer懸浮細(xì)胞500g，4 °C離心去除上清液，立即進(jìn)行轉(zhuǎn)座反應(yīng)[15]。

1.2.4.2 轉(zhuǎn)座反應(yīng)與純化

確保細(xì)胞始終置于冰上，配置轉(zhuǎn)座反應(yīng)體系懸浮細(xì)胞，37 °C孵育30 min，立即用Qiagen MinElute PCR Purification Kit純化DNA，之后用10μl elution buffer洗脫。

1.2.4.3 PCR擴(kuò)增

配置PCR反應(yīng)體系循環(huán)擴(kuò)增，之后用Qiagen MinElute PCR Purification Kit純化DNA[8]。

1.2.4.4 上機(jī)測(cè)試

庫(kù)檢合格后，把不同文庫(kù)按照有效濃度及目標(biāo)下機(jī)數(shù)據(jù)量的需求合并后進(jìn)行Illumina HiSeq測(cè)序。測(cè)序基于邊合成邊測(cè)序(Sequencing by Synthesis)的原理進(jìn)行，在序的流動(dòng)池中加入四種熒光標(biāo)記的dNTP、DNA聚合酶以及接頭引物進(jìn)行擴(kuò)增，在每一個(gè)測(cè)序簇延伸互補(bǔ)鏈時(shí)，每加入一個(gè)被熒光標(biāo)記的dNTP就能釋放出相對(duì)應(yīng)的熒光，測(cè)序儀通過(guò)捕獲熒光信號(hào)，并通過(guò)計(jì)算機(jī)軟件將光信號(hào)轉(zhuǎn)化為測(cè)序峰，從而獲得待測(cè)片段的序列信息。本測(cè)序?qū)嶒?yàn)由上海嘉因生物科技有限公司完成。

1.2.4.5 ATAC-seq的數(shù)據(jù)分析

首先對(duì)ATAC-seq下機(jī)數(shù)據(jù)進(jìn)行預(yù)處理，如FastQC質(zhì)量控制，去除duplicate sequences，去除blacklist reads[16]，原始序列比對(duì)，比對(duì)后去除重復(fù)序列和細(xì)胞器序列等，再利用生物信息學(xué)軟件BWA、Macs2、Homer、Deep tools等對(duì)ATAC-seq測(cè)序結(jié)果進(jìn)行參考序列比對(duì)分析、差異結(jié)合位點(diǎn)檢測(cè)、特征峰在全基因組上的分布注釋等，建立DLD1和SW480細(xì)胞FOXQ1敲低前后的染色質(zhì)易接近性變化圖譜。

1.2.5 RNA-seq和ATAC-seq數(shù)據(jù)的關(guān)聯(lián)分析

分別將DLD1和SW480兩個(gè)細(xì)胞系中RNA-seq的差異基因和ATAC-seq的差異表達(dá)峰進(jìn)行關(guān)聯(lián)分析，得到兩組關(guān)聯(lián)分析結(jié)果，以明確兩組細(xì)胞中染色質(zhì)易接近性變化區(qū)域?qū)ο掠位虻恼{(diào)控功能，并找出每種細(xì)胞實(shí)驗(yàn)組與對(duì)照組相比差異染色質(zhì)易接近性區(qū)域可能調(diào)控的下游基因。

1.2.6 兩個(gè)細(xì)胞系關(guān)聯(lián)分析數(shù)據(jù)的重疊分析

從DLD1和SW480的關(guān)聯(lián)分析結(jié)果中找出兩細(xì)胞共有的重疊部分，即獲得FOXQ1敲低后引起的染色質(zhì)易接近性改變所導(dǎo)致的表達(dá)改變基因。

1.2.7 代謝通路分析

將重疊分析結(jié)果中所獲基因基于代謝通路數(shù)據(jù)庫(kù)(Kyoto Encyclopedia of Genes and Genomes, KEGG)進(jìn)行代謝通路注釋，得到注釋基因參與的所有代謝通路名稱，采用Fisher檢驗(yàn)計(jì)算代謝通路的顯著性水平(P<0.05)，從而篩選出注釋基因富集的顯著性代謝通路[17]。

2 結(jié)果分析

2.1 在mRNA和蛋白水平驗(yàn)證 DLD1和SW480細(xì)胞中FOXQ1的敲低效率

qRT-PCR和WB實(shí)驗(yàn)結(jié)果(見(jiàn)圖1)表明成功構(gòu)建了FOXQ1基因敲低組DLD1-sh-FOXQ1、SW480-sh-FOXQ1 以及對(duì)照組DLD1-sh-Control和SW480-sh-Control。FOXQ1基因敲低組在mRNA和蛋白水平與對(duì)照組相比FOXQ1基因的表達(dá)量均顯著降低(***P<0.001)。

圖1 穩(wěn)定低表達(dá)FOXQ1的DLD1和SW480細(xì)胞驗(yàn)證(***P<0.001)Fig.1 Verification of the FOXQ1 expression quantity in DLD1 and SW480 cells which knock down FOXQ1 gene steadily (***P<0.001)

2.2 RNA-seq

RNA-seq結(jié)果表明，DLD1細(xì)胞FOXQ1基因敲低組與對(duì)照組相比表達(dá)顯著上調(diào)的基因有215個(gè)，表達(dá)顯著下調(diào)的基因有131個(gè)；SW480細(xì)胞FOXQ1敲低組與對(duì)照組相比表達(dá)顯著上調(diào)的基因有171個(gè)，表達(dá)顯著下調(diào)的基因有358個(gè)(表達(dá)差異基因的篩選閾值為FDR<0.05，log2FC>1)。對(duì)這些表達(dá)顯著差異基因進(jìn)行分析，發(fā)現(xiàn)在FOXQ1基因敲低后，與侵襲、自噬相關(guān)的EI24基因在表達(dá)下調(diào)基因中位居前列，與先天免疫應(yīng)答和炎癥反應(yīng)相關(guān)的TLR2基因以及與遷移相關(guān)的SMAD3基因在表達(dá)上調(diào)的基因中位居前列。推測(cè)這些基因的表達(dá)改變與CRC的發(fā)生、侵襲、發(fā)展等密切相關(guān)，將進(jìn)一步結(jié)合ATAC-seq分析這幾種基因在FOXQ1敲低后，染色質(zhì)開(kāi)放區(qū)域的變化情況(見(jiàn)圖2)。

2.3 ATAC-seq

ATAC-seq結(jié)果表明DLD1-shFOXQ1組有39 146個(gè)特征峰，DLD1-shControl組有49 381個(gè)特征峰，其中，基因敲低組與對(duì)照組相比染色質(zhì)易接近性增強(qiáng)的差異表達(dá)峰有2 385個(gè)，染色質(zhì)易接近性減弱的差異表達(dá)峰有6 205個(gè)。SW480-shFOXQ1組有38 962個(gè)特征峰，SW480-shControl組有42 244個(gè)特征峰，其中，基因敲低組與對(duì)照組相比染色質(zhì)易接近性增強(qiáng)的差異表達(dá)峰有4 563個(gè)，染色質(zhì)易接近性減弱的差異表達(dá)峰有3 733個(gè)(見(jiàn)圖3)。

圖3 ATAC-seq分析圖Fig.3 ATAC-seq analysis

2.4 RNA-seq與ATAC-seq的關(guān)聯(lián)分析

本實(shí)驗(yàn)將RNA-seq結(jié)果中的差異基因和ATAC-seq結(jié)果中差異染色質(zhì)開(kāi)放區(qū)域進(jìn)行關(guān)聯(lián)分析，預(yù)測(cè)由FOXQ1基因敲低所導(dǎo)致的染色質(zhì)開(kāi)放區(qū)域改變引起的轉(zhuǎn)錄因子的結(jié)合能力的改變，及最終導(dǎo)致下游基因的表達(dá)上調(diào)或下調(diào)。結(jié)果顯示，DLD1染色質(zhì)易接近性減弱區(qū)域?qū)ο掠位虮磉_(dá)具有抑制作用(見(jiàn)圖4(a))；染色質(zhì)易接近性增強(qiáng)區(qū)域?qū)ο掠位虮磉_(dá)具有促進(jìn)作用(見(jiàn)圖4(b))。在SW480細(xì)胞中也觀察到相同的結(jié)果(見(jiàn)圖4(c)、圖4(d))[4]。

圖4 染色質(zhì)功能預(yù)測(cè)Fig.4 Prediction of chromatin function

2.5 SW480與DLD1關(guān)聯(lián)分析結(jié)果的重疊分析

將兩個(gè)細(xì)胞系的關(guān)聯(lián)分析結(jié)果基因進(jìn)行重疊分析。重疊部分為在兩種細(xì)胞中染色質(zhì)易接近性變化相同且調(diào)控的表達(dá)量變化趨勢(shì)相同的基因。其中， DLD1細(xì)胞系中染色質(zhì)易接近性減弱區(qū)域調(diào)控的下調(diào)基因有1 436個(gè)，SW480細(xì)胞系中染色質(zhì)易接近性減弱區(qū)域調(diào)控的下調(diào)基因有407個(gè)(見(jiàn)圖5)，兩細(xì)胞系染色質(zhì)易接近性減弱區(qū)域共同調(diào)控的表達(dá)下調(diào)基因有70個(gè)。DLD1細(xì)胞系中染色質(zhì)易接近性增強(qiáng)區(qū)域調(diào)控的上調(diào)基因有792個(gè)，SW480細(xì)胞系中染色質(zhì)易接近性增強(qiáng)區(qū)域調(diào)控的上調(diào)基因有531個(gè)(見(jiàn)圖5)，兩細(xì)胞系染色質(zhì)易接近性增強(qiáng)區(qū)域共同調(diào)控的表達(dá)上調(diào)基因有61個(gè)。

圖5 SW480與DLD1關(guān)聯(lián)分析結(jié)果的重疊分析Fig.5 Overlapping analysis of correlation analysis results in SW480 and DLD1

在這兩個(gè)細(xì)胞的交集基因中發(fā)現(xiàn)EI24、TLR2、SMAD3基因也在其中，所以觀察這些基因在染色質(zhì)易接近性方面發(fā)生的變化(見(jiàn)圖6～圖8)。其中，TLR2、SMAD3基因在ATAC-seq數(shù)據(jù)中基因敲低組與對(duì)照組相比染色質(zhì)區(qū)域有較為明顯的變化，而EI24基因的染色質(zhì)區(qū)域變化很弱。

圖6 ATAC-seq EI24基因比對(duì)Fig.6 ATAC-seq EI24 genome comparison

圖7 ATAC-seq TLR2基因比對(duì)Fig.7 ATAC-seq TLR2 genome comparison

圖8 ATAC-seq SMAD3基因比對(duì)Fig.8 ATAC-seq SMAD3 genome comparison

2.6 代謝通路分析

將兩個(gè)細(xì)胞系的重疊基因進(jìn)行代謝通路分析，結(jié)果表明染色質(zhì)易接近性增強(qiáng)區(qū)域調(diào)控的表達(dá)上調(diào)基因顯著性富集的代謝通路包括炎癥性腸病(Inflammatory bowel disease, IBD)、甘油磷脂代謝通路(Glycerophospholipid metabolism)、細(xì)胞周期代謝通路(Cell cycle)等。其中SMAD3、TLR2基因顯著性富集在IBD代謝通路(見(jiàn)圖9(a))。染色質(zhì)易接近性減弱區(qū)域調(diào)控的表達(dá)下調(diào)基因顯著性富集的代謝通路包括p53信號(hào)通路(p53 signaling pathway)、TRP通道炎癥調(diào)控通路(Inflammatory mediator regulation of TRP channels)、半胱氨酸和蛋氨酸代謝通路(Cysteine and methionine metabolism)等。EI24基因顯著性富集在p53信號(hào)通路(圖9(b))。

圖9 代謝通路分析Fig.9 Metabolic pathway analysis

有研究表明，患有IBD的患者與健康人群相比患結(jié)直腸癌的風(fēng)險(xiǎn)更高[18-19]，當(dāng)腸道系統(tǒng)與其微生物群之間的關(guān)系(包括屏障功能、免疫信號(hào)和代謝物)受到干擾后引起的慢性炎癥是其發(fā)病的主要潛在原因[20]。其次，在CRC中P53信號(hào)通路也是一種重要的信號(hào)通路，p53抑癌基因突變是導(dǎo)致CRC發(fā)生的最主要原因之一，同時(shí)也是結(jié)直腸癌侵襲和轉(zhuǎn)移的原因之一，還有研究認(rèn)為p53突變?cè)谙倭?癌轉(zhuǎn)移過(guò)程中也發(fā)揮重要作用[21]。進(jìn)一步研究IBD、P53信號(hào)通路在CRC發(fā)生發(fā)展中的作用將具有重要意義。

3 討 論

一般來(lái)說(shuō)，染色質(zhì)有“關(guān)閉”“開(kāi)放”兩種狀態(tài)，處于“關(guān)閉”狀態(tài)的染色質(zhì)，在異染色質(zhì)蛋白以及修飾酶的作用下，被包裝成致密、緊湊的結(jié)構(gòu)，阻遏轉(zhuǎn)錄因子等蛋白的結(jié)合，此時(shí)染色質(zhì)處于沉默失去生物功能的階段；而處于“開(kāi)放”狀態(tài)的染色質(zhì)，具有不太緊致的結(jié)構(gòu)，可招募轉(zhuǎn)錄因子等蛋白的結(jié)合，進(jìn)而調(diào)控基因的表達(dá)水平。ATAC-seq作為一種繪制全基因組染色質(zhì)可及性圖譜的方法[8]，利用超活性Tn5轉(zhuǎn)座酶檢測(cè)染色質(zhì)的可接近性，是本實(shí)驗(yàn)的重要研究手段之一[7]。通過(guò)對(duì)ATAC-seq和RNA-seq數(shù)據(jù)的分析，初步確定了FOXQ1基因敲低后發(fā)生差異表達(dá)的基因，豐富了FOXQ1轉(zhuǎn)錄因子的下游調(diào)控網(wǎng)絡(luò)。