前列腺癌差異基因篩選、功能預(yù)測及臨床意義分析

楊 智 米 軍

蘭州大學(xué)第二醫(yī)院泌尿外科研究所/甘肅省泌尿系疾病重點實驗室/甘肅省泌尿系疾病臨床醫(yī)學(xué)中心(甘肅 蘭州 73000)

前列腺癌的發(fā)病率在全球范圍內(nèi)不斷升高。 在美國,前列腺癌位于惡性腫瘤及癌癥相關(guān)住院率的第二位[1]。 有研究表明,多種mRNA 參與前列腺癌發(fā)生、發(fā)展過程。 本研究基于美國基因表達(dá)綜合數(shù)據(jù)庫NCBI Gene Expression Omnibus(GEO)和美國癌癥腫瘤基因圖譜(the cancer genome atlas, TCGA)數(shù)據(jù)庫對前列腺癌組織及良性組織中的差異表達(dá)基因進(jìn)行篩選,并探討其相關(guān)分子機制和臨床預(yù)后意義。

材料與方法

一、數(shù)據(jù)提取

從美國基因表達(dá)綜合數(shù)據(jù)庫NCBI Gene Expression Omnibus(GEO)下載目標(biāo)前列腺癌組織及非癌組織mRNA 轉(zhuǎn)錄組數(shù)據(jù),形成三個數(shù)據(jù)集。 其中數(shù)據(jù)集GSE46602 包含36 例前列腺癌樣本,14 例前列腺良性組織樣本;數(shù)據(jù)集GSE3325 包含13 例前列腺癌樣本,6例前列腺良性增生組織樣本;數(shù)據(jù)集GSE104749 包含4例前列腺癌組織樣本,4 例前列腺良性增生組織樣本。

二、差異表達(dá)分析

在R 語言環(huán)境中,利用edgeR package 對GSE46602,GSE3325 和GSE104749 三個選擇的前列腺癌數(shù)據(jù)集進(jìn)行正常組織與癌癥組織的差異表達(dá)分析,篩選出差異表達(dá)基因。 edge R 工具包讀取文件,包括基因ID、log2轉(zhuǎn)化后的差異倍數(shù)(Fold Change)值、顯著性p值以及校正后p值(默認(rèn)偽發(fā)現(xiàn)率,FDR 校正)。 然后,設(shè)定篩選閾值為:FC 值＜1 且P＜0.05,用于比較篩選出前列腺癌組織和非癌組織的差異表達(dá)的基因。 定義|log2FC|＞1 以及FDR＜0.01 為差異表達(dá)基因。 其中l(wèi)og2FC＞1 設(shè)定為“up”,即差異表達(dá)基因中的上調(diào)基因;log2FC＜-1設(shè)定為“down”,即差異表達(dá)基因中的下調(diào)基因。 獲得識別了顯著差異表達(dá)的基因后,載入ggplot2,繪制數(shù)據(jù)集的火山圖。 載入pheatmap package,將獲得的基因差異表達(dá)結(jié)果,完成基因ID 轉(zhuǎn)換后,獲得相應(yīng)的基因表達(dá)熱圖。

三、關(guān)鍵基因篩選

利用韋恩圖,選定上調(diào)的差異表達(dá)基因中前六位的基因認(rèn)定為關(guān)鍵基因,進(jìn)一步進(jìn)行關(guān)鍵基因的功能分析和臨床預(yù)后預(yù)測。

四、GO 及KEGG 信號通路分析

在R 語言環(huán)境下,載入openxlsx R 包,用于讀取之前獲得的差異表達(dá)篩選的文件。 載入stringr R 包用于完成基因的ID 轉(zhuǎn)換;載入enrichplot R 包和clusterProfiler R 包用于對功能分析的基因進(jìn)行富集分析及結(jié)果可視化,最后利用barplot R 包和dotplot R 包完成圖片繪制。

五、關(guān)鍵基因的蛋白組學(xué)表達(dá)情況

在HPA(Human Protein Atlas,https:/ /www.proteinatlas.org/)直接檢索目標(biāo)基因可獲得目標(biāo)基因的組織特性、組織及單細(xì)胞表達(dá)聚類性質(zhì)及蛋白質(zhì)功能等一系列信息。

六、基因表達(dá)水平及生存分析

