操利超，巴 穎，盧曉萍，張核子

（深圳市核子基因科技有限公司，廣東 深圳 518071）

結(jié)腸癌（colorectal cancer，CRC）是一種常見的惡性腫瘤，是世界上第二大致死原因[1]。盡管結(jié)腸癌的診斷和治療已經(jīng)取得了很大的進(jìn)展，但結(jié)腸癌患者通常會出現(xiàn)復(fù)發(fā)和轉(zhuǎn)移，導(dǎo)致5 年生存率顯著下降[2]。因此，迫切需要改善結(jié)腸癌患者的診斷、治療和預(yù)后。近些年來，分子診斷技術(shù)已廣泛應(yīng)用于腫瘤的治療、預(yù)后領(lǐng)域[3-5]。生物信息學(xué)和機(jī)器學(xué)習(xí)技術(shù)已廣泛應(yīng)用于腫瘤診斷或預(yù)后分子標(biāo)志物的識別，這種分子標(biāo)志物類型多種多樣，如microRNAs[6]、長鏈非編碼RNA[7]、差異表達(dá)基因[8]、DNA 甲基化[9]等。其中，差異表達(dá)基因作為潛在的腫瘤診斷或預(yù)后標(biāo)志物應(yīng)用最為廣泛。為得到廣泛驗證的結(jié)腸癌相關(guān)的差異表達(dá)基因，本文利用生物信息學(xué)方法，從多個數(shù)據(jù)集、不同的數(shù)據(jù)庫中尋找共同的結(jié)腸癌相關(guān)的差異表達(dá)基因，并進(jìn)一步利用機(jī)器學(xué)習(xí)的方法，從這些差異基因中挑選出結(jié)腸癌預(yù)后相關(guān)的預(yù)測因子，并建立預(yù)后風(fēng)險評估模型。

1 材料與方法

1.1 數(shù)據(jù)下載和獲取 通過GEO 數(shù)據(jù)庫（https://www.ncbi.nlm.nih.gov/geo/）下載基因芯片表達(dá)數(shù)據(jù)集GSE44076、GSE28000 和GSE39582，每個參考數(shù)據(jù)集的正常和腫瘤樣本情況見表1。TCGA 中mRNA表達(dá)數(shù)據(jù)集和對應(yīng)的臨床信息從UCSC Xena 平臺（https://xenabrowser.net/datapages/）下載，選擇隊列為GDC TCGA Colon Cancer（COAD），樣本信息見表2。

表1 3 個GEO 數(shù)據(jù)集的樣本量情況

表2 TCGA 數(shù)據(jù)集的樣本信息[n（%）]

1.2 差異基因分析和統(tǒng)計分析 利用R 包分別對3個GEO 數(shù)據(jù)集和TCGA 數(shù)據(jù)集進(jìn)行差異基因分析，過濾標(biāo)準(zhǔn)為adjustedP-value＜0.05 和差異倍數(shù)1.5倍（|log2FC|＞0.585），然后取交集，得到共同的上調(diào)差異基因和下調(diào)差異基因。

1.3 構(gòu)建和評估預(yù)后風(fēng)險評分模型 為了確定與生存相關(guān)的差異表達(dá)基因，使用R 包Survival 進(jìn)行單變量Cox 比例風(fēng)險回歸模型（P＜0.05）。接著，使用LASSO 回歸分析進(jìn)一步縮減預(yù)后因子數(shù)量，通過多因子回歸分析確定每個預(yù)后因子的回歸系數(shù)，建立預(yù)后風(fēng)險評估模型，預(yù)測患者生存率。公式為：

1.4 繪制生存曲線和ROC 曲線 根據(jù)風(fēng)險評分預(yù)后模型，計算每個腫瘤樣本的預(yù)后因子風(fēng)險評分。利用R 包survivalROC 繪制ROC 曲線，用以展示構(gòu)建的風(fēng)險評估模型的敏感性和特異性。在ROC 曲線的轉(zhuǎn)折點選擇最佳風(fēng)險評分臨界值，轉(zhuǎn)折點處真陽性和假陽性之間的差異最大。高于臨界值的患者屬于高危組，低于臨界值的患者屬于低危組。使用未配對t檢驗估計兩組正態(tài)分布變量的統(tǒng)計顯著性，并使用R 包Survminer 繪制兩組的生存曲線。

