陳 思，劉春良，趙 倩，孫海鵬，劉云霞

上海交通大學基礎(chǔ)醫(yī)學院病理生理學系，細胞分化與凋亡教育部重點實驗室，上海200025

乳腺癌是最為常見的惡性腫瘤之一，也是女性中發(fā)病率最高的惡性腫瘤[1]。根據(jù)雌激素受體（estrogen receptor，ER）、孕激素受體（progesterone receptor，PR）和人類表皮生長因子受體2（human epidermal growth factor receptor-2，HER2）的表達情況，可以將乳腺癌分為4種亞型（Luminal A、Luminal B、HER2陽性和三陰乳腺癌）。根據(jù)不同乳腺癌亞型的生物學特征和臨床病理分期制定相應的個體化治療策略，可以使乳腺癌5年生存率達90%以上[2]。然而，乳腺癌的復發(fā)和轉(zhuǎn)移仍舊是一大難題。除此之外，一些乳腺癌亞型，例如三陰乳腺癌，由于缺乏有效的治療靶點，一直以來是臨床治療的一個瓶頸。因此，明確乳腺癌發(fā)生發(fā)展相關(guān)的關(guān)鍵基因和通路，有助于認識乳腺癌潛在的發(fā)病機制，或?qū)榕R床尋找更多的診斷和治療靶點提供參考。

自基因芯片技術(shù)和高通量測序技術(shù)問世以來，生物信息學迅速發(fā)展，目前已發(fā)現(xiàn)了許多疾病的生物學標志物[3-4]。從信息豐富的公共數(shù)據(jù)庫如GEO數(shù)據(jù)庫（Gene Expression Omnibus，基因表達匯編）和TCGA數(shù)據(jù)庫（The Cancer Genome Atlas，癌癥基因組圖譜）中可獲得基因表達數(shù)據(jù)，通過生物信息學方法，對數(shù)據(jù)庫中的基因表達數(shù)據(jù)進行聚類分析、統(tǒng)計分析、通路分析和可視化作圖等，能夠預測基因的功能以及基因間的相互作用，了解疾病基因?qū)用娴陌l(fā)病機制，發(fā)現(xiàn)潛在的生物學標志物，從而為疾病的分子靶向藥物研發(fā)和精準治療提供理論依據(jù)。本研究將乳腺癌基因表達數(shù)據(jù)與臨床生存分析相結(jié)合，以篩選樞紐基因和關(guān)鍵信號通路；基于生存期篩選出來的基因可能更具有臨床意義，或能為乳腺癌的診斷和治療提供新的思路。

1 材料與方法

1.1 數(shù)據(jù)采集

為了獲取乳腺癌的基因表達數(shù)據(jù)集，本研究在GEO數(shù)據(jù)庫下載了3個乳腺癌數(shù)據(jù)集（http://www.ncbi.nlm.nih.gov/geo/），分別為GSE54002、GSE29431和GSE61304，這3個數(shù)據(jù)集都是基于GPL570平臺。GSE54002包含417例乳腺癌樣本和16例正常樣本。GSE29431包含了54例乳腺癌樣本和12例正常樣本：54例乳腺癌樣本中有15例HER2免疫組織化學評分為3+，且伴有HER2基因擴增；26例評分為2+，其中13例伴有HER2基因擴增，13例不伴有HER2基因擴增；13例評分為0/1+，且不伴有HER2基因擴增。GSE61304包含了58例乳腺癌樣本和4例正常樣本：58例乳腺癌樣本中18例為ER+PR+，19例為ER-PR-，4例為ER+PR-，1例為ER-PR+，其他16例樣本未說明。

1.2 篩選差異表達基因

在Rstudio軟件中（版本3.4.0）加載并下載Bioconductor網(wǎng)站上軟件包來分析上述3個乳腺癌數(shù)據(jù)集。首先使用affy包導入CEL文件，使用simpleaffy包評估微陣列數(shù)據(jù)質(zhì)量[5]，gcrma包中的RMA算法預處理原始數(shù)據(jù)[6]，genefilter包過濾非特異性結(jié)合的探針和數(shù)據(jù)質(zhì)量低的探針，limma包進行差異基因表達的統(tǒng)計學驗證[7]?；虮磉_變化倍數(shù)的對數(shù)值的絕對值（|log2FC|） ＞1且P＜0.05認定為差異表達基因。最后，為了提高差異表達基因的穩(wěn)健性，使用Funrich軟件（版本3.1.3）獲得3個數(shù)據(jù)集中均上調(diào)或下調(diào)的基因，用于下一步的分析。

1.3 篩選與總體生存期相關(guān)的差異表達基因

