冠狀動脈粥樣硬化性心臟病患者心外膜脂肪組織的生物信息學研究

柴晏，趙玉青，郭旭男，王東英，邊云飛

動脈粥樣硬化（atherosclerosis，AS）導致的心血管疾?。╟ardiovascular disease，CVD）是全球死亡的主要原因，世界衛(wèi)生組織預測2030年將有2 360萬人死于CVD［1］。AS形成機制復雜，主要涉及血管內(nèi)皮損傷、單核巨噬細胞黏附、脂質(zhì)沉積等［2］，但這些因素如何影響AS進程，其機制尚未完全清楚。因此找尋AS新的診斷靶點，預防AS的發(fā)生，延緩其進展尤為重要。利用基因芯片技術(shù)對臨床患者標本進行檢測分析，篩選出一些具有價值的基因，進而深入研究冠狀動脈粥樣硬化性心臟?。╟oronary artery disease，CAD）的發(fā)病機制具有重要的臨床意義。肥胖是CVD的獨立危險因素，但肥胖的主要評價指標體質(zhì)指數(shù)（BMI）與AS患者死亡率呈“U”型關(guān)聯(lián)，這種現(xiàn)象被稱為肥胖悖論，因此脂肪組織不再被認為是單純的“能量倉庫”，而是一種代謝活躍的內(nèi)分泌和旁分泌器官［3］。心外膜脂肪組織（epicardial adipose tissue，EAT）在解剖上與冠狀動脈緊密相連，同AS、心房顫動、心力衰竭等CVD密切相關(guān)［4-6］。多項研究發(fā)現(xiàn)功能失調(diào)的EAT通過分泌外泌體（exosome，EXO）和生物活性物質(zhì)促進AS疾病進展［7］，但其作用機制仍需進一步研究。本研究通過對美國國家生物技術(shù)信息中心（NCBI）的基因表達綜合數(shù)據(jù)庫（Gene Expression Omnibus database，GEO）中 CAD患者EAT和EXO的基因測序數(shù)據(jù)進行生物信息學分析，探討EAT參與AS的分子機制及其與EXO的聯(lián)系，旨在為CAD的診斷和治療提供新的理論基礎(chǔ)。

1 材料與方法

1.1 數(shù)據(jù)集的下載與預處理 以“epicardial adipose tissue OR EAT”為關(guān)鍵詞，選擇“Homo sapiens”為研究對象檢索GEO數(shù)據(jù)庫，根據(jù)數(shù)據(jù)集提供的相關(guān)信息，獲得所需芯片GSE64554、GSE120774，根據(jù)臨床信息將EAT的測序數(shù)據(jù)分為CAD組和健康對照組。下載series matrix的數(shù)據(jù)，采用“sva”包“ComBat”方法去除批次效應(yīng)，合并基因表達矩陣。通過exoRBase 2.0數(shù)據(jù)庫獲取CAD組與健康對照組血液EXO基因表達譜。

1.2 差異基因的篩選 通過“l(fā)imma”R語言包［8］對合并后數(shù)據(jù)集中CAD組與健康對照組EAT間差異表達基因（differentially expressed genes，DEGs）和CAD組與健康對照組EXO間DEGs進行篩選，以P＜0.05為篩選條件，并將結(jié)果可視化為火山圖和熱圖。

1.3 加權(quán)基因共表達網(wǎng)絡(luò)（WGCNA）的構(gòu)建與CAD相關(guān)模塊的識別 為了防止合并數(shù)據(jù)集對基因關(guān)聯(lián)度的影響，使用“WGCNA”R語言包構(gòu)建GSE64554數(shù)據(jù)集基因共表達網(wǎng)絡(luò)并識別出模塊。選擇滿足無尺度網(wǎng)絡(luò)的標準軟閾值β構(gòu)建基因共表達網(wǎng)絡(luò)，轉(zhuǎn)化拓撲重疊矩陣（topological overlap matrix，TOM），采用動態(tài)剪枝法對模塊進行初步劃分。根據(jù)GSE64554數(shù)據(jù)集樣本臨床信息，決定納入的臨床信息包括疾病和年齡，得到相關(guān)性最高的模塊進行后續(xù)分析。計算基因與模塊的相關(guān)系數(shù)（module membership，MM）來評估基因與模塊之間的相關(guān)性，以|MM|＞0.8篩選模塊內(nèi)樞紐基因（hub gene）。