1.5 構(gòu)建和驗證列線圖 為了提高預(yù)后模型的性能，通過整合風(fēng)險評分模型和臨床信息，包括年齡、性別和腫瘤分期，可視化不同患者特征的預(yù)后價值，構(gòu)建列線圖。該分析使用R 軟件包rms 繪制校準(zhǔn)曲線，以評估預(yù)測概率，并與理想預(yù)測線進(jìn)行比較。此外，基于單因子回歸分析的森林圖說明了臨床信息與OS 之間的關(guān)系。其中，一致性指數(shù)（C-index）表明了列線圖的預(yù)測準(zhǔn)確性。

2 結(jié)果

2.1 差異基因分析 通過比較測試腫瘤樣本組和正常樣本組，TCGA 數(shù)據(jù)集和三個GEO 數(shù)據(jù)集的差異基因數(shù)量分布見表3 和圖1，可以看到共同的上調(diào)差異基因為48 個，共同的下調(diào)差異基因數(shù)為77 個。

表3 差異表達(dá)基因的統(tǒng)計信息

圖1 三個GEO 數(shù)據(jù)集和TCGA 數(shù)據(jù)集的差異基因數(shù)量情況

2.2 結(jié)腸癌預(yù)后模型的建立 通過單因子回歸分析表明，有14 個DEGs 與總生存期（OS）有關(guān)（P＜0.05），見表4。進(jìn)一步LASSO 回歸將基因數(shù)量縮減為10 個，見圖2。根據(jù)逐步回歸模型，Akaike 信息標(biāo)準(zhǔn)（AIC）為995.94，C 指數(shù)為0.63。

圖2 LASSO 回歸分析結(jié)果

表4 與結(jié)腸癌預(yù)后相關(guān)的基因信息

2.3 結(jié)腸癌預(yù)后模型的性能評估 使用預(yù)后模型公式計算每個結(jié)腸癌患者的風(fēng)險評分，然后根據(jù)R軟件包survminer 中預(yù)后因子相關(guān)風(fēng)險評分的最佳臨界值（cut-off 為0.16），將結(jié)腸癌患者分為高評分組和低評分組。結(jié)果顯示，隨著風(fēng)險得分的增加，生存時間呈現(xiàn)縮短的趨勢，并且高危組的死亡比例比低危組高，高風(fēng)險評分患者的OS 比低評分患者預(yù)后更差，其中基因CILP 和C7 在低風(fēng)險組表達(dá)量低，在高風(fēng)險組表達(dá)量低，而其余8 個基因趨勢相反，見圖3。

圖3 預(yù)后風(fēng)險評分分組和評估

2.4 預(yù)后模型的統(tǒng)計分析 為了進(jìn)一步評估預(yù)后風(fēng)險評分模型的性能，繪制ROC 曲線和腫瘤分層分析。通過將風(fēng)險評分的預(yù)后準(zhǔn)確性作為一個連續(xù)變量進(jìn)行研究，OS 預(yù)后模型的ROC 曲線下面積（AUC）在3 年時為0.628，4 年時為0.678，5 年時為0.730，見圖4。Wilcoxon 檢驗表明，較高的風(fēng)險評分與較高的病理分期（P=0.0019）、T 分期（P=0.049）、M分期（P=0.003）、N 分期（P=0.0015）相關(guān)。

圖4 預(yù)后風(fēng)險評估模型性能評價

2.5 列線圖模型的構(gòu)建與驗證 在列線圖中，每個變量的得分可以在分?jǐn)?shù)表上找到，然后通過計算總分來估計3 年、4 年和5 年的生存概率，見圖5A。森林圖顯示患者特征，包括年齡（＞60）、腫瘤分期（Ⅲ和Ⅳ）和風(fēng)險評分與OS 相關(guān)（P＜0.05），見圖5B。為了驗證列線圖的性能，繪制校準(zhǔn)曲線，可觀察到預(yù)測曲線接近理想曲線，性能良好，見圖5C～圖5E。此外，該列線圖（C-index：0.74）的預(yù)測準(zhǔn)確性高于風(fēng)險評分模型（C-index：0.63）。

圖5 列線圖模型的構(gòu)建與驗證

3 討論

