心肌梗死潛在生物標志物的生物信息學研究

錢 晗, 曹夢菲, 呂書梅, 袁 偉

(江蘇大學附屬醫院 心血管內科, 江蘇 鎮江, 212001)

在人口老齡化和代謝危險因素的雙重壓力下，中國居民的心血管疾病風險持續增加，心肌梗死患者的病死率亦逐漸上升[1]。心肌梗死是指心肌缺血性壞死，通常是在冠狀動脈病變的基礎上發生冠狀動脈血供減少或中斷，心肌細胞因嚴重且持續的缺血缺氧而發生壞死[2]。心肌細胞的死亡會引發一系列免疫炎癥反應，這些免疫炎癥反應介導的受損心肌修復在梗死后心室重構、心功能和患者結局等方面起著關鍵作用[3]。急性心肌梗死的診斷主要依據臨床表現、血清生物學標志物和心電圖，其中臨床常用的生物學標志物包括肌紅蛋白、肌酸激酶-MB、心肌肌鈣蛋白I和心肌肌鈣蛋白T[4-5]。然而，這些標志物水平在心力衰竭、腎衰竭和甲狀腺疾病患者中也同樣升高，故結合多組學平臺探尋更多特異性的生物學標志物非常必要[6-7]。外周血單個核細胞(PBMCs)包括T細胞、B細胞、單核細胞、樹突狀細胞和自然殺傷細胞，在免疫應答、免疫監測和免疫治療中發揮著重要作用[8]。本研究基于PBMCs分析心肌梗死和梗死恢復過程中潛在的關鍵基因和生物學改變，旨在為心肌梗死的診斷、評估、管理和預后分析提供新的見解。

1 材料和方法

1.1 數據集獲取

從GEO數據庫中下載關于心肌梗死PBMCs的轉錄組測序芯片數據集GSE59867[9-10], 芯片注釋平臺為GPL6244, 芯片信息為Affymetrix Human Gene 1.0 ST Array。該數據集收集了111例ST段抬高型心肌梗死(STEMI)患者入院、出院、心肌梗死后1個月、心肌梗死后6個月時的外周血樣本數據。對照組為46例無心肌梗死病史的穩定型冠心病患者。基于芯片注釋平臺對數據集進行基因名轉換，用于后續共表達網絡的構建。

1.2 加權基因共表達網絡分析(WGCNA)

通過WGCNA R包篩選與心肌梗死有關的易感基因模塊[11]。WGCNA是一種常用的生物信息學方法，可以識別具有相似表達模式的基因模塊，分析基因模塊與樣本表型之間的關系，繪制基因模塊中的調控網絡，識別關鍵模塊和調控基因[12-13]。首先，選取數據集中方差變異程度最高的前3 000個基因構建共表達網絡; 其次，利用R函數pickSoftThreshold計算軟閾值1～20; 然后，對這些基因進行層次聚類和動態樹切割來劃分模塊，并合并相似的模塊(閾值為0.25); 最后，計算模塊與疾病之間的Spearman相關系數，進而確定易感基因模塊。為了探討所選模塊的潛在生物學功能，通過ClusterProfiler R包對模塊中的基因進行基因本體論(GO)分析和京都基因與基因組百科全書(KEGG)通路分析[14]。GO分析包括生物過程(BP)、分子功能(MF)和細胞成分(CC)。通過Ggplot2 R包進行結果展示。

1.3 基因差異分析和轉錄因子預測

使用Limma R包進行基因差異分析[15]。差異表達基因(DEGs)的篩選標準為矯正后P<0.05, |log2(倍數變化)|>0.5。基于ClusterProfiler R包對DEGs進行GO分析和KEGG通路分析,P<0.05表示顯著富集。對上述易感基因模塊和DEGs進行交集分析，通過ChEA3數據庫對交集基因進行轉錄因子預測。ChEA3是一種轉錄因子富集分析工具，整合了多個數據庫，如ENCODE、GTEx、ARCHS4等數據庫[16]。由于ChEA3數據庫并未對預測結果可靠度進行嚴格劃分，參考相關研究[16-17]結論，本研究對平均排名分≤30分的共同轉錄因子進行結果展示。基于STRING數據庫構建交集基因的蛋白相互作用網絡，并將交互評分>0.4分的數據導入Cytoscape 3.8.2軟件進行可視化處理[18-19]。

