亞洲人阿爾茲海默癥miRNA-mRNA網絡的生物信息學分析

楊澤若，張 燚，胡 柳，溫 軼

(浙江養生堂天然藥物研究所， 杭州 310024)

阿爾茲海默癥(Alzheimer’s disease，AD)是一種慢性神經退行性疾病，在臨床上主要表現為記憶喪失與日常行動功能障礙。據阿爾茲海默協會(Alzheimer’s Association)報道[1]，截至2018年， 阿爾茲海默癥是造成美國65歲以上人群死亡的第五大原因。根據中國疾病預防控制中心發布的數據顯示[2]，作為一種高致殘率、影響獨立生活能力的疾病，AD所導致的死亡人數在中國所有疾病中排名第五。AD的病因眾說紛紜，但確切的發病機制尚未闡明，迄今也無特效治療或逆轉疾病進展的藥物。

隨著高通量測序技術和全基因組關聯分析研究的發展，多個與AD患病相關的基因被發現，其中已證實有三個基因——淀粉樣前體蛋白(APP)、早老素1(PSEN1)和早老素2(PSEN2)是家族性阿爾茲海默癥(Familial Alzheimer disease，FAD)的致病基因。目前針對AD的研究主要集中在高加索人群，而對亞洲人群以及中國人群的研究相對較少。Jia等[3]招募了404個家系的1 330例AD或輕度認知障礙(Mild cognitive impairment，MCI)患者并檢測PSENs / APP突變，鑒定了11個新的突變位點，新的PSENs / APP突變表明中國人與其他種族之間AD發病機制可能存在異質性。Zheng等[4]對漢族FAD患者進行了全基因組測序，鑒定出了漢族人特有的基因外顯子錯義突變rs3792646，表明Complement C7是漢族人患AD的新風險基因。Han等人指出[5]，CLU rs11136000多態性與白種人AD發病顯著相關，而它與亞洲人群AD發病并沒有顯著性關聯。由此我們猜測，不同人種之間的AD遺傳/發病機制存在差異，因此有必要對研究相對匱乏的亞洲人群進行研究，以挖掘亞洲人群特有的AD致病基因、生物標志物以及風險因素。

microRNA(miRNA)是一類保守的非編碼RNA小分子，它們主要通過抑制靶mRNA翻譯或促使其降解實現對靶基因的轉錄后調控[6]。miRNA在神經發育、組織分化和突觸形成過程中起著重要作用[7]。此外，越來越多的證據表明血漿miRNA含量可能作為潛在的疾病標志物，比如血漿miRNA-206水平升高可預測認知能力下降和癡呆程度[8]。目前，有關miRNA 在亞洲人群AD發生過程中的作用研究甚少，而且關于亞洲人AD的表達譜研究往往只針對mRNA 或者miRNA表達譜，并沒有將二者進行關聯分析，結合miRNA 及mRNA表達譜的亞洲人AD系統分析尚未見報道。

因此，收集現有公開的亞洲人AD與正常人mRNA/miRNA表達譜數據，通過生物信息學手段，篩選核心差異基因與差異miRNAs，構建miRNA-mRNA的差異表達網絡。從轉錄組水平系統的對亞洲人群可能的AD發生機制進行分析，為亞洲人AD的診療提供新的潛在靶標。

1 材料與方法

分析流程主要包含以下環節：數據收集與下載，數據集預處理與單個數據集差異分析，可合并腦區的薈萃分析(Meta分析)和合并分析，核心差異基因選取，潛在靶基因預測，核心miRNA-mRNA網絡構建，差異基因/靶基因富集分析，靶基因相關蛋白的互作分析(見圖1)。

圖1 分析流程圖Fig.1 Analysis flow chart

1.1 基因表達數據收集

