候圣祥 ，張敬芳，胡 俊，劉 璐，陳天禹，龐 燊

雄激素性脫發(fā)（AGA）是導(dǎo)致臨床成年人脫發(fā)的最主要的疾病。在我國，男性患病率約為21.3%，女性患病率約為6.0%[1]。近年來隨著對AGA的流行病學(xué)調(diào)查的深入研究，發(fā)現(xiàn)遺傳因素在發(fā)病的時間及程度過程中起著重要作用。盡管AGA不是威脅生命的疾病，但它會嚴重影響患者的精神狀態(tài)和生活質(zhì)量。目前，2%～5%米諾地爾溶液可通過局部應(yīng)用促進頭發(fā)生長[2]，而5α-R抑制劑非那雄胺和度他雄胺是經(jīng)美國聯(lián)邦食品藥品監(jiān)督管理局（FDA）批準可用于治療AGA的口服合成藥物[3]。但是這些藥物均具有不良反應(yīng)、功效有限，且必須持續(xù)用藥的缺點[4]。因此需要更多方法來發(fā)現(xiàn)潛在的關(guān)鍵（樞紐、Hub）基因，進而揭示其分子機制。加權(quán)基因共表達網(wǎng)絡(luò)分析（WGCNA）又稱基因權(quán)重共表達網(wǎng)絡(luò)分析，常用于分析樣品基因與表型相關(guān)聯(lián)系。因此，本研究基于基因表達數(shù)據(jù)庫（GEO）中AGA數(shù)據(jù)，利用WGCNA挖掘AGA相關(guān)的基因網(wǎng)絡(luò)，發(fā)現(xiàn)脫發(fā)中特異表達的關(guān)鍵基因，進一步分析和AGA相關(guān)的蛋白功能及通路，為AGA建立相關(guān)的基因功能和信號通路分析，以期獲得治療AGA的潛在靶標。

1 材料與方法

1.1 數(shù)據(jù)收集與數(shù)據(jù)預(yù)處理

本項目的AGA基因表達矩陣（GSE90594）來自GEO數(shù)據(jù)庫。包含14份AGA患者頭皮樣本，和14份正常人樣本，平臺為GPL17077。使用R語言程序包（V4.0.1）對芯片的數(shù)據(jù)加以分析，對樣本進行聚類分析后去除離群樣本，利用WGCNA程序包分析芯片數(shù)據(jù)，選取與AGA正、負相關(guān)性最高的模塊，使用clusterprofiler軟件包完成基因組百科全書(KEGG)及基因本體論(GO)、limma軟件包分析基因的表達。

1.2 加權(quán)共表達模塊的計算

利用WGCNA程序包分析芯片數(shù)據(jù)，選取β值為5（圖1）用來估算網(wǎng)絡(luò)拓撲重疊TOM，通過層次聚類分析鏈接值的差異獲得的基因聚類樹，基于動態(tài)分層剪切數(shù)法運算來獲得模塊。計算基因模塊的特征值（ME），聯(lián)系臨床信息，對ME進行分段聚類并重新排列樹狀圖，設(shè)置高度值0.7為分割線，合并相似程度較高的基因模塊，再用剪切后的模塊替換新的聚類樹和模塊圖。

圖1 WGCNA中的軟閾值篩選圖

1.3 共表達模塊與臨床表現(xiàn)相關(guān)性的分析

選取與AGA正、負相關(guān)性最高的模塊各一個，通過對模塊內(nèi)基因分析基因臨床特征相關(guān)性(GS)和基因模塊表達水平(MM)的計算和取絕對值,篩選出高度正相關(guān)或負相關(guān)的基因,從而識別關(guān)鍵樞紐基因。

1.4 構(gòu)建基因網(wǎng)絡(luò)與尋找Hub基因

本項目選取與臨床顯著相關(guān)的基因模塊，計算模塊內(nèi)基因的MM和GS值，將|MM|＞0.8且|GS|＞0.6的基因?qū)С?，它們對具有模塊主要的調(diào)控作用。用Cytoscape軟件繪制模塊中共表達基因網(wǎng)絡(luò)關(guān)系[5]通過基因之間加權(quán)的共表達關(guān)系使用Cytoscape繪制網(wǎng)絡(luò)圖，通過MCODE[6]插件尋找Hub基因。

1.5 GO和KEGG富集分析