1.4 關鍵基因的確定和受試者工作特征(ROC)曲線分析

LASSO回歸分析是常用的篩選變量的壓縮估計方法，通過構造懲罰函數獲得較為精煉的模型，使得其壓縮一些回歸系數，同時設定一些回歸系數為0，因此保留了子集收縮的優點，是一種處理具有復共線性數據的有偏估計，常被用于協變量篩選[20]。基于glmnet R包對交集基因中存在蛋白相互作用的基因進行LASSO回歸分析，進一步篩選關鍵基因[21]。采用10倍交叉驗證法選擇懲罰項(λ), 并選擇最小二項偏差在1個標準誤差內的最簡單模型的λ值(lambda.1se)。通過pROC R包對關鍵基因進行ROC曲線分析，計算曲線下面積(AUC), 評估關鍵基因對心肌梗死的診斷價值[22]。

1.5 基因集富集分析(GSEA)和外部數據集驗證

根據每個關鍵基因的表達量將數據集劃分為高表達組和低表達組，通過GSEA探討高表達組和低表達組之間潛在的生物學功能變化，將特征基因集作為預先定義的基因集。作為一種常見的生物信息學方法, GSEA可評估一個預先定義的基因集的基因在與表型相關度排序的基因表中的分布趨勢，從而判斷其對表型的影響[23]。從GEO數據庫中選取另一個數據集GSE123342, 將其用于驗證關鍵基因的表達。GSE123342數據集包括急性心肌梗死(65例)、梗死后30 d(64例)、梗死后1年(37例)、穩定型冠心病(22例)樣本數據以及4例技術重復樣本數據。

2 結 果

2.1 易感基因模塊的確定

根據無標度網絡原理，選取10作為共表達網絡的軟閾值，見圖1A。通過層次聚類和動態樹切割，將前3 000個基因劃分為15個基因模塊，見圖1B。剔除聚類失敗的灰色模塊，棕色模塊的總體表達與心肌梗死的相關性最高(r=0.4), 且與出院至心肌梗死后6個月的演變過程呈負相關(r=-0.33), 見圖1C。棕色模塊中的共表達基因伴隨著心肌梗死的發生呈現高表達趨勢，而在心肌梗死恢復過程中則呈現低表達趨勢。由此提示，棕色模塊可能是易感基因模塊，其中包含243個基因，即ABHD5、ACER3、ACP3、ACSL1、ADAM17、ADAM9、ADM、AHR、ALCAM、ALDH2、ALPK1、ANKRD50、ANO5、AQP9、ARHGAP24、ARHGAP29、ARRDC4、ASPH、ATP6V0A1、ATP6V1A、B3GNT5、BACH1、BLVRB、BST1、C3AR1、C5、C9orf72、CALCRL、CAPG、CAPZA2、CARD6、CASP1、CASP5、CCDC88A、CCR2、CD14、CD163、CD1D、CD33、CD36、CD63、CLEC1A、CLEC4A、CLEC4D、CLEC4E、CLEC7A、CLMN、CNTLN、CPED1、CPM、CPNE8、CPVL、CR1、CR1L、CREG1、CRISPLD2、CSTA、CTNNA1、CTSH、CYBRD1、CYP1B1、CYP1B1-AS1、DACH1、DDIAS、DOCK4、DOCK5、DRAM1、DSC2、DSE、DUSP6、EDNRB、ENTPD1、ERLIN1、EVI5、F5、FAM114A1、FAM13A、FAM151B、FAM20A、FAM198B、FAR2、FBN2、FBP1、FCGR2A、FGD4、FGD6、FLT3、FLVCR2、FMO5、FPR1、FPR2、FRRS1、FUCA1、FUCA2、GAPT、GAS2L3、GCA、GIMAP8、GLT1D1、GPAT3、GPR141、HAL、HGF、HMGB2、HNMT、HORMAD1、HP、HPSE、IDH1、IFNGR1、IL13RA1、IL15、IL18、IMPA2、IRAK3、JAK2、KCNJ15、KCNJ2、KIF13A、KLHL8、KYNU、LHFPL2、LILRA5、LIN7A、LIPN、LMNB1、LPCAT2、LRMDA、LRRK2、LTA4H、LY86、LYVE1、MAP3K20、MCEMP1、ME1、METTL7A、METTL7B、MFSD1、MGST1、MGST2、MILR1、MNDA、MOSPD2、MS4A4A、MSRB1、MTARC1、MYOF、NAAA、NAIP、NETO2、NLN、NLRC4、NPL、OLFML2B、P2RX7、P2RY13、PGD、PIP4P2、PLA2G4A、PLBD1、PLD1、PLIN2、PLSCR1、PPARG、PPT1、PRRG4、PSTPIP2、PYGL、QPCT、RAB39A、RALB、RBP7、RGL1、RNASE2、RNASE6、RNF141、RNF217、RRAGD、S100A12、S100A9、S100Z、SAT2、SEMA3C、SEPTIN10、SERPINB1、SERPINB2、SESTD1、SGMS2、SH3PXD2B、SHTN1、SIPA1L2、SIRPD、SLC15A2、SLC1A3、SLC22A15、SLC22A4、SLC26A8、SLC31A1、SLC36A4、SLC7A7、SLC8A1、SLITRK4、SMAD1、SMPDL3A、SNX10、SOD2、SORT1、SPATA6、ST3GAL6、ST6GALNAC3、STEAP4、SULT1B1、TASL、TBC1D12、TCN2、TDRD9、TFEC、TGFBI、TLR1、TLR2、TLR4、TLR5、TLR6、TLR7、TLR8、TM6SF1、TMEM144、TMEM167A、TMTC2、TNFAIP6、TNFSF13B、TPST1、TSPO、UBE2D1、UGGT2、VNN2、VNN3、WASF1、WDFY3、WLS、ZC3H12C、ZFYVE16、ZNF438。

