基于心力衰竭HMOX2診斷標志物的機器學習模型

武淑琴，王遠涵，曹紅艷，張巖波,3

1.山西醫科大學基礎醫學院數學教研室，山西太原 030001；2.山西醫科大學公共衛生學院衛生統計教研室，山西太原 030001；3.重大疾病風險評估山西省重點實驗室，山西太原 030001

心力衰竭（heart failure，HF）是心臟病發展到嚴重階段導致心臟機械活動功能障礙的一種高發病率和高死亡率的臨床綜合征[1]。HF 通常是由其他心臟疾病引起，可以是原發性疾病，也可繼發于其他疾病。原發或繼發性的心臟疾病導致心臟功能減退，使心臟無法有效泵血，最終引發心力衰竭。近年來，免疫和炎癥被認為是影響HF 的常見病理生理學因素，基于生物信息學發現潛在的特征基因已成為一種更新穎、更可靠的疾病診斷方法。對疾病采用免疫治療已成為一種較新穎的臨床治療策略。有研究表明，SERCA2a 等基因已成為最有可能治療HF 的基因之一[2-4]。本研究選取基因表達綜合（gene expression omnibus，GEO）數據庫中的微陣列表達數據集，尋找HF 與正常樣本之間的差異表達基因（differentially expressed genes，DEGs），并通過生物富集分析探索潛在的生物學功能。同時采用最小絕對收縮和選擇算子（the least absolute shrinkage and selection operator，LASSO）、支持向量機遞歸特征消除（ support vector machine recursive feature elimination，SVM-RFE）和隨機森林（random forest，RF）三種高效機器學習算法對差異基因進行篩選，并進行特征基因外部數據集驗證[5]。使用R 軟件包“CIBERSORT”探討HF 組織免疫浸潤的結果，并綜合分析診斷基因與免疫細胞之間的相關性。

1 資料與方法

1.1 資料來源與處理

在GEO 數據庫中使用關鍵詞“HF、正常”或“HF、健康”，通過美國國家生物技術信息中心（National Center for Biotechnology Information，NCBI）數據庫平臺進行數據搜索。選擇的數據集類型是數組表達譜，生物類型是智人，數據集的樣本量大于20。使用R 包sva 消除來自不同平臺的數據批次效應。對該數據集的log2 轉換分位數歸一化信號強度進行校正，并輸出校正結果。篩選后獲取人類HF 樣本基因表達譜公共數據集GSE5406 和GSE84796 作為聯合芯片數據集，而GSE21610 作為獨立驗證數據集，下載各矩陣數據集文件及相應的平臺文件，見表1。

表1 訓練/驗證數據集的信息

1.2 方法

1.2.1 差異基因鑒定和基因功能富集分析 合并兩個數據集為訓練集，使用R 包“limma”對標準化的基因表達譜進行差異分析，調整后P<0.05 和logFC>0.5被確定為DEGs 的顯著性標準。使用ggplot2 繪制DEGs 火山圖，用pheatmap 繪制DEGs 熱圖。使用clusterProfiler 包對差異基因進行基因本體論（gene ontology，GO）、京都基因和基因組數據庫（Kyoto Encyclopedia of Genes and Genomes，KEGG）、疾病本體論（disease ontology，DO）和基因集富集分析（gene set enrichment analysis，GSEA）來探索和分析差異基因的生物學功能[11]。GO、KEGG 和DO 的截止標準設置為校正后P<0.05。GSEA 分析中，ggplot2 包用于繪圖，GSEA 的截止值設置為|NES|>1.0、校正后P<0.05。

1.2.2 機器學習算法模型構建 LASSO 回歸的目的是得到預測誤差最小的變量結果及其對應的回歸系數。通過約束回歸系數，得到最優結果[6]。用“glmnet”包獲取最佳值，具體計算方法：Cost()；其中，為懲罰系數，為向量。SVM-RFE 是基于SVM 的最大間隔原理的序列后向選擇其算法，通過模型訓練樣本，然后對每個特征得分進行排序，去掉最小得分的特征，用剩余特征繼續訓練模型，進行下次迭代，直到選出最佳的特征數，其算法用“e1071”包實現[7]。

1.2.3 特征基因的驗證 為進一步測試機器學習基因信號篩查效果，使用外部數據集（GSE21610）作為驗證數據集進行驗證。其受試者操作特征曲線（receiver operating characteristic curve，ROC 曲線）是使用R 包“pROC”生成，并使用曲線下面積（area under the curve，AUC）和95%CI來驗證特征基因的診斷效率。