比較關(guān)鍵基因在前列腺癌組織及非癌組織中的表達(dá)水平。 利用GEPIA(http:/ /gepia. cancer-pku. cn/),結(jié)合TCGA 和GTEx 相關(guān)數(shù)據(jù),可直接分析出單一基因在腫瘤組織和非腫瘤組織的表達(dá)差異和相關(guān)性結(jié)果。利用UCSC Xena(https:/ /xena. ucsc. edu),將目標(biāo)基因選擇“Kaplan Meier plot”,輸出目標(biāo)基因表達(dá)量和無進(jìn)展間隔期(progression-free interval,PFI)的生存曲線結(jié)果。

七、統(tǒng)計學(xué)方法

本研究采用SPSS 22、R 語言(R 4.04)及R studio(2021.09. 1 Build 372)處理數(shù)據(jù)。 Cytoscape(version 3.5.1 java 1.8.2_291)、GraphPad prism 8 繪制圖像。 計量資料若符合正態(tài)分布,以均數(shù)±標(biāo)準(zhǔn)差(±s) 表示,比較用t檢驗,若非正態(tài)性分布則采用非參數(shù)檢驗;計數(shù)資料以率(%)表示,比較用卡方檢驗。 UCSC Xena采用Kaplan-Meier 法繪制目標(biāo)基因表達(dá)量的生存曲線,分為兩組或三組,比較采用Log rankχ2檢驗。P＜0.05為差異有統(tǒng)計學(xué)意義。

結(jié) 果

一、差異表達(dá)基因

根據(jù)差異基因的篩選條件,GSE46602 共讀取229個差異表達(dá)基因,包含68 個上調(diào)的差異表達(dá)基因和161 個下調(diào)的差異表達(dá)基因;GSE3325 共讀取468 個差異表達(dá)基因,包含232 個上調(diào)的差異表達(dá)基因和236 個下調(diào)的差異表達(dá)基因;GSE104749 共讀取1007 個差異表達(dá)基因,包含402 個上調(diào)的差異表達(dá)基因和605 個下調(diào)的差異表達(dá)基因。 篩選出上調(diào)基因中,靠前六位的基因作為差異表達(dá)基因的關(guān)鍵基因,關(guān)鍵基因為:PCA3(前列腺癌抗原3), SIM2(Single-Mind 2), AMACR(異構(gòu)梅2-甲酰輔酶A 消旋梅), HOXC6(同源異型盒-6 基因), ERG(ETS 轉(zhuǎn)錄調(diào)節(jié)因子)及TOP2A(拓?fù)洚悩?gòu)梅ⅢA 型)。 (見圖1A-圖1H)

圖1 篩選差異表達(dá)基因的熱圖、火山圖及前六位關(guān)鍵基因的韋恩圖

二、GO 及KEGG 富集分析結(jié)果

GO 結(jié)果顯示,差異表達(dá)基因的生物功能主要與減數(shù)分裂核分裂、減數(shù)分裂細(xì)胞周期、三價鐵結(jié)合和氧化還原酶活性相關(guān)。 KEGG 數(shù)據(jù)庫結(jié)果顯示,差異表達(dá)基因主要富集的代謝途徑包括:原代膽汁酸生物合成、嘧啶代謝、谷胱甘肽代謝、鉑類耐藥性及p53 信號通路(見圖2)。

圖2 GO+KEGG 可視化結(jié)果

三、蛋白互作網(wǎng)絡(luò)結(jié)果

將數(shù)據(jù)導(dǎo)入STRING 數(shù)據(jù)庫獲得蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)(見圖3),PPI 共有17 點節(jié)點,35 個面;從類聚結(jié)果分析,主要類聚為三個部分,分別以TOP2A、TACC3、AMACR、TP53、SIM2、ERG、HOXC6 等為中心。

圖3 蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò)

四、關(guān)鍵基因表達(dá)水平驗證

借助人類蛋白質(zhì)表達(dá)圖集(The Human Protein Atlas,HPA),能查看和分析關(guān)鍵基因編碼蛋白在不同腫瘤組織/器官中的表達(dá)情況(見圖4,由于PCA3為非編碼基因,故未納入PCA3 在組織中的染色標(biāo)本)。

圖4 HPA 數(shù)據(jù)庫中關(guān)鍵基因的組織學(xué)染色

五、關(guān)鍵基因生存分析

