丙戊酸鈉抗癲癇療效相關(guān)樞紐基因鑒定

楊 炯，張明華，耿 淼，艾倫娜，田 磊，范 皎

(1.中國(guó)人民解放軍總醫(yī)院醫(yī)療保障中心，2.中國(guó)人民解放軍總醫(yī)院第二醫(yī)學(xué)中心老年醫(yī)學(xué)研究所，國(guó)家老年疾病臨床醫(yī)學(xué)研究中心，北京 100853)

癲癇是一種多因素參與的神經(jīng)系統(tǒng)疾病，全世界約有6 000多萬(wàn)人受到影響。當(dāng)前主要依靠藥物對(duì)癥治療，通過(guò)控制患者癲癇發(fā)作頻率改善患者的生活質(zhì)量。患者對(duì)于抗癲癇藥物(antiepileptic drug，AED)的反應(yīng)存在個(gè)體差異，約有40%～50%的患者對(duì)第一次AED單藥治療無(wú)效[1]，其中30%的患者有抗藥性[2]。有報(bào)道認(rèn)為，早期發(fā)病、治療前高發(fā)作頻率、隱原性癲癇、腦神經(jīng)解剖學(xué)異常等臨床因素是藥物反應(yīng)性差的原因。也有研究發(fā)現(xiàn)，AED藥物的代謝酶、轉(zhuǎn)運(yùn)體和靶點(diǎn)等分子的遺傳多態(tài)性可用于預(yù)測(cè)藥物治療的效果[3]。神經(jīng)環(huán)路由眾多未知細(xì)胞和物質(zhì)類(lèi)型構(gòu)成，癲癇的多因素特性涉及基因組和環(huán)境的復(fù)雜相互作用，任何遺傳或環(huán)境因素的改變最終都會(huì)影響到基因的表達(dá)和功能，需要超越單一通路的范疇去探究疾病的信號(hào)傳導(dǎo)機(jī)制，進(jìn)一步為開(kāi)發(fā)有效的藥物靶標(biāo)提供必要依據(jù)[4]。

目前，用于治療癲癇的藥物有二十多種，包括苯妥英鈉、卡馬西平和丙戊酸鈉等。其中丙戊酸鈉無(wú)論單用還是多藥聯(lián)用對(duì)多種類(lèi)型的癲癇都有較好的療效，目前是抗癲癇的一線用藥。丙戊酸鈉緩釋片屬于廣譜抗癲癇藥物，抗癲癇作用可能是通過(guò)競(jìng)爭(zhēng)性抑制γ-氨基丁酸轉(zhuǎn)移酶，提高腦內(nèi)γ-氨基丁酸的含量，抑制神經(jīng)元異常放電實(shí)現(xiàn)的。該藥的安全性比較高，能夠在維持藥物濃度的情況下，防止患者出現(xiàn)藥物不良反應(yīng)。但是丙戊酸鈉的抗癲癇與其他抗癲癇藥物一樣，其藥效和藥代動(dòng)力學(xué)都存在個(gè)體差異，原因尚不清楚。

加權(quán)基因共表達(dá)網(wǎng)絡(luò)分析(weighted gene co-expression network analysis, WGCNA)[5-6]依據(jù)基因的表達(dá)模式相似性構(gòu)建權(quán)重網(wǎng)絡(luò)，進(jìn)一步通過(guò)分析基因模塊與性狀數(shù)據(jù)之間的相關(guān)性得到與表型高度相關(guān)的基因集合，最后通過(guò)模塊內(nèi)部基因關(guān)聯(lián)分析得到模塊內(nèi)部的關(guān)鍵基因。WGCNA通過(guò)選擇合適的加權(quán)系數(shù)對(duì)基因之間的相關(guān)系數(shù)加權(quán)，使之滿足基因網(wǎng)絡(luò)近似服從無(wú)尺度網(wǎng)絡(luò)分布，擁有更好的統(tǒng)計(jì)功效，避免了大量的多重校正導(dǎo)致的假陰性結(jié)果。本研究通過(guò)對(duì)未接受過(guò)藥物治療的癲癇患者及接受丙戊酸鈉治療并出現(xiàn)不同療效患者的外周血樣本表達(dá)譜數(shù)據(jù)進(jìn)行WGCNA分析，為探尋造成丙戊酸鈉療效差異相關(guān)的基因提供理論依據(jù)。