1.4 基因功能的富集分析 通過“clusterProfiler”R語言包［9］將所獲基因進行富集分析，以FDR＜0.05為篩選標準，進行通路富集分析。通過GO富集分析，得到所選基因所富集的生物過程（biological process，BP）、細胞組分（cellular components，CC）及分子功能（molecular function，MF）。通過KEGG通路富集分析，得到所選基因所富集的信號通路。將所獲基因?qū)隡etascape在線數(shù)據(jù)庫［10］進行基因的富集分析以及轉(zhuǎn)錄因子預測。

1.5 蛋白質(zhì)-蛋白質(zhì)相互作用（protein-protein interaction，PPI）網(wǎng)絡(luò)的構(gòu)建與hub 基因篩選 將所選基因?qū)隨TRING數(shù)據(jù)庫（https://string-db.org）構(gòu)建PPI網(wǎng)絡(luò)，將PPI網(wǎng)絡(luò)導入Cytoscape（Version 3.9.0）軟件［11］通過CytoHubba插件MCC算法獲得連接度最高的10個基因作為關(guān)鍵基因，并進行繪圖。

1.6 免疫細胞浸潤分析 在Cibersort網(wǎng)站（https://cibersort.stanford.edu/）下載22種免疫細胞的基因表達矩陣，采用Cibersort反卷積算法，計算GSE64554芯片中46個樣本的免疫細胞浸潤情況，并作圖對組織的免疫細胞浸潤情況進行比較。

1.7 臨床樣本的采集及實時熒光定量PCR（quantitative real-time PCR，qRT-PCR） 選取山西醫(yī)科大學第二醫(yī)院2021年6—10月行冠狀動脈造影患者20例，其中CAD患者及健康者各10例，留取外周血儲存在-80 ℃冰箱中備用。引物由武漢賽維爾生物科技有限公司合成，見表1。按照總RNA快速提取試劑盒（RP4001，百泰克公司）的說明提取血漿中的總RNA，cDNA的合成按照HiScriptⅢ RT SuperMix for qPCR（+gDNA wiper）試劑盒的步驟進行，第一步去除總RNA中殘余的DNA，第二步反轉(zhuǎn)錄cDNA，得到的cDNA保存于-20 ℃冰箱。以cDNA為模板，進行qRT-PCR，反應(yīng)程序為：95 ℃3 min；95 ℃ 5 s；60 ℃ 10 s，72 ℃ 30 s，40 個循環(huán)；95 ℃15 s，65 ℃ 5 s，95 ℃ 50 s。結(jié)果用2-ΔΔCT進行計算分析。

表1 PCR引物Table 1 Primer sequences for qRT-PCR

1.8 統(tǒng)計學方法 采用Excel對錄入整理數(shù)據(jù)，應(yīng)用SPSS 26.0軟件對數(shù)據(jù)進行統(tǒng)計學分析，符合正態(tài)分布的連續(xù)變量采用（±s）表示，兩組間比較采用獨立樣本t檢驗，以P＜0.05為差異有統(tǒng)計學意義。

2 結(jié)果

2.1 CAD組與健康對照組EAT間的DEGs

2.1.1 DEGs的獲取 采用“l(fā)imma”包獲取CAD組與健康對照組EAT間DEGs，以P＜0.05作為篩選差異基因的閾值，得到差異表達的1 511個mRNA，其中表達上調(diào)956個，表達下調(diào)555個，見圖1、2。

圖1 CAD組與健康對照組EAT間DEGs火山圖Figure 1 Volcano plot of differentially expressed genes in epicardial adipose tissue between CAD group and healthy control group

圖2 CAD組與健康對照組EAT間DEGs熱圖Figure 2 Heatmap of differentially expressed genes in epicardial adipose tissue between CAD group and healthy control group