為了分析差異表達基因?qū)θ橄侔┗颊呖傮w生存期的影響，在Kaplan-Meier plotter數(shù)據(jù)庫（www.kmplot.com）中根據(jù)基因的中位表達值將患者樣本分為高表達組和低表達組。使用默認參數(shù)，計算每個基因高表達組和低表達組的中位生存期；若log-rank P＜0.05，則該基因被視為與總體生存期相關(guān)的差異表達基因。

1.4 分析與總體生存期相關(guān)的差異表達基因的功能

基因本體數(shù)據(jù)庫（Gene Ontology，GO）富集分析和京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genome，KEGG）通路富集分析被廣泛用于識別基因的功能和通路。本研究使用Rstudio軟件中clustProfiler包對與乳腺癌患者總體生存期相關(guān)的差異表達基因進行GO分析和KEGG分析，P＜0.05認為具有統(tǒng)計學意義[8]。

1.5 蛋白 - 蛋白相互作用網(wǎng)絡(luò)構(gòu)建及篩選樞紐基因

蛋白與蛋白相互作用在調(diào)節(jié)生物學過程中起著至關(guān)重要的作用。這種關(guān)系可以通過蛋白 - 蛋白相互作用（proteinprotein interaction，PPI）網(wǎng)絡(luò)表示，每個節(jié)點代表一個蛋白，邊代表蛋白質(zhì)之間的相互作用。緊密相連的區(qū)域可以作為富集功能群。STRING數(shù)據(jù)庫（https://string-db.org/）包含了豐富的蛋白質(zhì)之間相互作用的信息[9]。為了評估與總體生存期相關(guān)的差異表達基因之間的相互關(guān)系，將差異表達基因列表導入STRING數(shù)據(jù)庫，并設(shè)定信度為0.4；接著將PPI表格數(shù)據(jù)導入Cytoscape軟件中構(gòu)建PPI網(wǎng)絡(luò)，使用軟件中的插入式分子復合物檢測（MCODE評分＞4分，節(jié)點數(shù)＞5個）篩選出PPI網(wǎng)絡(luò)中的樞紐模塊，最后通過CytoHubba插件計算網(wǎng)絡(luò)中每一個基因的最大團中心性（maximal clique centrality，MCC）分數(shù)，將得分前10的基因作為樞紐基因[10-12]。

1.6 樞紐基因驗證

1.6.1 數(shù)據(jù)庫驗證樞紐基因表達 利用Oncomine數(shù)據(jù)庫（https://www.oncomine.org/resource/main.html/）和人類蛋白質(zhì)圖譜（Human Protein Atlas，HPA）數(shù)據(jù)庫（https://www.proteinatlas.org/）對樞紐基因在乳腺癌腫瘤組織和正常組織間的mRNA和蛋白水平的表達進行驗證[13]。

1.6.2 RNA抽提及實時熒光定量PCR 人乳腺癌細胞MDA-MB-231和人正常乳腺上皮細胞MCF-10A（購自中國科學院干細胞庫）的總RNA抽提按照Invitrogen公司TRIzol試劑盒提供的方法。反轉(zhuǎn)錄以1 000 ng RNA為模板，按TaKaRa公司M-MLV Reverse Transcription Kit試劑盒說明配制反轉(zhuǎn)錄反應液，進行反轉(zhuǎn)錄反應。實時熒光定量PCR（quantitative real-time PCR，qPCR）使用ABI公司7500 Real-Time PCR System試劑盒。每份采用10 μL體系，cDNA模板1 μL，每份樣品做3個復孔，求平均值；以18S rRNA為內(nèi)參，引物序列見表1。

表1 qPCR引物序列Tab 1 Primer sequences for qPCR

1.7 統(tǒng)計學分析

應用OriginPro 2017C和Adobe Illustraor CS6軟件進行統(tǒng)計學分析和作圖，qPCR數(shù)值用x—±s表示，組間比較采用獨立樣本t檢驗。P＜0.05認為差異具有統(tǒng)計學意義。

2 結(jié)果

2.1 篩選差異表達基因

在|log2FC|＞1且P＜0.05的篩選條件下，從GSE54002中得到差異表達基因3 389個，其中上調(diào)基因1 561個，下調(diào)基因1 828個；GSE29431中得到差異表達基因3 660個，其中上調(diào)基因1 097個，下調(diào)基因2 563個；GSE61304中得到差異表達基因1 828個，其中上調(diào)基因821個，下調(diào)基因1 007個。然后用Funrich軟件的Vene圖篩選得到了211個共同上調(diào)和374個共同下調(diào)的差異表達基因，如圖1所示。