1 材料與方法

1.1 數(shù)據(jù)獲取與預(yù)處理從GEO數(shù)據(jù)庫(kù)(Gene Expression Omnibus)中下載表達(dá)譜芯片數(shù)據(jù)，數(shù)據(jù)編號(hào)為GSE143272，測(cè)序平臺(tái)為Illumina Human HT-12 V4.0 expression beadchip。從該樣本集中選取34例此前未接受過(guò)藥物治療的癲癇患者(drug-free組)以及25例服用丙戊酸鈉并隨訪一年以上的患者(Valproate組)的外周血RNA樣本。其中又根據(jù)藥物治療隨訪過(guò)程中癲癇是否發(fā)作將25例患者分為丙戊酸鈉治療有效組患者(VA responders組)16例，無(wú)效組患者(VA non-responders組)9例。有效組指治療過(guò)程中無(wú)癲癇發(fā)作，無(wú)效組指用藥期間至少經(jīng)歷3次以上癲癇發(fā)作。數(shù)據(jù)集共包括13 165個(gè)基因探針表達(dá)數(shù)據(jù)，經(jīng)過(guò)合并和轉(zhuǎn)換篩選生成 9 924 個(gè)基因和59個(gè)樣本的表達(dá)矩陣，從9 924個(gè)基因的表達(dá)數(shù)據(jù)中篩選出SD值前5 000的基因進(jìn)行后續(xù)分析。

1.2 加權(quán)基因共表達(dá)網(wǎng)絡(luò)的構(gòu)建使用R語(yǔ)言中的WGCNA程序包構(gòu)建共表達(dá)網(wǎng)絡(luò)。對(duì)樣本進(jìn)行聚類(lèi)分析，根據(jù)樣本聚類(lèi)的距離鑒定是否存在離群樣本。共表達(dá)網(wǎng)絡(luò)符合無(wú)尺度網(wǎng)絡(luò)，即出現(xiàn)連接度為k的節(jié)點(diǎn)的對(duì)數(shù)lgk與該節(jié)點(diǎn)出現(xiàn)的概率的對(duì)數(shù)lg [P(k)]呈負(fù)相關(guān)，且相關(guān)系數(shù)應(yīng)>0.8。按照無(wú)尺度網(wǎng)絡(luò)標(biāo)準(zhǔn)篩選合適的連接函數(shù)加權(quán)參數(shù)，即軟閾值β。由軟閾值計(jì)算網(wǎng)絡(luò)拓?fù)渲丿B矩陣(topological overlap matrix, TOM)，通過(guò)動(dòng)態(tài)剪切樹(shù)法進(jìn)行模塊識(shí)別，設(shè)置每個(gè)基因網(wǎng)絡(luò)模塊最少的基因數(shù)目為30。

1.3 共表達(dá)模塊與臨床表型的相關(guān)性分析采用分層聚類(lèi)方法識(shí)別基因模塊，并用不同顏色來(lái)表示。計(jì)算每個(gè)模塊的特征向量基因ME(module eigengene)，降維處理模塊聚類(lèi)，并合并相似模塊，grey模塊是無(wú)法聚集到其他模塊的基因集合。通過(guò)計(jì)算模塊和臨床表型性狀相關(guān)系數(shù)，給出模塊和性狀之間的相關(guān)系數(shù)熱圖。本研究中關(guān)注的臨床表型有年齡(age)，性別(sex)，癲癇種類(lèi)(epilepsy type)、是否接受丙戊酸鈉治療(treatment)、丙戊酸鈉治療是否有效(VA response)。通過(guò)模塊的特征向量與性狀的相關(guān)系數(shù)以及模塊顯著性P<0.05篩選出與臨床表型顯著相關(guān)的基因模塊。計(jì)算模塊內(nèi)的基因表達(dá)與性狀的相關(guān)性(gene significance，GS)和某個(gè)基因表達(dá)與模塊內(nèi)基因主成分表達(dá)的相關(guān)系數(shù)(module membership，MM)，通過(guò)設(shè)置GS、MM、q.weighted的取值范圍對(duì)networkScreening函數(shù)計(jì)算得到的基因列表進(jìn)行篩選，從而識(shí)別和鑒定出關(guān)鍵樞紐基因(hub genes)。