棕色模塊GO分析最顯著的結果分別是骨髓白細胞激活(BP)、囊泡(CC)和模式識別受體活性(MF), 見圖1D; KEGG分析結果則主要為感染、免疫炎癥等相關生物學過程的改變，包括Toll樣受體信號通路等，見圖1E。

2.2 基因差異分析

進一步對數據集進行基因差異分析，根據篩選標準共得到142個DEGs(圖2A), 其中上調的DEGs有77個，下調的DEGs有65個，見表1。BP分析結果主要包括免疫反應、免疫系統過程和免疫效應過程, CC分析結果包括分泌顆粒、分泌小泡和細胞質小泡部分等, MF分析結果主要包括信號受體活性、分子傳感器活性和碳水化合物結合，見圖2B。KEGG分析提示，這些DEGs主要涉及免疫、炎癥和感染的相關生物學改變，見圖2C。DEGs與上述棕色模塊存在35個交集基因，見圖2D。基于ChEA3數據庫對35個交集基因進行共同轉錄因子預測，其中平均排名分≤30分的共同轉錄因子分別是CREB5、MTF1、NFE4、SPI1、ZNF467、NFE2、TFEC、MXD1、NR1H3、BORCS8MEF2B、NFIL3和CEBPE, 見表2。構建蛋白互作網絡，其中有18個基因存在蛋白質相互作用，見圖2E。

表1 差異基因與交集基因

表2 共同轉錄因子預測結果(平均排名分≤30分)

2.3 關鍵基因的確定

對18個存在蛋白相互作用的基因進一步行LASSO回歸分析，這些基因不同懲罰參數值所對應的系數見圖3A。本研究選取最小二項偏差在1個標準誤差內的λ值(lambda.1se, 6個基因)，該λ值提供更精簡的模型，見圖3B。ROC曲線分析顯示,CD163、RNASE2、HP、FAM20A、MCEMP1和FAM198B基因表達水平對心肌梗死的發生均具有良好的診斷價值,AUC分別為0.831、0.798、0.775、0763、0.866和0.829, 見圖3C。與對照組穩定型冠心病患者相比，這些關鍵基因的表達水平在入院時心肌梗死患者中顯著上升，并在出院時、心肌梗死后1個月、心肌梗死后6個月逐步下降，見圖3D。