圖1 GSE29431、GSE54002和GSE61304中篩選得到的差異表達基因的Venn圖Fig 1 Venn diagrams of the differentially expressed genes screened in GSE29431, GSE54002 and GSE61304

2.2 乳腺癌中與總體生存期相關(guān)的差異表達基因

為了篩選與總體生存期相關(guān)的基因，使用Kaplan-Meier plotter數(shù)據(jù)庫對這585個差異表達基因進行總體生存分析，計算每個基因的高表達組和低表達組的中位生存期和log-rankP值。結(jié)果顯示，有262個基因的高表達組和低表達組的總體生存期之間的差異有統(tǒng)計學意義。表2列出了與生存相關(guān)最顯著的前10個差異表達基因。這262個與總體生存期相關(guān)的差異表達基因?qū)⒂糜谙乱徊降幕蚬δ芊治觥?/p>

表2 與總體生存期相關(guān)最顯著的前10個差異表達基因Tab 2 List of the top 10 overall survival-related differentially expressed genes

2.3 乳腺癌中與總體生存期相關(guān)的差異表達基因的功能分析

為進一步了解這些基因的功能，利用Rstudio軟件中的clusterProfiler包對得到的生存期相關(guān)的差異表達基因進行功能分析，P＜0.05認為具有統(tǒng)計學意義。GO功能分析結(jié)果顯示，這些基因的功能主要與細胞分裂（cell division）、核分裂（nuclear division）、細胞器分裂（organelle fission）、細胞周期的調(diào)控（regulation of mitotic cell cycle）、負向調(diào)控細胞有絲分裂（mitotic nuclear division）以及染色體分離（chromosome segregation）等生物學過程有關(guān)（圖2A）。KEGG通路分析顯示這些基因參與細胞周期（cell cycle）、人T細胞白血病病毒1感染（human T-cell leukemia virus 1 infection）、FoxO信號通路（FoxO signaling pathway）、卵母細胞減數(shù)分裂（oocyte meiosis）、細胞衰老（cellular senescence）以及孕酮介導的卵母細胞成熟（progesterone-mediated oocyte maturation）等信號通路（圖2B）。

圖2 與總體生存期相關(guān)的差異表達基因GO功能分析和KEGG通路分析Fig 2 GO function analysis and KEGG pathway analysis of the overall survival-related differentially expressed genes

2.4 PPI網(wǎng)絡(luò)構(gòu)建及模塊分析

將與總體生存期相關(guān)的差異表達基因列表上傳至STRING，并設(shè)定信度0.4作為判斷相互作用是否有意義的標準，構(gòu)建了PPI網(wǎng)絡(luò)（圖3）。Cytoscape軟件根據(jù)MCODE評分排序篩選出了3個樞紐模塊（圖4），并對樞紐模塊中的基因進行GO功能分析，結(jié)果顯示模塊1和模塊2的基因主要富集在細胞分裂（cell division）、細胞核分裂（nuclear division）和細胞器分裂（organelle fission）等生物學過程，模塊3的基因主要集中在血小板脫顆粒（platelet degranulation）、負向調(diào)控線粒體細胞色素c的釋 放（negative regulation of release of cytochrome c from mitochondria）等生物學過程（表3）。

圖3 與總體生存期相關(guān)的差異表達基因的PPI網(wǎng)絡(luò)Fig 3 PPI network of overall survival-related differentially expressed genes

圖4 MCODE評分排序篩選出的3個樞紐模塊Fig 4 Three modules selected by MCODE scoring sorting

表3 3個模塊基因的GO功能分析Tab 3 GO functional analysis of the differentially expressed genes in three modules

2.5 乳腺癌中樞紐基因的鑒定

樞紐基因是一類在生物學過程中發(fā)揮至關(guān)重要作用的基因，在相關(guān)通路中其他非樞紐基因的調(diào)控往往要受到這類基因的影響，因此樞紐基因有可能成為乳腺癌的生物學標志物和治療靶標。使用Cytoscape中的插件Cytohubba，通過MCC法得到了得分前10的樞紐基因，分別是NDC80、BUB1、CDCA8、BUB1B、BIRC5、CCNB1、KIF2C、CENPF、MAD2L1和CDC20（圖5）。

圖5 通過MCC法得到的10個樞紐基因及其相互作用Fig 5 Ten hub genes and their interactions by MCC

2.6 樞紐基因的驗證