1.4 網(wǎng)絡(luò)可視化及功能富集分析選取與臨床表型顯著相關(guān)的基因模塊，利用Cytoscape軟件繪制網(wǎng)絡(luò)圖。基于基因本體論(gene ontology, GO)數(shù)據(jù)庫(kù)對(duì)與性狀高度相關(guān)的模塊分別進(jìn)行基因功能注釋分析，采用Fisher檢驗(yàn)篩選出模塊基因顯著性富集的GO條目；基于Kyoto Encyclopedia of Genes and Genomes(KEGG)數(shù)據(jù)庫(kù)對(duì)核心模塊中基因進(jìn)行信號(hào)通路(pathway)富集分析。

2 結(jié)果

2.1 加權(quán)基因共表達(dá)網(wǎng)絡(luò)的構(gòu)建首先對(duì)所有樣本的基因表達(dá)值進(jìn)行聚類(lèi)分析并作圖(Fig 1A)。通過(guò)WGCNA算法，根據(jù)無(wú)尺度網(wǎng)絡(luò)擬合指數(shù)和平均連接度，計(jì)算并選取β=6作為本數(shù)據(jù)集的軟閾值(Fig 1B)，計(jì)算基因間的鄰接矩陣和TOM矩陣，并根據(jù)TOM矩陣構(gòu)建基因間的分層聚類(lèi)樹(shù)，基于動(dòng)態(tài)剪切樹(shù)的方法把基因分成12個(gè)模塊，分別用12種顏色矩形表示，縱坐標(biāo)為基因占比(Fig 1C)。各個(gè)模塊包含基因個(gè)數(shù)如下：black(177)，blue(675)，brown(540)，green(199)，greenyellow(59)，magenta(136)，pink(175)，purple(86)，red(181)，turquoise(1897)，yellow(322)，grey(553)，其中無(wú)法聚類(lèi)到其他任何模塊中的基因放入grey模塊，在后續(xù)分析中將grey模塊移除。

Fig 1 Construction of co-expression network

2.2 共表達(dá)模塊與臨床表型的相關(guān)性分析計(jì)算不同基因模塊與臨床表型之間的關(guān)系，繪制共表達(dá)模塊與表型的相關(guān)性熱圖(Fig 2A)。與是否接受丙戊酸鈉治療(treatment)表型相關(guān)性最強(qiáng)的模塊是yellow模塊(r=0.57，P<0.000 1)，yellow模塊中的基因總體上與治療方式表型正相關(guān)，即yellow模塊整體上在服用丙戊酸鈉治療組患者(valproate組)中高表達(dá)，在此前未接受過(guò)藥物治療的患者(drug-free 組)中低表達(dá)。與丙戊酸鈉治療有效(VA response)表型相關(guān)性最強(qiáng)的模塊是blue模塊(r=-0.53，P<0.000 1)，blue模塊中的基因總體上與丙戊酸鈉治療有效表型負(fù)相關(guān)，即blue模塊整體上在丙戊酸鈉治療有效組患者(VA responders組)中低表達(dá)，在無(wú)效組患者(VA non-responders組)中高表達(dá)。對(duì)各個(gè)模塊進(jìn)行層次聚類(lèi)和用熱圖分析各模塊之間的相關(guān)性(Fig 2B)。