亞洲人AD患者與正常人的表達譜數據來自基因表達綜合數據庫(Gene Expression Omnibus, GEO)[9](https://www.ncbi.nlm.nih.gov/geo/)。通過限定關鍵詞“Alzheimer AND HAN”或者“Alzheimer AND ASIAN”， 選取了包含額葉 (Frontal Cortex, FC) 、顳葉 (Temporal Cortex，TC)、海馬體 (Hippocampus，HP) 與內嗅皮質 (Entorhinal Cortex，EC) 等AD發病相關腦區的數據。截止2020年6月，一共收集4個表達譜數據集，分別是GSE131617[10]、GSE36980[11]、GSE139384[12]和GSE120584[13]。為確保數據質量，剔除了小于50歲的樣本以及在主成分分析中明顯離群的樣本。最終一共獲得了212個AD患者腦組織mRNA數據(76個FC，71個TC，58個EC以及7個HP)與92個正常人(Control，CT)腦組織mRNA數據(34個FC，35個TC，13個EC以及10個HP)。此外，miRNA血清數據集中包含1 021個AD以及288個CT。4個數據集的具體信息(見表1)。

表1 從GEO收集的亞洲人AD患者與正常對照的表達譜數據集Table 1 Expression data sets for AD patients and normal controls in Asian populations collected from GEO

1.2 數據集下載、預處理以及單個數據集差異分析

首先用R包GEOquery[14]下載數據集，利用探針注釋的文件對表達矩陣進行探針ID轉換。接著利用R包limma[15]的線性模型，對單個表達譜數據進行歸一化處理、性別校正與差異分析：其中miRNA數據集進行Benjamini-Hochberg(BH)法[16]校正，P校正值 < 0.05為顯著差異；mRNA數據集取P值 < 0.05 為顯著差異。由于GSE131617數據集根據布拉克分期(Braak Staging)[17]將AD樣本分為多組，將該數據集中布拉克分期為0的樣本劃分為正常人，將其余布拉克分期樣本統一合并為AD患者。

1.3 FC/TC腦區多個mRNA數據集整合分析

由于多個數據集中均包含FC和TC兩個腦區的數據，同時采用以下兩種算法進行整合分析：1)合并同腦區多個數據集，去除批次效應后進行差異分析；2)meta分析中的魯棒排名聚集算法(Robust Rank Aggregation，RRA)[18]。

1.3.1 FC/TC腦區mRNA數據集的合并分析

將不同研究的芯片數據合并進行歸一化，可以移除不同研究的批次效應，同時保留真實的生物學差異[19]。通過R包sva中的ComBat函數[20]，將三個mRNA數據集的表達矩陣進行合并，并且去除批次效應，然后用limma的線性模型對合并之后的表達矩陣進行性別校正，最后對合并之后的表達矩陣進行差異分析。

1.3.2 FC/TC腦區mRNA數據集的RRA分析

對單個數據集中FC或TC腦區得到的差異基因(P< 0.05)，按照絕對值表達差異倍數(|logFC|)從大到小進行排列，并進行RRA運算，以RRAscore< 0.05 作為篩選閾值，獲得三個mRNA數據集數據中排名相對靠前的差異基因。

1.4 核心差異基因選取

對于存在于多個數據集中的腦區(FC/TC)，通過1.3中兩種方法的交集來獲得核心差異基因。對于僅存在于單個數據集的腦區(EC/HP)，通過進一步降低P值與增加絕對值表達差異倍數(|logFC|)的篩選閾值來獲得核心差異基因。

1.5 差異miRNA潛在靶基因預測