通過Oncomine數(shù)據(jù)庫和HPA數(shù)據(jù)庫對篩選得到的10個樞紐基因進行表達驗證。Oncomine數(shù)據(jù)庫結(jié)果顯示，BIRC5、CDC20、NDC80、CENPF、MAD2L1、CDCA8、KIF2C、BUB1、CCNB1和BUB1B的mRNA水平在乳腺癌組織中明顯上調(diào)（圖6）。HPA數(shù)據(jù)庫免疫組織化學檢測結(jié)果顯示，腫瘤 組 織 中CDCA8、BIRC5、CDC20、CENPF、MAD2L1和CCNB1的蛋白表達水平高于正常乳腺組織，另外4個基因NDC80、KIF2C、BUB1、BUB1B未被HPA數(shù)據(jù)庫收錄。以上結(jié)果提示，篩選得到的樞紐基因具有較好的穩(wěn)健性。

圖6 通過Oncomine數(shù)據(jù)庫驗證樞紐基因轉(zhuǎn)錄水平的表達Fig 6 Validation of the expression of hub genes at mRNA level using Oncomine database

此外，通過提取人乳腺癌細胞MDA-MB-231和人正常乳腺上皮細胞MCF-10A的RNA進行反轉(zhuǎn)錄和qPCR，驗證了這些樞紐基因的表達水平。如圖7結(jié)果所示，BIRC5、NDC80、MAD2L1、CENPF、CCNB1、BUB1、BUB1B、KIF2C、CDC20和CDCA8在乳腺癌細胞中的表達均高于人正常乳腺上皮細胞；同時，雖然2個樞紐基因MAD2L1與BUB1的表達差異無統(tǒng)計學意義，但乳腺癌細胞中水平仍略高于正常乳腺上皮細胞。以上結(jié)果在細胞水平驗證了各樞紐基因表達。

圖7 qPCR檢測乳腺癌細胞樞紐基因的表達Fig 7 qPCR analysis of hub genes in breast cancer cells

3 討論

本研究利用基于生存分析的生物信息學方法，將基因表達數(shù)據(jù)與臨床生存分析相結(jié)合，篩選乳腺癌中的樞紐基因和關(guān)鍵通路。首先，我們從GEO數(shù)據(jù)庫下載了3個乳腺癌組織和癌旁組織的基因表達數(shù)據(jù)集，篩選出了585個差異表達基因，其中上調(diào)的基因有211個，下調(diào)的基因有374個。然后用Kaplan-Meier plotter數(shù)據(jù)庫篩選出了262個與乳腺癌患者總體生存期相關(guān)的差異表達基因。接著對這262個基因進行了GO功能分析，結(jié)果顯示這些基因主要與細胞分裂、細胞周期的調(diào)控以及染色體分離等生物學過程相關(guān)；KEGG通路分析結(jié)果表明，這些與總體生存期相關(guān)的差異表達基因主要富集在細胞周期、FoxO信號通路和卵母細胞減數(shù)分裂等通路上。此外我們發(fā)現(xiàn)人類T細胞白血病病毒1感染信號通路在乳腺癌中失調(diào)，但目前并無關(guān)于人類T細胞白血病病毒1感染信號通路與乳腺癌的報道；由于其在乳腺癌中的分子機制仍不明確，因此需要進一步研究。以上篩選出來的通路可以為乳腺癌發(fā)生發(fā)展的分子機制研究提供依據(jù)。

此外，我們通過構(gòu)建PPI網(wǎng)絡(luò)篩選出了10個乳腺癌的樞紐基因，分別為BIRC5、CDC20、NDC80、CENPF、MAD2L1、CDCA8、KIF2C、BUB1、CCNB1和BUB1B，經(jīng)Oncomine數(shù)據(jù)庫和HPA數(shù)據(jù)庫以及qPCR驗證，它們在乳腺癌中均高表達，并與乳腺癌患者較差的生存期相關(guān)。通過查閱文獻發(fā)現(xiàn)這些基因主要參與了細胞有絲分裂染色體的分離和細胞周期的調(diào)控。我們還嘗試在Oncomine數(shù)據(jù)庫中研究這些樞紐基因的表達是否與乳腺癌分子分型、分期及惡性程度有關(guān)；但由于數(shù)據(jù)庫中數(shù)據(jù)的一些局限，無法直接分析樞紐基因在不同分子分型乳腺癌的表達情況，因此無法確認樞紐基因的表達與分子分型的相關(guān)性。但我們發(fā)現(xiàn)樞紐基因在ER陰性乳腺癌和PR陰性乳腺癌的表達量分別高于ER陽性乳腺癌和PR陽性乳腺癌的表達量；在HER2陽性或陰性的乳腺癌中，這些基因的表達卻沒有明顯差別（數(shù)據(jù)未展示），提示樞紐基因的表達可能與分子分型有關(guān)。通過分析樞紐基因在不同分期的乳腺癌中的表達情況，我們發(fā)現(xiàn)這些基因在ⅠB期的表達量明顯低于其他分期的表達量，而BUB1B在Ⅲ期的表達量顯著高于其他分期（數(shù)據(jù)未展示）。MAD2L1、CDC20、CENPF、KIF2C、CCNB1、NDC80、BUB1、CDCA8和KIF2C在各分期的表達量沒有明顯差異。上述發(fā)現(xiàn)說明樞紐基因可能可以作為特定的分子分型和分期治療的判斷依據(jù)。