Fig 2 Correlation between co-expression modules and clinical traits

將yellow和blue模塊作為關(guān)鍵模塊進(jìn)行GS和MM分析，yellow模塊的GS與MM的相關(guān)系數(shù)r=0.42，P<0.000 1(Fig 3A)，blue模塊的GS與MM的相關(guān)系數(shù)r=0.29，P<0.000 1(Fig 3B)。進(jìn)一步分析yellow模塊和blue模塊的eigengene表達(dá)熱圖(Fig 3C、3D)及yellow模塊和blue模塊基因與各表型之間的熱圖和聚類(lèi)圖(Fig 3E、3F)。

Fig 3 Correlation between yellow and blue modules with clinical traits

在兩個(gè)模塊中通過(guò)基因顯著性(GS)和模塊成員關(guān)系(MM)篩選yellow和blue模塊中關(guān)鍵樞紐基因。設(shè)置3個(gè)篩選標(biāo)準(zhǔn)：|GS|>0.4，|MM|>0.2，q.weighted<0.05。yellow模塊篩選出14個(gè)樞紐基因(Tab 1)，blue模塊篩選出10個(gè)樞紐基因(Tab 2)。

Tab 1 Hub genes screened in yellow module

Tab 2 Hub genes screened in blue module

2.3 網(wǎng)絡(luò)可視化及功能富集分析分別根據(jù)兩個(gè)模塊中篩選出的樞紐基因權(quán)重做出共表達(dá)網(wǎng)絡(luò)基因間的相互作用關(guān)系圖(Fig 4)，其中節(jié)點(diǎn)形狀的大小代表與該節(jié)點(diǎn)連接的邊的數(shù)量，邊的寬度代表兩節(jié)點(diǎn)間連接的權(quán)重大小。

Fig 4 Hub gene mapA: Hub gene map of yellow module; B： Hub gene map of blue module.

分別對(duì)兩個(gè)模塊中的樞紐基因進(jìn)行GO富集分析(Fig 5A、5B)，發(fā)現(xiàn)yellow模塊的基因功能主要富集于：免疫反應(yīng)調(diào)控(GO:0050776)、細(xì)胞表面受體信號(hào)通路(GO:0007166)、T細(xì)胞受體信號(hào)通路(GO:0050852)、MAPK信號(hào)轉(zhuǎn)導(dǎo)(GO:0000188)；blue模塊的基因功能主要富集于：細(xì)胞增殖調(diào)控(GO:0042127)、脂質(zhì)糖反應(yīng)(GO:0032496)、GTPase活性激活(GO:0090630)、吞噬作用(GO:0006909)、神經(jīng)元凋亡(GO:0051402)等生物學(xué)過(guò)程。此外，KEGG通路分析顯示，yellow模塊基因主要富集于：自然殺傷細(xì)胞介導(dǎo)的細(xì)胞毒性(hsa04650)、Ras信號(hào)通路(hsa04014)、細(xì)胞黏附分子(hsa04514)、初級(jí)免疫缺陷(hsa05340)(Fig 5C)；blue模塊基因主要富集于：NF-κB信號(hào)通路(hsa04064)、細(xì)胞周期(hsa04110)、抗生素生物合成(hsa01130)等通路(Fig 5D)。

Fig 5 GO and KEGG enrichment analysis of yellow and blue modules

3 討論