利用miRWalk在線工具[21](http://zmf.umm.uniheidelberg.de/apps/zmf/mirwalk/)進行差異miRNA的靶基因預測。選擇靶基因3端UTR作為miRNA的作用區域，并選擇PicTar5[22]、 miRWalk2.0、 TargetScan[22]、DIANA[23]、miRDB[24]以及miRanda[25]等6個靶基因數據庫，選取至少被其中2個靶基因數據庫收錄的靶基因作為潛在靶基因。

1.6 靶基因選取、miRNA-mRNA網絡構建、靶基因蛋白互作分析

將1.5得到的潛在靶基因與1.4得到的各腦區核心差異基因進行交集，得到四個腦區的靶基因。將靶基因與對應的差異miRNA以及差異表達倍數等導入Cytoscape軟件[26]，構建四個腦區的miRNA-mRNA差異表達網絡。用geneMANIA[27](https://genemania.org/)數據庫，對靶基因進行蛋白互作分析。

1.7 基因本體論(GO)以及京都基因組百科全書數據庫(KEGG)通路富集分析

通過R包clusterprofiler[28]對四個腦區的差異基因與靶基因進行GO以及KEGG通路富集分析。

2 結 果

2.1 AD患者與正常對照組miRNA/mRNA差異表達結果

2.1.1 miRNA差異表達結果

以P校正值< 0.05和|logFC| > 0.5作為差異miRNA的選擇標準，共篩選出5個差異miRNAs， 其中4個(hsa-mi-22、hsa-miR-24、hsa-miR-125b-1-3p、hsa-mi-125-3p)在AD組下調，1個(hsa-miR-208a-5p)在AD組上調。

2.1.2 mRNA差異表達結果(差異基因與核心差異基因)

差異基因與核心差異基因的具體篩選要求及差異個數統計數據(見表2)。對于存在于多個數據集中的腦區(FC/TC)，采用方法1.31的合并分析進行差異基因篩選，以上差異基因再與1.32中RRA分析結果取得交集獲得核心差異基因；對于僅存在于單個數據集的腦區(EC/HP)，通過常規篩選條件獲得差異基因，并通過提高篩選條件來獲得核心差異基因。

表2 差異基因與核心差異基因的篩選條件與個數統計Table 2 Thresholds for DEGs and KDEGs selections and summary statistics

如表2所示，額葉、顳葉、海馬體與內嗅皮質中分別存在346、662、461與234個差異基因(見圖2)，研究發現NPTX2、SLC14A1與GJA1在四個腦區中均為差異基因，且上述基因在AD患者(與正常對照相比)的四個腦區中呈現出相同的上調或下調表達趨勢。額葉、顳葉、海馬體與內嗅皮質中分別存在15、29、63與20個核心差異基因(見表2)。SLC14A1是唯一一個在四個腦區中都被選入核心差異基因的基因。此外，NPAS4同時出現在FC、TC與EC腦區；GMPR、ITGB4、APLNR與NPTX2同時出現在FC、TC與HP腦區。

2.2 miRNA潛在靶基因預測和miRNA-mRNA網絡構建

利用miRWalk等6個數據庫，對5個差異miRNAs進行潛在靶基因預測，并與前文所述四個腦區的核心差異基因進行交集，獲得靶基因。共有4個差異miRNAs的潛在靶基因與核心差異基因存在交集，分別是hsa-mi-22、hsa-miR-24、hsa-miR-125b-1-3p、hsa-mi-125-3p。4個差異miRNAs與四個腦區靶基因的調控關系見miRNA-mRNA網絡(見圖3)。

2.3 四個腦區靶基因統計及GeneMANIA的蛋白互作分析

額葉、顳葉、海馬體與內嗅皮質中分別存在6、10、25與7個靶基因(見表3)。其中，SLC14A1同時出現在四個腦區的miRNA-mRNA網絡中，NPTX2出現在FC、TC與HP腦區的網絡中，AQP1與ANTXR1出現在FC與TC腦區的網絡中，TXNIP出現在TC與HP腦區的網絡中。此外，DDIT3、FOSB為FC腦區特有靶基因；HLA-DMA、TPT1、BRMS1、DLG4、ITPKA為TC腦區特有靶基因；CD44、EMP1、PRRX1、HCN1、KCNH5、HTR3B、CBLN4、RAB15、GFRA2、GABRA3、HPCA、PTPN3、IL12RB2、CACNG3、CTXN3、LRRC2、NCALD、HTR2A、GABRA1、RGS4、PCSK1、MET為HP腦區特有靶基因；FOS、SOD2、PDYN、NAMPT、PLAC8、MT1M為EC腦區特有靶基因。靶基因的GeneMANIA結果表明，四個腦區中的靶基因間均存在復雜的相互作用，包括共表達相關性、共享蛋白結構域、基因相互作用、物理相互作用與共同通路等等。值得注意的是，HP靶基因數目最多、差異倍數更大、具有更復雜的蛋白相互作用關系、且多在HP中特異表達；EC靶基因與其他腦區靶基因交集最少，SLC14A1為唯一交集。

表3 靶基因統計以及geneMANIA蛋白互作分析結果Table 3 Target genes and geneMANIA results of the PPI network analyses

2.4 四個腦區中差異基因/靶基因GO與KEGG通路富集結果

四個腦區差異基因的GO與KEGG富集結果(見圖4)。結果表明FC、TC與HP腦區的GO富集結果比較接近，主要集中在神經遞質傳導、化學突觸傳遞等神經信號傳遞相關通路，而EC腦區的GO富集結果更偏向于免疫相關，如中性粒細胞脫顆粒等。KEGG分析中， FC與TC腦區都富集到了神經性退行性疾病通路，如亨廷頓氏病、漸凍癥、帕金森癥等，HP富集到了神經遞質傳導相關通路，如神經活性配體受體相互作用信號通路、鈣離子信號通路、GABA能突觸信號、谷氨酸能突觸信號，而EC主要富集在炎癥、免疫相關通路。

四個腦區靶基因的GO與KEGG富集結果(見圖5a,5b)。從GO分析結果來看，四個腦區各自的靶基因可被同時富集到轉運相關的生物學過程中；FC、TC與HP的靶基因被同時富集到神經遞質受體活性調控的生物學過程；HP的靶基因可被特異性地富集到離子跨膜運輸調控過程中；EC的靶基因可被特異性地富集到低溫應激、衰老、神經元死亡、活性氧應激等4個生物學過程中。從KEGG分析結果來看，四個腦區的靶基因未能富集到同一通路上，但FC、TC與EC腦區的靶基因能被同時富集到可卡因成癮通路。

3 討 論

本研究收集了GEO數據庫中亞洲人AD患者與正常人mRNA/miRNA數據，通過薈萃分析與合并分析，從轉錄組水平系統地對亞洲人群可能的AD發病機制進行了分析，首次構建了亞洲人四個腦區AD的miRNA-mRNA差異表達網絡。

目前有多項生信研究利用已有AD表達譜數據對AD機制進行了系統性研究。Xu等[29]運用合并分析方法，整合了20個表達譜數據，構建了AD患者四個腦區的mRNA轉錄譜，發現了YAP1等AD上游調控因子。Shi等[30]通過分析1個RNA-seq數據集，鑒定了AD患者上游調控的lncRNA，構建了lncRNA-mRNA共表達網絡。Moradifard等[31]分析了1個miRNA以及6個mRNA芯片數據，通過薈萃分析發現了AD中重要的mRNA/miRNA相互作用關系。以上工作往往只采用一種整合分析方法，如薈萃分析或合并分析，或者只進行了單組學分析。本研究對多數據集腦區的數據同時應用了兩種整合分析方法，可以更加準確地得到多數據集腦區中的差異基因。此外，上述生信研究的表達譜數據均為高加索人群，而不同人種之間的AD遺傳/發病機制可能存在差異，因此本研究旨在收集亞洲人AD表達譜數據，構建亞洲人AD的miRNA-mRNA差異表達網絡。

發現5個差異miRNAs，已有文獻報道miR-22與miR-24在AD病人中低表達，與本研究中的發現一致：Jovicic等[32]指出miR-22在阿爾茲海默癥患者中低表達；Hu等[33]通過薈萃分析發現miR-24在AD患者的腦脊液中低表達。盡管Ylmaz等[34]指出hsa-mi-125-3p和hsa-miR-125b-5p在AD患者中(相比于正常人)并沒有顯著差異，但是其研究樣本較少(AD 172個，CT 109個)，且亞洲人群占比較少。因此，這兩個差異miRNAs及其靶基因在亞洲人群AD發生發展中的作用值得進一步的研究。此外，本研究發現 NTPX2、GJA1與SLC14A1在四個腦區中均為差異基因。Xiao等[35]表明NPTX2在 AD患者皮質中大量減少，且與患者認知能力和海馬體積密切相關。Kajiwara等[36]指出星形膠質細胞特異性表達的GJA1在AD病人大腦上調，并發現GJA1與AD淀粉樣蛋白、tau病理以及認知功能密切相關。Kerstin等[37]指出，SLC14A1的表達在AD病人和小鼠模型中均顯著上調。上述基因在本研究的表達趨勢與文獻報道一致，進一步證明了本研究分析策略和結果的可靠性。

SLC14A1(尿素轉運蛋白1)是唯一在四個腦區中共有的靶基因。通過蛋白互作分析，發現SLC14A1和AQP1(水通道蛋白1)存在蛋白相互作用，同時AQP1在AD患者中表達上調，且四個腦區靶基因的GO都富集到液體轉運等生物學過程。有研究表明SLC14A1在亨廷頓舞蹈病中顯著增高，從而破壞大腦中的尿素平衡，導致紋狀體神經元大量死亡[38]。此外，Hansmannel等[39]通過AD患者額葉皮層的轉錄組測序，指出腦中尿素平衡與AD發病相關。因此猜測SLC14A1與AQP1在AD患者中的表達上調與相互作用，可能會引起腦中尿素循環與水循環紊亂，從而造成神經損傷并誘導AD發病。此外，通過艾倫大腦圖譜(Allen Brain Map)[40]中的單細胞轉錄組數據，發現AQP1、SLC14A1都在星形膠質細胞(Astro L1-6 FGFR3)中特異表達。 Misawa等[41]也指出，星形膠質細胞中AQP1上調可能與AD患者腦中Aβ沉積相關。因此本研究進一步推測，AD患者星形膠質細胞SLC14A1、AQP1的差異表達，將破壞星形膠質細胞的水平衡與尿素平衡，從而引起細胞腫脹破裂，最終影響神經元功能。

通過整合四個腦區靶基因和miRNA-mRNA的差異網絡構建，對單腦區的特異性靶基因進行了分析，特別是HP腦區，其靶基因數目最多、差異倍數更大、多數在海馬中特異表達、并且主要是配體門控性離子通道相關基因，包括GABRA1(A型γ-氨基丁酸受體亞基Alpha1)，GABRA3(A型γ-氨基丁酸受體亞基Beta3)，HTR3B(5-羥色胺受體3B)與HTR2A(5-羥色胺受體2A)等。Agenor等[42]表明GABA受體 α1、β3亞型在AD患者的海馬中表達下調。Garcia-Alloza等[43]發現膽堿能與5-羥色胺能系統之間的不平衡，可能與AD的認知障礙相關。此外，Pedro 等[44]證明Cbln4參與形成并維持GABA能連接，敲除Cbln4后神經元GABA能連接大大減少，且Cbln4在AD小鼠的海馬區中表達顯著下調。上述基因在本研究中的HP腦區中均表達下調，因此，本研究推斷亞洲人群AD患者的海馬體GABA受體亞型與5-羥色胺能受體亞型表達下調，可能破壞興奮性能與抑制性能神經元之間信號平衡，最終造成海馬體中神經元的功能障礙。

此外，EC腦區靶基因與其他腦區靶基因交集只有SLC14A1，且其GO/KEGG富集通路與其他三個腦區有較大區別。Flynn等[45]表明SOD2(線粒體超氧化物歧化酶2)可以將超氧化物還原為過氧化氫，保護神經元細胞免受谷氨酸誘導的氧化應激和細胞毒性。Xing等[46]指出，NAMPT(煙酰胺磷酸核糖基轉移酶)減少導致NAD +減少，從而導致NAD+ / NADH的比率降低與線粒體功能障礙。在本研究中，SOD2、NAMPT等基因在EC中下調，因此本研究猜測亞洲人群AD患者EC中SOD2與NAMPT的下調可能引起線粒體功能障礙，繼而造成突觸受損、神經細胞凋亡，最終促使AD的病情發展。

本研究主要基于GEO數據庫中的微陣列芯片數據，因此具有一些局限性。首先，微陣列芯片數據受制于其技術的局限性，例如技術動態范圍較窄、靈敏度較低等。其次，本研究使用的芯片數據來自不同條件下的多個研究，樣本量不均并且許多患者臨床細節未知(如用藥情況等)，增加了本研究不確定性。盡管本研究采用兩種整合分析方法，也無法完全消除各項研究的異質性。此外，本研究受限于亞洲人數據集的樣本量，某些腦區(如海馬體)樣本量較少，其結果可能不具備普適性。值得注意的是，上述芯片數據的基因表達量為單腦區所有細胞的平均值，并沒有達到單細胞精確度，因此我們無法直接證實SLC14A1是否在AD患者的星形膠質細胞特異差異表達，并影響AD患者星形膠質細胞的功能。最后，本研究沒有詳細探究FC與TC腦區特異性變化的靶基因，可能忽略了這些靶基因在上述腦區的特異性作用。本研究雖然對SLC14A1在亞洲人AD中的潛在作用進行了推測，但仍然需要構建相應的細胞和動物模型進一步研究其確切的作用機制。此外，本研究希望通過單細胞轉錄組技術，核實本研究的推測，并進一步研究星形膠質細胞在亞洲人AD中的作用機制。

4 結 論

收集現有公開的亞洲人AD與正常人mRNA/miRNA表達譜數據，通過生物信息學手段，從轉錄組水平系統地對亞洲人群可能的AD發生機制進行了分析，最終篩選出差異miRNAs與核心差異基因，首次構建了亞洲人四個腦區的關鍵miRNA-mRNA差異表達網絡，發現了SLC14A1基因在四個腦區均差異表達。本研究推測SLC14A1和AQP1蛋白相互作用，影響了AD患者星形膠質細胞的體液轉運功能，最終影響神經元細胞，從而誘導AD的發病。此外，研究發現亞洲人AD海馬體與內嗅皮質中的特異表達靶基因，分別與神經元間信號傳遞與線粒體功能障礙相關。以上研究為亞洲人AD的診療提供新的潛在靶標，也為組學分析提供了新的思路。