2.4 GSEA結果和驗證

GSEA結果提示，這些關鍵基因涉及的生物學改變主要與糖脂代謝、活性氧、免疫炎癥等有關，見圖4。這些基因在外部數據集中大多也存在差異表達。基于外部驗證數據集GSE123342, 本研究同樣發現，與穩定型冠心病相比,HP、FAM198B、CD163、FAM20A、MCEMP1表達水平在急性心肌梗死發生時顯著上升，而RNASE2僅表現出上升趨勢，見圖5。

3 討 論

本研究基于PBMCs轉錄組學的變化探討心肌梗死及梗死恢復過程中可能存在的生物學改變和潛在標志物, WGCNA和DEGs分析結果顯示，心肌梗死過程伴隨著免疫炎癥紊亂。基于LASSO分析，本研究從上述2種分析方法得到的共同基因中鑒定出6個關鍵基因，即CD163、RNASE2、HP、FAM20A、MCEMP1和FAM198B。ROC曲線分析結果提示，這些關鍵基因對心肌梗死的發生均具有較高的診斷價值。本研究還發現，在梗死恢復過程中，這些基因的表達呈現明顯下降趨勢。GSEA分析結果表明，這些關鍵基因涉及的生物學改變主要與糖脂代謝、活性氧、免疫炎癥等有關。

CD163分子是Ⅰ型膜蛋白，蛋白表達限于單核細胞/巨噬細胞系, CD163抗原特異性釋放機制可能在炎癥調節過程中起重要作用。在動脈粥樣硬化過程中, CD163+巨噬細胞能促進血管生成和增加血管通透性，并伴隨炎癥反應；此外，破裂的冠狀動脈斑塊中CD163的表達增加了心肌梗死和冠心病發生風險[24]。RNASE2編碼的蛋白質是非分泌型核糖核酸酶，屬于胰核糖核酸酶家族[25]。FAM20A是激酶編碼基因家族的成員，本身不具有激酶活性，通過與家族成員FAM20C形成復合物，增強FAM20C的激酶活性，從而使分泌通路內的蛋白磷酸化[26]。MCEMP1又稱C19ORF59, 可編碼表達于肥大細胞、巨噬細胞等的跨膜蛋白[27]。研究[27]表明,MCEMP1基因可能是卒中診斷和預后評估的新生物標志物。生物信息學研究[28-29]提示,RNASE2、FAM20A和MCEMP1可能是心肌梗死的關鍵基因，然而這些基因在心肌梗死及梗死恢復過程中的具體作用機制有待進一步研究。HP基因編碼一種前蛋白，經處理后產生α鏈和β鏈，結合為四聚體產生觸珠蛋白，其基因型與急性心肌梗死的發病風險密切相關，并且能決定心肌梗死面積[30-31]。FAM198B(GASK1B)是目前未知功能的新基因，可能編碼定位在高爾基體上的膜結合糖蛋白，被認為參與癌癥的轉移與進展[32]。本研究分析結果顯示,FAM198B可能是心肌梗死和反映心肌梗死恢復的關鍵基因。

本研究尚存在一定局限性: 首先，本研究基于GEO數據庫中下載的心肌梗死患者PBMCs的測序數據集進行分析，還需要進一步開展分子生物學實驗來驗證這些基因在心肌梗死中的表達及其可能的生物學機制; 其次，本研究僅分析了現有的樣本數據，未來還需基于更大的樣本規模和更詳細的樣本類型進一步深入研究; 最后，由于缺乏臨床資料，本研究未分析這些關鍵基因的表達與心肌梗死患者年齡、性別、射血分數等臨床指標的相關性。

綜上所述，本研究基于PBMCs確定了6個可能與心肌梗死密切相關的關鍵基因，即CD163、RNASE2、HP、FAM20A、MCEMP1和FAM198B, 為心肌梗死的診斷、評估、管理和預后分析提供了新的思路。