結(jié)合TCGA 數(shù)據(jù)庫中前列腺癌的臨床相關(guān)信息,將篩選的關(guān)鍵基因做臨床預(yù)后分析,結(jié)果見圖5-7。 對于PCA3,正常非癌組織中的表達(dá)量和T2、T3 和T4 期前列腺癌組織比較,其表達(dá)量均有統(tǒng)計學(xué)差異。 PCA3 在前列腺癌N 分期中,隨著疾病的進(jìn)展,PCA3 的表達(dá)量階梯式上升,但是PCA3 和前列腺癌的遠(yuǎn)處轉(zhuǎn)移無明顯統(tǒng)計學(xué)差異。 SIM2 和AMACR 在前列腺癌中表達(dá)水平和T 分期基本無關(guān)系,非癌組織中的表達(dá)量明顯低于前列腺癌組織;N 分期和M 分期的結(jié)果類似,SIM2 和AMACR 的表達(dá)量隨著疾病的進(jìn)展并未明顯變化。HOXC6 和ERG 表達(dá)量在TNM 分期的各組間基本無明顯的統(tǒng)計學(xué)差異,多為良性前列腺組織的表達(dá)量低于前列腺癌組織。 TOP2A 的表達(dá)量與PCA3 的趨勢相似,隨著癌癥的進(jìn)展,TOP2A 的表達(dá)量在TNM 分期中呈現(xiàn)階梯式上升。 進(jìn)一步分析關(guān)鍵基因與前列腺癌無進(jìn)展間隔期(progression-free interval, PFI)的關(guān)系。 PCA3 和TOP2A的表達(dá)情況對前列腺癌的無進(jìn)展間隔期(progression-free interval,PFI)時間存在統(tǒng)計學(xué)差異(p＜0.05),PCA3 的高表達(dá)量組和低表達(dá)量組的PFI 存在統(tǒng)計學(xué)差異(p=0.005), HR 為0.55(0.36-0.84);TOP2A 的高表達(dá)量組和低表達(dá)量組的PFI 也存在統(tǒng)計學(xué)差異(p＜0.001),HR 為2. 33(1. 51-3. 62),其余關(guān)鍵基因(SIM2、AMACR、HOXC6 和ERG)未顯著影響前列腺癌PFI。

圖5 關(guān)鍵基因表達(dá)量對TNM 分期的影響

圖6 關(guān)鍵基因表達(dá)量對TNM 分期的影響

圖7 關(guān)鍵基因表達(dá)量對PFI 的生存分析

討 論

前列腺癌作為美國男性最常見的惡性腫瘤,全球約有21% 的癌癥病例為前列腺癌[1]。 本研究結(jié)合GEO、TCGA 等相關(guān)數(shù)據(jù)庫,利用生物信息學(xué)技術(shù),篩選出與前列腺癌相關(guān)的六個關(guān)鍵基因,分別為PCA3、SIM2、AMACR、HOXC6、ERG 和TOP2A。 篩選的關(guān)鍵基因中,PCA3 和TOP2A 在前列腺癌中表達(dá)較顯著,與前列腺癌TNM 分期,腫瘤的無進(jìn)展間隔期顯著現(xiàn)關(guān),PCA3 和TOP2A 有望成為預(yù)測前列腺癌預(yù)后的相關(guān)指標(biāo)及新的腫瘤標(biāo)志物。

PCA3 是一種長非編碼RNA(lncRNA),最初通過差異顯示實驗法( differential display experimental approach)被命名為DD3[2]。 PCA3 通過調(diào)節(jié)雄激素受體(AR)信號傳導(dǎo)進(jìn)而影響前列腺癌(PCa)細(xì)胞功能,同時控制多個雄激素應(yīng)答和癌癥相關(guān)基因的表達(dá),包括上皮-間質(zhì)轉(zhuǎn)化(EMT)過程。 此外,尿液中的PCA3已成功應(yīng)用于PCa 的臨床診斷[3]。 Qin Z 等[4]對8139名病例和14116 名對照病例進(jìn)行薈萃分析,以評估PCA3 的診斷價值,PCA3 的總診斷優(yōu)勢比(DOR)較高,95%置信區(qū)間為5.44(4.53-6.53)。 Ye LF等[5]發(fā)明了一種新方法,即LBXexo 評分法,用于測量尿液中的外泌體PCA3/PRAC 表達(dá)水平,結(jié)果證實LBXexo 可以改進(jìn)中國人群前列腺癌的診斷。 對于臨床應(yīng)用,Alshalafa M 等[6]發(fā)現(xiàn)PCA3 在初始活檢中預(yù)測高級別前列腺癌(Gleason Score,GS≥8)的表現(xiàn)較差,假陰性率較高。Kotova ES 等[7]在PCa 患者尿液樣本中檢測PCA3 和AMACR 的mRNA 水平。 在其研究中,PCA3 評分AUC為0.632(95%置信區(qū)間:0.511-0.752),AMACR 評分AUC 為0.711(95%置信區(qū)間:0.617-0806),較高的預(yù)測價值使其成為潛在的前列腺癌生物標(biāo)記物。