2.1.2 EAT間DEGs的GO/KEGG富集分析與PPI網(wǎng)絡(luò)對合并后CAD組較健康對照組EAT間DEGs進行GO富集分析，共獲得符合篩選條件的GO術(shù)語共800條，其中BP 567條，主要涉及細胞壓力反應(yīng)、細胞DNA損傷刺激反應(yīng)等。獲得CC 127條，主要涉及細胞質(zhì)基質(zhì)、MHC蛋白復合物、內(nèi)質(zhì)網(wǎng)腔側(cè)的組成部分等。獲得MF 106條，主要涉及MHCⅡ類蛋白復合物、特殊核糖核酸酶激活等，見圖3。共富集到KEGG通路20個，主要涉及RNA降解、病毒性心肌炎、移植物抗宿主病、1型糖尿病等信號通路，見圖4。

圖3 EAT間DEGs的GO富集分析Figure 3 GO enrichment analysis of differentially expressed genes in epicardial adipose tissue

圖4 EAT間DEGs的KEGG富集分析Figure 4 KEGG enrichment analysis of differentially expressed genes in epicardial adipose tissue

將所獲EAT間DEGs輸入STRING數(shù)據(jù)庫，構(gòu)建PPI網(wǎng)絡(luò)，其中實際上有互作關(guān)系的節(jié)點有1 507個，有1 258條邊，每個節(jié)點的平均Degree為1.67分。將PPI網(wǎng)絡(luò)圖導入Cytoscape軟件，使用CytoHubba插件根據(jù)MCC算法得到連接度前10的基因并作圖，分別為 RPS27A、PSME4、PSMA1、PSMA5、GLI3、GLI1、PSME2、PSMB9、SUFU和SHH，見圖5。

圖5 EAT間DEGs的PPI網(wǎng)絡(luò)與關(guān)鍵基因Figure 5 Differentially expressed genes in epicardial adipose tissue and the key genes identified in the PPI network

2.1.3 EAT間DEGs的Metascape富集分析 對合并后CAD組與健康對照組EAT間DEGs進行Metascape富集分析，發(fā)現(xiàn)DEGs主要富集于嗅覺轉(zhuǎn)導、細胞對DNA損傷刺激反應(yīng)、單純皰疹病毒1型感染、RNA代謝、調(diào)節(jié)細胞對壓力的反應(yīng)、適應(yīng)性免疫系統(tǒng)等，TRRUST數(shù)據(jù)庫預測MHCⅡ反式激活蛋白（CIITA）轉(zhuǎn)錄因子可能參與了EAT間DEGs的調(diào)控，見圖6、7。

圖6 EAT間DEGs的Metascape富集分析Figure 6 Metascape enrichment analysis of differentially expressed genes in epicardial adipose tissue

圖7 EAT間DEGs的TRRUST數(shù)據(jù)庫預測轉(zhuǎn)錄因子Figure 7 TRRUST database-based analysis of differentially expressed genes in epicardial adipose tissue to predict transcription factors

2.2 EAT組織的WGCNA分析

2.2.1 EAT組織WGCNA軟閾值的選擇與模塊的初步劃分 本研究最終以β=5作為軟閾值來構(gòu)建共表達網(wǎng)絡(luò)，設(shè)置每個基因模塊中基因數(shù)目最小為30，設(shè)定基因的聚類高度上限為0.995。最終得到156個網(wǎng)絡(luò)模塊，模塊中的基因個數(shù)為30～2 169。根據(jù)拓撲重疊程度，制作模塊相關(guān)性圖，見圖8、9。

圖8 軟閾值參數(shù)網(wǎng)絡(luò)拓撲結(jié)構(gòu)的分析Figure 8 Network topology for soft thresholding powers

圖9 基因聚類樹狀圖Figure 9 Dendrogram of gene clustering

2.2.2 臨床信息關(guān)聯(lián)與hub基因的獲取 將各個劃分出的模塊與臨床信息進行關(guān)聯(lián)分析，得到與AS患者的EAT樣本臨床信息最為相關(guān)的模塊進行后續(xù)研究。結(jié)果發(fā)現(xiàn)，AS患者EAT與green4（r=0.50，P=0.02）、slateblue（r=0.54，P=8.0×10-3）、darkseagreen4（r=0.55，P=6.4×10-3）、darkolivegreen1（r=0.55，P=7.0×10-3）、sienna4（r=0.56，P=5.6×10-3）、brown3（r=0.57，P=4.2×10-3）模塊呈正相關(guān)，與 deeppink（r=0.50，P=0.01）、lavender（r=0.60，P=2.6×10-3）、lavenderblush3（r=0.56，P=5.9×10-3）呈負相關(guān)，以|MM|＞0.8篩選所選模塊內(nèi)hub基因，見圖10。