有 研 究 報 道，NDC80、MAD2L1、CDC20、BUB1、BIRC5和CCNB1在乳腺癌的發(fā)生發(fā)展中發(fā)揮了重要功能。NDC80基因編碼的NDC80蛋白參與構(gòu)成微管與動粒連接復合體，為染色體正常分離所必需，在細胞有絲分裂過程中起著至關(guān)重要的作用；NDC80的異常表達會造成染色體的異常分離，從而使染色體不穩(wěn)定，最終導致腫瘤的發(fā)生[14]。多項研究[15-16]表明NDC80在多種腫瘤中高表達并參與腫瘤的發(fā)生發(fā)展。NDC80抑制劑TAI-95可以在體外和體內(nèi)抑制乳腺癌腫瘤生長，NDC80有望成為乳腺癌的治療靶點[17]。MAD2L1和CDC20蛋白是紡錘體檢查點（spindle assembly checkpoint，SAC）復合體的組成成分，SAC負責確保姐妹染色單體的動粒和紡錘體結(jié)合并準確分離到2個子細胞中[18]。研究[19]顯示通過shRNA敲除MAD2L1可以抑制乳腺癌細胞的生長和侵襲能力。CDC20的高表達與乳腺癌患者較差的預后相關(guān)，且由它介導的BTG3相關(guān)核蛋白 （BTG3 associated nuclear protein，SMAR1）降解可以促進乳腺癌細胞遷移和侵襲[20-21]。BUB1是一類絲氨酸/蘇氨酸蛋白激酶，在有絲分裂和DNA損傷應答中發(fā)揮著重要的功能[22]。Han等[23]發(fā)現(xiàn)，敲低BUB1會降低腫瘤干細胞的潛能，BUB1可能成為針對乳腺癌干細胞的治療靶點。BIRC5編碼的蛋白是凋亡抑制（inhibitor of apoptosis，IAP）家族的成員之一，BIRC5具有明顯的抑制細胞凋亡的作用。Li等[24]認為BIRC5與乳腺癌患者進展及不良預后相關(guān)；研究[25-27]發(fā)現(xiàn)，BIRC5是miR-485-5p、ZIC家族成員1（Zic family member 1，ZIC1）和半乳凝素1（galectin-1）的靶點，可能在乳腺癌中發(fā)揮重要作用。CCNB1蛋白是調(diào)控細胞周期G2/M過渡階段所必需的；研究[28]表明，CCNB1蛋白與侵襲程度（腫瘤分級、腫瘤體積和淋巴結(jié)狀態(tài)）相關(guān)，是乳腺癌的預后因素之一。

目前KIF2C、CENPF、CDCA8和BUB1B這4個樞紐基因與乳腺癌相關(guān)的研究還比較少，具體分子機制有待進一步研究。KIF2C和CENPF蛋白參與有絲分裂的染色體分離。研究[29]報道，KIF2C在多種上皮型腫瘤中高表達，并與腫瘤的分級、分期和預后相關(guān)。CENPF基因編碼一種與著絲粒 - 動粒復合體相關(guān)的蛋白，其表達水平在有絲分裂前期增加，并在后期開始時發(fā)生蛋白水解[30-31]；CENPF的上調(diào)與乳腺癌的不良預后相關(guān)，可能是有臨床意義的治療靶點[32]。CDCA8是細胞周期分裂相關(guān)蛋白家族的一員，此蛋白受細胞周期調(diào)控，對染色質(zhì)誘導的微管穩(wěn)定和紡錘體形成是必需的；CDCA8的高表達與乳腺癌患者的低生存率和較差預后相關(guān)[33]。BUB1B基因編碼的蛋白是SAC復合體的重要成員，BUB1B表達增加是惡性程度高的乳腺癌的特征之一[34]。

綜上所述，本研究基于與總生存期相關(guān)的差異表達基因，篩選了乳腺癌中的樞紐基因和關(guān)鍵通路；基于生存期的樞紐基因可能更具有臨床意義，但這些樞紐基因的功能仍需進一步研究。