SIM2(single-Mind 2)基因是具有基本螺旋-環(huán)-螺旋/果蠅Per 蛋白-芳香烴受體核轉(zhuǎn)運蛋白-果蠅Sim 蛋白(helix-loop-helix/per-Arnt-Sim,bHLH/PAS)結(jié)構(gòu)域的轉(zhuǎn)錄因子家族的成員,影響實體瘤的發(fā)生過程。 在乳腺癌中,SIM2s 直接下調(diào)SNAI2(Snail Family Transcriptional Repressor 2)的表達(dá),抑制EMT,從而抑制腫瘤生長和侵襲[8]。 Wyatt GL 等[9]發(fā)現(xiàn)SIM2s 和NFκB(核因子κB,nuclear factor kappa-B)之間存在串?dāng)_,對乳腺腫瘤的生長和轉(zhuǎn)移具有抑制作用。

AMACR 最為人所知的生理作用是催化支鏈脂肪酸的α-甲基質(zhì)子在線粒體和過氧化物酶體中的β-氧化[10]。 前列腺癌細(xì)針穿刺活檢中,AMACR 表達(dá)對PCa檢測具有97%的敏感性和100%的特異性[11]。 Xie H等[12]使用miR200c 抑制前列腺癌中的AMACR 活性,研究證實抑制AMACR 可以阻礙細(xì)胞增殖和遷移。

Homeobox C6,HOXC6 屬于同源盒家族(homeobox,HOX),其成員編碼高度保守的轉(zhuǎn)錄因子,其家族有39個HOX 基因位于4 個染色體位點。 HOXC6 已被證實為侵襲性前列腺癌中可靠的生物標(biāo)記物,可提高早期診斷率并預(yù)測治療后癌癥復(fù)發(fā)[14]-[17]。 Luo Z 等[13]發(fā)現(xiàn),大多數(shù)HOXC4 結(jié)合位點與HOXC6 結(jié)合位點重疊。通過全基因組分析,發(fā)現(xiàn)HOXC6 與侵襲性前列腺癌的臨床預(yù)后指標(biāo)相關(guān)。

ERG 是屬于ETS 轉(zhuǎn)錄因子家族,是紅細(xì)胞轉(zhuǎn)化特異性(ETS)家族成員,與胚胎發(fā)育、細(xì)胞增殖、分化、血管生成、炎癥和凋亡相關(guān)[18]。 ERG 家族基因與致癌基因融合有關(guān),是前列腺癌在內(nèi)的多種癌癥的典型特征[19]-[21]。 Eerola SK 等[22]通過分析mRNA 表達(dá)和染色質(zhì)免疫沉淀測序(ChIP-Seq)數(shù)據(jù)集,發(fā)現(xiàn)PIM(Proviras Integration of Maloney Kinase)與MYC 原癌基因和ERG癌蛋白在前列腺癌的發(fā)生和進(jìn)展中具有協(xié)同作用。

TOP2A 作為拓?fù)洚悩?gòu)酶Ⅱ(TOP2)家族的兩個成員之一,僅在循環(huán)細(xì)胞中表達(dá),并負(fù)責(zé)DNA 鏈復(fù)制過程中的酶解偶聯(lián)[23]。 據(jù)報道,在結(jié)腸癌組織中TOP2A 過度表達(dá)[24]。 然而,TOP2A 在組織發(fā)育和腫瘤進(jìn)展中的作用和機制仍然未知[25]。 在根治性前列腺切除術(shù)組織中,TOP2A 與細(xì)胞周期進(jìn)程(cell cycle progression,mCCP)的相關(guān)性最顯著(r=0.7)[26]。 Cattrini C 等[27]還發(fā)現(xiàn)TOP2A mRNA 過度表達(dá)與患者預(yù)后不良相關(guān)。

綜上,本研究通過對前列腺癌差異表達(dá)基因的篩選,建立蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò),探究關(guān)鍵基因參與的生物學(xué)過程和代謝途徑,結(jié)合前列腺癌的臨床數(shù)據(jù),評價PCA3 和TOP2A 對前列腺癌臨床預(yù)后的關(guān)系。 由于基因測序結(jié)果、篩選閾值、實驗方法、數(shù)據(jù)庫選擇等多種因素的影響,得出的結(jié)論存在偏差,需要進(jìn)一步細(xì)胞、動物以及大樣本的人體實驗的驗證。