癲癇的發(fā)作是一種由多因素共同參與的復(fù)雜過(guò)程，機(jī)制較為復(fù)雜，遺傳因素、腦部疾病等均可誘發(fā)該疾病。從該病臨床特點(diǎn)上分析，發(fā)作一般沒(méi)有征兆，有可能導(dǎo)致患者溺水、摔倒和燙傷等意外情況，威脅患者生命安全，影響生活質(zhì)量。目前盡管已有20多種藥物應(yīng)用于抗癲癇臨床治療，但仍有約1/3的患者藥物治療無(wú)效。分析單個(gè)基因研究癲癇發(fā)病機(jī)制難以實(shí)現(xiàn)突破，基于網(wǎng)絡(luò)的分析方法可能有助于發(fā)現(xiàn)疾病相關(guān)的基因網(wǎng)絡(luò)。在傳統(tǒng)的基因水平研究更多關(guān)注于強(qiáng)效應(yīng)基因，本研究采用的WGCNA算法補(bǔ)充了對(duì)弱效應(yīng)基因的分析。WGCNA算法構(gòu)建的基因網(wǎng)絡(luò)關(guān)系服從近似無(wú)尺度網(wǎng)絡(luò)分布，因此相較于常規(guī)的聚類(lèi)方法，更具有生物學(xué)數(shù)據(jù)特性，很好還原基因在生物學(xué)過(guò)程中的作用，有助于鑒定出與特定臨床表型相關(guān)的重要基因模塊和關(guān)鍵樞紐基因。

在本研究中，對(duì)癲癇患者外周血樣本的表達(dá)譜數(shù)據(jù)集進(jìn)行了WGCNA分析，根據(jù)丙戊酸鈉的治療效果，將患者分為丙戊酸鈉治療有效組患者和無(wú)效組患者，同時(shí)和未經(jīng)藥物治療的癲癇患者一起構(gòu)建共表達(dá)網(wǎng)路分析基因表達(dá)差異。在共表達(dá)網(wǎng)絡(luò)分析中，識(shí)別并聚類(lèi)成12個(gè)顏色模塊，并對(duì)各模塊進(jìn)行基因與臨床表型的相關(guān)性分析，其中yellow模塊中的基因與治療方式表型的相關(guān)性最為明顯，也即yellow模塊中的基因在接受丙戊酸鈉治療后的患者中會(huì)發(fā)現(xiàn)明顯變化；blue模塊中的基因與丙戊酸鈉療效表型的相關(guān)性最為顯著，也即blue模塊的基因在丙戊酸鈉有效患者和無(wú)效患者中表達(dá)存在明顯差異。有研究報(bào)道，MAPK信號(hào)通路調(diào)控神經(jīng)元的基因表達(dá)調(diào)控癲癇和認(rèn)知障礙中的突觸興奮，與癲癇的發(fā)作有密切聯(lián)系[7]。p38 MAPK主要由炎癥因子和環(huán)境壓力激活，抑制p38 MAPK可以減少癲癇誘導(dǎo)的錐體細(xì)胞缺失、細(xì)胞變性、神經(jīng)元損傷和海馬神經(jīng)元退化[8]。本研究發(fā)現(xiàn)，yellow模塊中的基因部分富集到MAPK信號(hào)失活通路，提示丙戊酸鈉可能通過(guò)抑制MAPK信號(hào)通路來(lái)發(fā)揮抗癲癇的作用。NF-κB廣泛存在于神經(jīng)細(xì)胞中，參與免疫應(yīng)激反應(yīng)、炎癥反應(yīng)及細(xì)胞的增殖與凋亡。在癲癇患兒的外周血單個(gè)紅細(xì)胞中NF-κB活化與癲癇發(fā)作的嚴(yán)重程度成正相關(guān)[9]。也有研究顯示NF-κB的活化與癲癇易感性有關(guān)。本研究分析顯示與丙戊酸鈉療效明顯相關(guān)的blue模塊基因富集在NF-κB通路，癲癇患者體內(nèi)NF-κB通路的活化可能影響丙戊酸鈉的藥物敏感性。此外，yellow模塊和blue模塊基因明顯富集在T細(xì)胞受體信號(hào)通路和T細(xì)胞激活機(jī)制，這可能與丙戊酸鈉對(duì)組蛋白乙酰化的修飾能力有關(guān)；有研究在腫瘤細(xì)胞治療過(guò)程中，利用丙戊酸鈉的這種活性增強(qiáng)治療效果[10]。綜合來(lái)看，丙戊酸鈉抗癲癇的藥理作用與其對(duì)免疫反應(yīng)的調(diào)控機(jī)制存在密切的關(guān)聯(lián)。