圖10 基因模塊與臨床信息關(guān)聯(lián)性Figure 10 Correlation between gene modules and clinical information

2.3 關(guān)鍵基因的獲得 將DEGs同模塊內(nèi)hub基因取交集，獲得差異表達的模塊內(nèi)hub基因，將其作為EAT組織參與AS病理生理過程的關(guān)鍵基因，其中上調(diào)hub基因為DDX47、FEM1C、NOL11、SRP54、ABI1、PATL1、BNIP2，下調(diào)hub基因為C1orf159、CHCHD4，見圖11。

圖11 差異性表達基因與模塊內(nèi)hub基因Venn圖Figure 11 Venn diagram of differentially expressed genes and hub genes in the module

2.4 免疫細胞浸潤分析 通過Cibersort反卷積算法，計算GSE64554芯片中CAD組和健康對照組EAT的免疫細胞浸潤情況。研究發(fā)現(xiàn)，CAD組EAT組織較健康對照組中初始CD4+T細胞表達豐度升高（P=0.048），靜息樹突細胞表達豐度減低（P=0.044），見圖12、13。

圖12 GSE64554數(shù)據(jù)集EAT組織免疫細胞浸潤柱狀圖Figure 12 Bar graph of immune cell infiltration in epicardial adipose tissue between CAD patients and controls included in dataset GSE64554

圖13 GSE64554數(shù)據(jù)集EAT組織免疫細胞浸潤小提琴圖Figure 13 Violin plot of immune cell infiltration in epicardial adipose tissue between CAD patients and controls included in dataset GSE64554

2.5 CAD組與健康對照組EXO間DEGs 采用“l(fā)imma”包選獲取CAD組與健康對照組EXO間DEGs，以P＜0.05，|logFoldchange|＞1作為篩選差異基因的閾值，共得到DEGs 1 658個，其中上調(diào)278個，下調(diào)1 380個，見圖14、15。

圖14 CAD組與健康對照組EXO間火山圖Figure 14 Volcano plot of differentially expressed genes in exosomes between CAD patients and heathy controls

圖15 CAD組與健康對照組EXO間DEGs熱圖Figure 15 Heatmap of differentially expressed genes in exosomes between CAD patients and heathy controls

2.6 EAT來源EXO 將CAD組與健康對照組EAT、EXO間DEGs取交集并做Venn圖，共獲得129個基因，其中上調(diào)16個，下調(diào)113個，經(jīng)篩選，將BPI、BIRC5、CXCL12、RNASE1、F2R作為關(guān)鍵基因并通過qRT-PCR進行驗證，見圖16。

圖16 CAD患者EAT組織、EXO間DEGs Venn圖Figure 16 Venn diagram of differentially expressed genes in epicardial adipose tissue and exosomes from CAD patients

2.7 EAT來源EXO間DEGs的驗證 收集CAD組和健康對照組外周血，運用qRT-PCR驗證生物信息學分析所得的EAT來源EXO基因BPI、BIRC5、CXCL12、RNASE1、F2R的mRNA水平。結(jié)果顯示，與健康對照組相比，CAD患者的5個所選基因的表達水平差異有統(tǒng)計學意義（P＜0.05），其中BPI、BIRC5、CXCL12和RNASE1的mRNA水平升高，F(xiàn)2R基因的mRNA水平下降（P＜0.05），見表2。

表2 CAD組與健康對照組EXO差異性基因表達（±s）Table 2 Differentially expressed genes in exosomes between CAD patients and healthy controls

表2 CAD組與健康對照組EXO差異性基因表達（±s）Table 2 Differentially expressed genes in exosomes between CAD patients and healthy controls

注：CAD=冠狀動脈粥樣硬化性心臟??；BPI、BIRC5、CXCL12、RNASE1、F2R為心外膜、EXO間DEGs