1.2.4 免疫滲透評估 本研究探索免疫細胞在HF樣本中的浸潤，并分析診斷基因與免疫細胞之間的相關性，以尋找可能的病理生理過程。反卷積工具CIBERSORT 用于計算組織基因表達中的定量免疫細胞成分。比較HF 和正常人22 種免疫細胞的表達水平及免疫細胞之間的相關性。通過Spearman 秩相關分析得到診斷基因與免疫細胞亞型之間的關系，圖形用“ggplot2”包繪制。

1.3 統計學方法

數據預處理、模型的構建、驗證及評價均使用R 語言（版本4.2.2）及相關的R 程序包。所有統計學檢驗均為雙側檢驗，P<0.05 為差異有統計學意義。

2 結果

2.1 差異基因的識別

使用Affymetrix Human Genome U133 Plus 2.0 Array 分別分析GSE5406、GSE21610 數據集基因表達；使用Agilent-028004 SurePrint G3 Human GE 8x60K 微陣列分析GSE84796 數據集基因表達（表1）。根據篩選標準logFC>0.5 和調整后的P<0.05，獲得在正常樣本和HF 中差異表達的177個基因。采用火山圖描述數據集中所有差異基因的表達狀態，見圖1。上調和下調基因之間的差異明顯。

圖1 訓練數據集中所有差異基因的火山圖

圖3 三種機器學習算法示意圖

圖4 特征基因的相關信息

圖5 免疫細胞浸潤與通路活性分析

2.2 差異基因的功能富集分析

GO 分析3 個層面：生物過程（biological processes，BP）分析，差異基因主要參與細胞外基質組織、細胞外結構組織、外部封裝結構組織、細胞對細胞外刺激的反應、多機體繁殖等過程；細胞成分（cell component，CC）分析，差異基因主要富集在含膠原蛋白的細胞外基質、膠原蛋白三聚體、纖維狀膠原蛋白三聚體、帶狀膠原蛋白纖維和膠原蛋白三聚體復合體中；分子功能（molecular function，MF）分析，差異基因主要在細胞外基質結構成分、糖胺聚糖結合、肝素結合、膠原蛋白結合和細胞外基質結構成分賦予抗拉強度等方面發揮作用，見圖2A。在DO富集分析中，發現差異基因主要集中在一些心血管系統和炎癥性疾病中，包括心肌梗死、動脈硬化癥、肝炎、生殖細胞癌、腎功能不全、動脈硬化性心血管疾病等，見圖2B。KEGG 術語表明，差異基因主要參與PI3K-Akt 信號傳導、細胞因子–細胞因子受體相互作用、蛋白質的消化和吸收、cGMP-PKG 信號傳導等多條通路，從而在HF 的發生和發展中發揮重要作用，見圖2C。GSEA 富集分析發現，Th1-Th2細胞分化、MAPK 信號通路和B 細胞受體信號通路等與炎癥和免疫相關，且在HF 疾病和正常樣本中表現出顯著差異，見圖2D、2E。

2.3 機器學習確定HF 的診斷遺傳信號

LASSO 采用10 折交叉驗證確定16 個診斷基因信號，見圖3A；SVM-RFE 算法在特征數為6 時誤差最小，因此篩選出6 個最佳診斷信號，見圖3B；在RF 算法中，最終選擇26 棵樹作為參數，見圖3C，并把重要性評分>2 的5 個基因作為最佳診斷基因信號，見圖3D。整合三種機器學習算法確定HMOX2為HF 的診斷基因信號，見圖3E，其AUC 值為94.3%，見圖4A。

2.4 診斷基因信號的驗證

將基因信號HMOX2在數據集GSE21610中進行驗證，得出其在正常對照組中的表達高于HF 組，見圖4B，其AUC 值為83.3%，見圖4C。表明HMOX2具有較高的診斷效能，可作為潛在的生物標志物。

2.5 免疫細胞浸潤與通路活性分析

組織免疫細胞浸潤分析揭示227 個樣本的基因集合中的22 種免疫細胞亞型，見圖5A，發現3 種免疫細胞的差異：HF 疾病組CD8T 細胞（P<0.001）和中性粒細胞（P=0.008）顯著升高，而幼稚B 細胞（P=0.002）則為一組在HF 組織中低表達的細胞亞型，見圖5B。就免疫細胞亞型的相關性而言，巨噬細胞M0 與巨噬細胞M1 的正相關最顯著（r=0.35），活化的肥大細胞與靜息的肥大細胞負相關最顯著（r=–0.58），見圖5C。HF 特征基因與免疫細胞亞型的相關性分析顯示，HMOX2 與幼稚B 細胞呈最顯著正相關（r=0.24，P<0.001），與γ-δ 型T 細胞呈最顯著負相關（r=–0.15，P=0.027）。此外，HMOX2 還與漿細胞和巨噬細胞M1 相關，見圖5D。