本研究中利用3 個GEO 數(shù)據(jù)集和TCGA 數(shù)據(jù)集來挖掘共同的差異表達(dá)基因。其中，GEO 數(shù)據(jù)集是基于基因芯片平臺，TCGA 數(shù)據(jù)集是基于二代測序平臺，不同數(shù)據(jù)集交叉驗證，使得得到的差異表達(dá)基因具有相對廣泛適用性。回歸分析研究表明，有10 個差異表達(dá)基因與結(jié)腸癌預(yù)后顯著相關(guān)。其中，SLC4A4 全稱Solute Carrier Family 4 Member 4，已有研究報道該基因在結(jié)腸癌患者中低表達(dá)，與結(jié)腸癌癌較差的預(yù)后相關(guān)，也與淋巴結(jié)浸潤和遠(yuǎn)處轉(zhuǎn)移有關(guān)[10，11]。CD177 被認(rèn)為是一種干細(xì)胞因子受體，有實驗證明CD177 可作為結(jié)直腸癌患者對含貝伐單抗的抗癌治療反應(yīng)的潛在預(yù)測生物標(biāo)記物[12]。另有研究[13]提出CD177 調(diào)節(jié)胃癌中的腫瘤細(xì)胞粘附和遷移，是生存預(yù)后的因素。有報道表明[14]，CD177 的異常表達(dá)與結(jié)腸癌的發(fā)生發(fā)展相關(guān)。C7 作為一種潛在的腫瘤抑制因子，被報道與前列腺癌的免疫相關(guān)預(yù)后生物標(biāo)志物[15]。此外，細(xì)胞膜上表達(dá)的C7 是過度促炎反應(yīng)的調(diào)節(jié)因子[16]，而非小細(xì)胞肺癌（NSCLC）中C7 的低表達(dá)可能是腫瘤抑制劑，與腫瘤進(jìn)展和預(yù)后相關(guān)[16]。UGT2A3 是葡萄糖醛酸轉(zhuǎn)移酶（UGT）家族成員之一，UDP-葡萄糖醛酸轉(zhuǎn)移酶負(fù)責(zé)外源和內(nèi)源性化合物的葡萄糖醛酸化，包括藥物、環(huán)境麻醉劑、類固醇、神經(jīng)遞質(zhì)、膽汁酸和其他激素。據(jù)報道，UGT2A3 主要在與藥物清除有關(guān)的組織中表達(dá)水平最高，其中肝臟是表達(dá)最多的器官，其次是胃腸道和腎臟。研究表明[17]，原發(fā)性結(jié)腸癌肝轉(zhuǎn)移患者的UGT2A3 水平明顯高于無肝轉(zhuǎn)移患者。DNASE1L3 與自身免疫性疾病相關(guān)，它可以消化凋亡細(xì)胞釋放的微粒中的染色質(zhì)，這是一種潛在的自身抗原，它的過度積累將導(dǎo)致身體產(chǎn)生自身免疫反應(yīng)[18，19]。有研究表明，DNASE1L3 可能是結(jié)腸癌免疫浸潤的重要生物標(biāo)志物，并將為結(jié)腸癌免疫治療靶點的選擇提供理論依據(jù)[20]。HEPACAM2 是粘附基因免疫球蛋白家族的成員，已有研究表明該基因與結(jié)腸癌的發(fā)生發(fā)展有關(guān)[21]。ITLN1 可作為多種癌癥的腫瘤抑制因子，如胃癌[22]、卵巢癌[23]、神經(jīng)母細(xì)胞瘤[24]和結(jié)腸癌[25]。有研究通過免疫組織化學(xué)發(fā)現(xiàn)[25]，148例大腸癌中87 例（59%）ITLN1 蛋白低表達(dá)，ITLN1表達(dá)低的結(jié)腸癌患者的M 分級高于ITLN1 表達(dá)高的結(jié)腸癌患者（P=0.0017），且ITLN1 表達(dá)高的患者比ITLN1 表達(dá)低的患者預(yù)后更為良好。MMP3 和MMP10 均屬于基質(zhì)金屬蛋白酶（matrix metalloproteinase，MMP）家族，該家族成員參與正常生理過程中細(xì)胞外基質(zhì)的分解，如胚胎發(fā)育、生殖和組織重構(gòu)，以及疾病的發(fā)生發(fā)展。研究表明[26，27]，MMP3 和MMP10 均可作為結(jié)腸癌的預(yù)后相關(guān)。CILP 基因編碼軟骨中間層蛋白，在早期骨關(guān)節(jié)病軟骨中增加，目前還未見該基因與腫瘤的相關(guān)性報道。

綜上所述，本研究基于以上10 個預(yù)后相關(guān)的基因，構(gòu)建了結(jié)腸癌風(fēng)險評分模型和列線圖模型，結(jié)果表明，構(gòu)建的模型在結(jié)腸癌預(yù)后和腫瘤分層中表現(xiàn)良好，具有一定的應(yīng)用價值。