組別 樣本數(shù) BPI BIRC5 CXCL12 RNASE1 F2R CAD 組 10 2.14±0.32 1.98±0.23 2.91±0.89 4.47±1.56 0.45±0.06健康對照組 10 0.99±0.17 1.02±0.11 1.14±0.44 1.28±0.56 1.00±0.12 t值 -9.945 -12.070 -5.630 -6.001 12.637 P值 0.011 0.010 ＜0.001 0.009 ＜0.001

3 討論

EAT作為一種特殊的脂肪組織，解剖結(jié)構(gòu)上緊貼冠狀動脈，在病理環(huán)境中其向炎癥表型轉(zhuǎn)化，通過分泌促炎細胞因子和EXO參與血管內(nèi)皮損傷、血管重塑、免疫細胞黏附等AS病理生理過程。EAT的炎性反應(yīng)是AS的全新診斷治療靶點［5］，OIKONOMOU等［12］通過影像組學技術(shù)將心外膜脂肪衰減系數(shù)（fat attenuation index，F(xiàn)AI）應(yīng)用于CAD患者的早期篩查，而目前已廣泛應(yīng)用于臨床的鈉-葡萄糖協(xié)同轉(zhuǎn)運蛋白2（sodiumdependent glucose transporters 2，SGLT-2）抑制劑達格列凈則可以通過減輕EAT炎癥而減少糖尿病患者心血管事件的發(fā)生率［13］，但EAT在CAD中的作用機制仍需進一步研究。

此次通過生物信息學分析的方法從數(shù)據(jù)集GSE64554和GSE120774中篩選出CAD組和健康對照組EAT間DEGs 1 511個，通過富集分析發(fā)現(xiàn)DEGs主要富集于細胞質(zhì)基質(zhì)、MHC蛋白復合物、MHCⅡ蛋白復合物、RNA降解、抗原加工和呈遞等。目前AS被認為是一種慢性炎癥性疾病，MCHⅡ作為T細胞的呈遞抗原，通過抗原識別、細胞活化、產(chǎn)生細胞因子等方式在炎癥過程和CAD進展中發(fā)揮作用［14］。通過Meatscape富集分析中的TRRUST數(shù)據(jù)庫，筆者預測到CIITA轉(zhuǎn)錄因子可能在EAT中發(fā)揮了重要的調(diào)節(jié)作用。研究發(fā)現(xiàn)，CIITA通過其轉(zhuǎn)錄后修飾甲基化［15］、去乙?；?6-17］調(diào)節(jié)MHCⅡ的轉(zhuǎn)錄而參與AS。通過免疫細胞浸潤分析本研究發(fā)現(xiàn)CAD患者EAT中幼稚CD4+T細胞豐度升高，早在2015年IBORRA等［18］便指出可以通過觀察幼稚CD4+T細胞的體外分化程度評估AS程度，GADDIS等［19］指出AS通過影響糖酵解而損害幼稚CD4+T細胞功能。綜上所述，MHCⅡ在AS中發(fā)揮了重要作用，但其在CAD患者EAT中的作用和調(diào)控機制仍需進一步研究。

通過篩選CAD患者與健康對照組EAT間DEGs，構(gòu)建EAT組織WGCNA，獲得了DDX47、FEM1C、NOL11、SRP54、ABI1、PATL1、BNIP2、C1orf159、CHCHD4、RPS27A共10個基因作為EAT參與CAD病理生理過程的關(guān)鍵基因。目前僅有關(guān)鍵基因被發(fā)現(xiàn)在CVD中發(fā)揮作用，XIA等［20］通過生物信息學技術(shù)發(fā)現(xiàn)免疫相關(guān)基因RPS27A可能參與AS的病理過程，HUANG等［21］指出RPS27A基因在糖尿病視網(wǎng)膜病變的關(guān)鍵基因且同幼稚CD4+T細胞密切相關(guān)。LI等［22］發(fā)現(xiàn)MBNL1通過對ABI1的剪切調(diào)控AS中平滑肌細胞的表型轉(zhuǎn)換。因此，關(guān)鍵基因仍需進一步進行細胞和動物實驗證明其在EAT中的差異表達以及相關(guān)功能。