選取與臨床表型顯著相關(guān)的基因模塊,使用clusterProfiler軟件包進行GO分析和KEGG富集分析[7]。

1.6 LASSO模型的構(gòu)建和ROC曲線分析

LASSO具有較強的預(yù)測值和相關(guān)性高,并且可以用于高維數(shù)據(jù)等功能[6,7]。為了區(qū)分AGA與對照組，通過glmnet包提取了基因的表達譜以構(gòu)建LASSO模型。結(jié)果使用來自LASSO分析的回歸系數(shù)為每個樣本創(chuàng)建模型索引，以下列公式加權(quán)所選基因的表達值：索引= ExpGene1×Coef1 + ExpGene2 × Coef2 +ExpGene3×Coef3?！癈oef”是基因的回歸系數(shù)，它是從LASSO Cox回歸得出的，“Exp”表示基因的表達值。為了評估LASSO模型識別AGA的能力，使用pROC軟件包[8]在自身進行ROC曲線分析。

2 結(jié)果

2.1 芯片數(shù)據(jù)處理及加權(quán)基因共表達模塊的構(gòu)建

GSE90594芯片數(shù)據(jù)及臨床信息的下載及預(yù)處理,過濾其中探針信息注釋不全和重復(fù)的基因，最終獲得29 499個基因?qū)?yīng)27個樣本的表達矩陣，繪制樣本的分層聚類圖與對應(yīng)的臨床信息的熱圖（圖2）。

圖2 樣本對應(yīng)性狀熱圖

通過WGCNA包的算法,根據(jù)無尺度網(wǎng)絡(luò)分布擬合，選取5作為軟閾值，并計算基因間的相關(guān)性矩陣和TOM，使用TOM構(gòu)建基因間分層聚類樹，同時使用動態(tài)剪切樹的方法把基因分成20個模塊（圖3）

圖3 模塊的動態(tài)切割圖

2.2 共表達模塊與臨床表型的相關(guān)性分析

通過計算各個模塊與臨床表型之間的關(guān)系，繪制基因模塊與臨床表型熱圖，從基因模塊與臨床表型的熱圖中可以看出brown模塊與turquoise模塊與臨床分型顯著相關(guān)。能體現(xiàn)出該評分高低在正常人和AGA患者間基因表達的差異有非常強的聯(lián)系。

2.3 共表達網(wǎng)絡(luò)的可視化

通 過 |MM|＞ 0.8且 |GS|＞ 0.6的標準篩選brown模塊和turquoise模塊，分別在brown模塊和turquoise模塊中篩選出73、138個符合標準的基因。分別將篩選出的基因上傳至STRING數(shù)據(jù)庫，并通過Cytoscape對共表達網(wǎng)絡(luò)基因間的相互作用關(guān)系的可視化，使用Cytoscape繪制網(wǎng)絡(luò)圖（圖4），通過MCODE篩選樞紐基因(hubgenes)。在brown模塊尋找出PDGFRA、PMP22、ZCCHC24、COL6A1、ISLR、PRRX1 和turquoise模 塊 尋找 出PKP1、CALN1、PNMAL2、PPP5D1、GJB6、DSC2、GJA3共計13個樞紐基因（圖5）。

圖5 brown和turquoise模塊基因共表達網(wǎng)絡(luò)圖

2.4 GO和KEGG富集分析

選取與臨床高度相關(guān)的brown模塊和turquoise模塊進行GO與KEGG富集分析（圖6）。GO分析發(fā)現(xiàn)brown模塊的基因主要參與淋巴細胞活化的調(diào)節(jié)、免疫效應(yīng)過程的調(diào)節(jié)、T淋巴細胞活化、細胞粘附的正向調(diào)節(jié)、細胞激活的正調(diào)節(jié)、B淋巴細胞活化、B淋巴細胞活化的調(diào)節(jié)等生物功能；turquoise模塊的基因主要參與表皮發(fā)育、皮膚發(fā)育、細胞器分裂、分化形成角質(zhì)形成細胞、染色體分離、角化、減數(shù)分裂細胞周期、翻譯起始等生物功能。KEGG分析發(fā)現(xiàn)brown模塊的基因參與內(nèi)吞作用、造血細胞系、弓形蟲病、T淋巴細胞受體信號通路、補體與凝血級聯(lián)、利什曼病、原發(fā)性免疫缺陷等重要通路；turquoise模塊的基因參與河馬（hippo）信號通路、Wnt信號通路、細胞周期、黑素生成、基底細胞癌等通路。