根據(jù)共表達(dá)網(wǎng)絡(luò)的顯著性，本研究在yellow模塊中篩選出包括S1PR5、GNLY、CD160在內(nèi)的14個(gè)樞紐基因，在blue模塊中篩選出包括MAGED1、FBXO31在內(nèi)的10個(gè)下調(diào)樞紐基因。S1PR5是一種淋巴細(xì)胞表面的G-蛋白偶聯(lián)受體，主要表達(dá)于NK細(xì)胞表面，S1PR5表達(dá)缺失影響穩(wěn)態(tài)及炎癥條件下小鼠體內(nèi)自然殺傷細(xì)胞(natural killer, NK)的分布[11]。有研究發(fā)現(xiàn)，癲癇患者發(fā)作期間外周血中NK細(xì)胞活性低下，但發(fā)作后外周血中NK細(xì)胞水平增高[12]。早期研究報(bào)道，部分癲癇患者服用苯妥英鈉后NK細(xì)胞活性降低，而服用卡馬西平的癲癇患者NK細(xì)胞活性增加[13]。目前，癲癇及抗癲癇藥物對(duì)NK細(xì)胞活性的研究結(jié)果尚不完全清楚。本研究提示，服用丙戊酸鈉的癲癇患者相比于未經(jīng)藥物治療的患者S1PR5基因水平出現(xiàn)明顯變化，有可能會(huì)調(diào)控NK細(xì)胞的分布或功能。

MAGED1是黑色素瘤相關(guān)抗原家族成員之一，廣泛表達(dá)于各組織中，參與細(xì)胞周期調(diào)控、細(xì)胞凋亡、分化、細(xì)胞黏附等多個(gè)生物學(xué)過(guò)程[14]。已有研究發(fā)現(xiàn)MAGED1調(diào)節(jié)脂肪前體細(xì)胞發(fā)育進(jìn)而參與糖脂代謝調(diào)控途徑[15]。FBXO31基因編碼的蛋白屬于F-box蛋白家族的一員，可以通過(guò)靶向介導(dǎo)細(xì)胞周期蛋白CyclinD1的泛素化和降解來(lái)抑制細(xì)胞增殖和細(xì)胞周期進(jìn)程[16]。此外，F(xiàn)BXO31通過(guò)降解MDM2在DNA損傷修復(fù)過(guò)程中發(fā)揮重要作用[17]。

丙戊酸鈉是抑制癲癇患者發(fā)作常用藥物，但是患者之間存在較大的個(gè)體差異。臨床上通常從藥代動(dòng)力學(xué)的角度，利用治療藥物檢測(cè)技術(shù)，檢測(cè)其穩(wěn)態(tài)血藥濃度調(diào)整給藥劑量的方式，實(shí)現(xiàn)個(gè)體化的用藥。如果仍然不能有效控制癲癇發(fā)作，一般聯(lián)合使用拉莫三嗪、加巴噴丁、卡馬西平、托吡酯、苯巴比妥等藥物。目前較少?gòu)乃幮?dòng)力學(xué)角度考慮提高個(gè)體化治療水平。本研究通過(guò)WGCNA分析，從藥效學(xué)的角度發(fā)現(xiàn)了與丙戊酸鈉療效發(fā)生高度相關(guān)的模塊，包含MAGED1、FBXO31等關(guān)鍵樞紐基因，可能通過(guò)調(diào)節(jié)NF-κB信號(hào)通路、免疫調(diào)控等生物學(xué)過(guò)程影響癲癇患者治療效果。本研究為揭示癲癇的耐藥機(jī)制、進(jìn)一步提升丙戊酸鈉的臨床療效、減少患者的不良反應(yīng)提供了藥效學(xué)理論依據(jù)。