內(nèi)皮細胞、平滑肌細胞、心肌細胞、脂肪細胞和血小板來源的EXO通過轉(zhuǎn)運蛋白質(zhì)、RNA、DNA和其他生物活性物質(zhì)在缺血再灌注損傷、血管鈣化、AS和心臟重塑中發(fā)揮了重要作用，已成為CAD的全新診斷治療靶點［23］。脂肪組織是循環(huán)中EXO的來源之一［24］，冠狀動脈血管周圍脂肪組織是一種特殊的EAT［25］，LIU等［26］發(fā)現(xiàn)血管周圍脂肪組織來源EXO中miR-382-5p通過ABCA1/ABCG1信號通路調(diào)節(jié)泡沫巨噬細胞形成。LI等［27］發(fā)現(xiàn)血管周圍脂肪組織來源EXO中miR-221-3p通過介導血管重塑參與AS。ZHAO等［7］發(fā)現(xiàn)芒果苷誘導血管周圍脂肪組織衍生的EXO可以改善AS中的內(nèi)皮功能障礙，因此EAT來源EXO在AS中的作用具有重要的研究價值。

為進一步探討EAT在CAD中的作用機制，本研究將CAD患者與健康對照組的EAT間DEGs、EXO間DEGs取交集，共獲得129個基因，并篩選表達程度較高、有意義的基因RNASE1、BPI、BIRC5、CXCL12、F2R。其中前4個基因為CAD患者中上調(diào)基因。RNASE1是一種可分泌的耐熱蛋白，具有強大的心臟保護功能［28］，其通過保護內(nèi)皮細胞在炎癥中的損傷而成為血管穩(wěn)態(tài)的關(guān)鍵調(diào)控者［29］，在CAD中，EAT可能通過EXO釋放RNASE1進而調(diào)節(jié)冠狀動脈內(nèi)皮細胞的炎性損傷。血管生成趨化因子介導細胞趨化、白細胞脫粒和血管生成，CXCL-12同AS密切相關(guān)［30］，其通過CXCL-12-CXCR4軸參與平滑肌細胞表型轉(zhuǎn)化保護血管穩(wěn)態(tài)［31］，特異性過表達巨噬細胞IGF-1通過減少CXCL-12介導的單核細胞募集和增加ABCA1依賴性巨噬細胞脂質(zhì)外流增加AS斑塊的穩(wěn)定性［32］，但目前尚無關(guān)于EXO內(nèi)CXCL-12的相關(guān)研究。BIRC5是凋亡抑制蛋白（IAP）家族的成員，LI等［33］通過對小鼠頸動脈單細胞測序發(fā)現(xiàn)BIRC5可能與血流紊亂所致血管巨噬細胞聚集相關(guān)。早在2002年便有報道指出BPI的單核苷酸多態(tài)性同心肌梗死密切相關(guān)［34］，BLEIJERVELD等［35］通過循環(huán)中粒細胞的深度蛋白質(zhì)組分析發(fā)現(xiàn)了BPI可以作為嚴重AS性冠狀動脈狹窄的生物標志物，這表明BPI可以作為晚期CAD乃至預測急性心血管事件的標志物。F2R是血管壁炎癥和血栓形成的關(guān)鍵調(diào)節(jié)劑，其多態(tài)性同氯吡格雷療效密切相關(guān)［36］，但目前少有F2R在AS中的研究。

綜上所述，本研究共確定了DDX47、FEM1C、NOL11、SRP54、ABI1、PATL1、BNIP2、C1orf159、CHCHD4、RPS27A共10個關(guān)鍵基因，RNASE1、BPI、BIRC5、CXCL12、F2R共5個可能源自EAT EXO的基因進行了驗證，證明其在CAD中的診療意義，此外還探討了CAD患者EAT免疫浸潤情況，為CAD的基礎(chǔ)研究提供理論支持。本研究有其局限性，本研究主要針對mRNA，但蛋白質(zhì)才是生命活動的執(zhí)行者，mRNA層面的改變并不能等同蛋白質(zhì)層面的改變，因此本研究所得結(jié)論仍需動物、細胞實驗進一步驗證。

作者貢獻：柴晏、趙玉青、郭旭男負責課題設(shè)計、文章撰寫；柴晏、郭旭男負責數(shù)據(jù)庫檢索及數(shù)據(jù)整理；趙玉青、王東英負責實驗實施及實驗數(shù)據(jù)處理；邊云飛負責文章的質(zhì)量控制及審校。

本文無利益沖突。