2.5 構(gòu)建LASSO模型

LASSO模型是AGA的潛在預(yù)測指標。本文提取了中樞基因的表達譜以構(gòu)建LASSO模型（圖7）。使用LASSO方法，鑒定了5個具有非零回歸系數(shù)的基因，并且lambda.min的值= 0.02903?；诨虻哪Ｐ退饕匆韵鹿絼?chuàng)建：索引=PRRX1 *(注:表示相乘)（1.9166078）+GJA3 * （-1.0668399）+DSC2*（-2.9185647）+COL6A1*（4.0211617）+CALN1*（0.3254326）。ROC曲線分析表明基于5基因的模型的訓(xùn)練集的AUC在訓(xùn)練集中為0.98，這表明LASSO模型可以用作AGA的生物標志物。

圖7 5基因lasso模型預(yù)測自身ROC曲線

3 討論

本項目對GSE90594數(shù)據(jù)集轉(zhuǎn)進行權(quán)重基因共表達網(wǎng)絡(luò)分析，成功構(gòu)建20個模塊，其中brown模塊和turquoise模塊與AGA密切相關(guān)的新基因共表達網(wǎng)絡(luò)模塊。筆者進一步探索與AGA相關(guān)的生物學(xué)過程。

功能富集分析表明與AGA高度相關(guān)的模塊參與了細胞分化、表皮形成等生物過程，并顯著參與WNT等通路。Wnt/β-catenin信號通路對毛囊干細胞的生長發(fā)育有著密切相關(guān)性，且在AGA中Wnt/β-catenin信號通路活性常常出現(xiàn)異常。Wnt/β-catenin信號通路在毛發(fā)周期期間對毛囊形態(tài)、發(fā)育和再生中起關(guān)鍵作用[8]。在毛發(fā)周期中，毛囊干細胞的端粒期-生長期轉(zhuǎn)變是由Wnt/β-catenin信號通路與BMP途徑之間的平衡相互作用驅(qū)動的[9]。在毛囊隆起和繼發(fā)性發(fā)芽中，表皮干細胞和黑素細胞干細胞共同經(jīng)歷Wnt/β-catenin信號通路激活參與生發(fā)[10]?？紤]到 Wnt/β-catenin 信號通路的改變會導(dǎo)致AGA中毛囊的狀態(tài)的改變，因此可以通過以下方式調(diào)節(jié)Wnt/β-catenin信號通路，促進頭發(fā)再生及發(fā)育[11]：①調(diào)節(jié)Wnt信號通路配體的分泌；②改變配體-受體的結(jié)合；③促進β-catenin向細胞核移位。對于AGA，由于β-catenin表達不穩(wěn)定，導(dǎo)致處于生長期的毛囊干細胞減少。通過上調(diào)Wnt/β-catenin信號通路中β-catenin蛋白的活性可以延長發(fā)囊生長期。AGA患者的頭皮在很大程度上保留了靜止的毛囊干細胞，但缺乏祖細胞[12]。顯然，毛囊干細胞向祖細胞轉(zhuǎn)化的減少歸因于雄激素介導(dǎo)的Wnt/β-catenin信號通路活性的異常[13]。調(diào)節(jié)雄激素受體（AR）反式激活并穩(wěn)定Wnt/β-catenin信號通路可以治愈，逆轉(zhuǎn)并預(yù)防AGA中的脫發(fā)，并將不良反應(yīng)降至最低。除新療法外，靶向激活Wnt/β-catenin信號通路途徑對頭發(fā)的生長也顯示出促進的效果。

綜上所述，PDGFRA、PMP22、ZCCHC24、COL6A1、ISLR、PRRX1的 表 達 與AGA呈 正 相 關(guān)，PKP1、CALN1、PNMAL2、PPP5D1、GJB6、DSC2、GJA3的表達與AGA呈負相關(guān)。本研究不足之處在于GSE90594選取的對象為發(fā)病樣本，缺少發(fā)病前期樣本，無法分析基因在疾病各個時間及空間的轉(zhuǎn)錄水平；且為頭皮樣本，無法進一步分析在不同種類細胞中轉(zhuǎn)錄組的表達。今后可做不同發(fā)病時間段的單細胞轉(zhuǎn)錄組數(shù)據(jù)處理，可以更加清晰的表現(xiàn)各個基因在不同時間段，不同種類細胞中表達的差異。