3 討論

本研究表明慢性HF是一個全身性T細胞集落激活過程。在HF 進展中，T 細胞亞群（Th1-Th2）時空分布的變化對心室重構有較大的影響。Th2 細胞可能代表慢性HF 的治療靶點，有助于減少組織炎癥。此外，MAPK 信號通路通過血管緊張素Ⅱ與HF 發生關聯。通常認為抑制MAPK 信號通路可通過阻斷血管緊張素Ⅱ的分泌減輕HF 的發展。有研究發現B細胞通過依賴和獨立于抗體產生的機制在HF 的進展中發揮重要作用[9]。這也為B 細胞反應通路在HF中的作用提供新的見解。

HMOX2 又稱血紅素加氧酶-2，是血紅素加氧酶家族的組成性亞型成員，主要在大腦和睪丸中表達[10-12]。與長期以來一直是心血管研究重點的血紅素加氧酶-1（HMOX1）相比，對HMOX2 的研究仍處于起步階段。有報道稱，HMOX2 通過頸動脈中的BKCa2+通道在氧感知中發揮重要作用[13]。同時，HMOX2還通過調節細胞中的血紅素濃度及CO 和H2S 的水平影響多種生物過程。CO 可激活cGMP 信號通路、抑制強血管收縮劑內皮素-1 的產生，還可通過調節H2S 途徑調節腦血流量。這些生物過程大多是通過HMOX2 首先感知O2濃度而實現的[14]。據推測，HMOX2 能抑制缺氧疾病中的全身反應，但具體機制仍不清楚[15]。一般來說，氧化應激存在于大多數心血管疾病中[16]；其具體機制是大量的心肌細胞（心肌細胞、內皮細胞和中性粒細胞）可產生活性氧（reactive oxygen species，ROS）。在正常的生理條件下，心臟發揮防御性抗氧化功能，與ROS 的產生保持動態平衡。然而，在病理因素的刺激下，該平衡很快被改變，大量的ROS 被釋放，引起功能蛋白和脂質的過氧化及DNA 損傷，導致心肌收縮功能受損和細胞外基質重塑[17]。HMOX2 對清除細胞內的ROS 很重要，且在保護細胞免受ROS 引起的損傷方面也發揮重要作用[18]。本研究表明，HMOX2 的表達水平在HF 患者中降低，因此，根據HMOX2 活性表達情況，可對HF 診斷發揮一定的潛在作用。在HF治療方面，HMOX2 靶點的作用還在研究中，尚未得到充分的證據支持，需要更多的實驗和臨床研究確定HMOX2 靶點的潛在治療效果及其在HF管理中的具體作用機制。HMOX2 對HF 靶點基因的確定，可為HF 免疫治療提供新的思路。通過調控HMOX2 活性，可對HF 的分子治療具有潛在的作用。

本研究進一步探索HF 組織中免疫細胞的分布。其疾病組CD8T 細胞和中性粒細胞顯著升高，低表達的細胞亞型為幼稚B 細胞。研究發現HF 患者的CD8T 細胞、巨噬細胞、中性粒細胞的表達與正常人有顯著差異[19]。HF 后的組織損傷和纖維化伴隨復雜的免疫細胞反應，免疫細胞功能是HF 心肌損傷病理過程的中心環節，激活CD8T 細胞可抑制心室重塑。有研究發現，中性粒細胞在壓力超負荷引起的HF 實驗模型中發揮有害作用[20]。這一發現與中性粒細胞通過釋放脫顆粒和募集微泡參與多種心血管疾病發生、發展的基本作用機制相吻合。通過免疫細胞相關性棒棒糖圖探討HMOX2 與免疫細胞之間的相關性；免疫細胞相關性棒棒糖圖是一種熱圖的變種，其中每個細胞類型都用一列表示，每個基因用一行表示。通過觀察該圖，可判斷特征基因在免疫細胞中是否具有特異性表達，或某些免疫細胞類型中是否存在共同的表達模式。結果顯示：HMOX2 與幼稚B 細胞呈最顯著正相關，與γ-δ 型T 細胞呈最顯著負相關。此外，HMOX2 還與漿細胞和巨噬細胞M1 相關。這些結果表明，HMOX2 可通過調節免疫浸潤在HF 中發揮關鍵作用。

綜上所述，本研究通過LASSO 回歸、SVM-RFE和RF 三種機器學習算法篩選得到HMOX2診斷基因信號，可用作HF 診斷的特異生物標志物，為探索HF 的發病機制提供科學依據。