焦趙爽，張懷念

0 引言

甲狀腺癌是人類內分泌系統(tǒng)和頭頸部最常見的惡性腫瘤之一，其發(fā)病率在國內外均呈上升趨勢[1]。許多危險因素與甲狀腺癌的發(fā)病相關，包括電離輻射、雌激素、碘攝入和遺傳因素等。研究表明，手術聯(lián)合放化療以及分子靶向治療有助于控制甲狀腺癌進展，延長患者的無進展生存期，提高總體生存率[2-4]。然而由于早期缺乏典型的臨床表現(xiàn)，甲狀腺癌很容易被患者自身忽略，從而錯過最佳治療時間[5]。分化型甲狀腺癌一般預后較好，但當其失分化為碘難治性甲狀腺癌時便進展迅速，病死率較高，預后較差[6]。乳頭狀甲狀腺癌（papillary thyroid carcinoma,PTC）是最常見的甲狀腺惡性腫瘤，約占甲狀腺癌的80%[7]。據(jù)報道，在發(fā)達國家PTC的發(fā)病率呈上升趨勢[8-10]，PTC的不良預后因素包括年齡、性別、腫瘤大小等[7]。Grogan等[11]對269例PTC患者進行平均27年的長期隨訪，結果發(fā)現(xiàn)超過25%的患者復發(fā)，并且在11%的病例中，PTC復發(fā)發(fā)生在治療后20多年。因此，評估患者的疾病進展情況，尋找影響PTC患者預后的相關因素值得重視[12]。

微小RNA（miRNA）是一類短的非編碼RNA分子，僅含有19～22個核苷酸。它們是高度保守的分子，通過以完全或不完全互補的方式與靶基因的3’端非翻譯區(qū)（3’-UTR）結合，參與靶基因的調節(jié)[13]。越來越多的研究證實，miRNA廣泛參與人類癌癥的發(fā)生和進展，包括甲狀腺癌。miRNA在腫瘤發(fā)生和發(fā)展過程中參與細胞增殖、凋亡和腫瘤侵襲等多種惡性生物學行為[14]。雖有研究報道某些miRNA可以作為PTC的分子標志物，但目前仍缺乏敏感的miRNA組合標志物以評估PTC患者的預后。因此，本研究通過使用TCGA數(shù)據(jù)庫，結合生物信息學分析對PTC組織中差異表達的miRNA進行篩選，利用Cox和Lasso回歸分析出與PTC患者預后相關的miRNA，構建基于miRNA表達的預測PTC患者預后的生存模型。

1 資料與方法

1.1 原始數(shù)據(jù)的下載與處理

在TCGA官方網(wǎng)站上下載PTC患者的miRNA測序數(shù)據(jù)和患者的臨床數(shù)據(jù)。其中癌組織514例，正常甲狀腺組織59例。使用Perl5.24.3軟件將原始miRNA測序數(shù)據(jù)轉換成miRNA表達矩陣。在分析過程中，若有多個探針檢測同一miRNA表達量，則取該miRNA表達量的平均值作為該miRNA的表達值。對于患者臨床資料的分析，刪除生存時間未知和生存時間為0的患者臨床信息。

1.2 差異miRNA的篩選

以|log foldchange|≥2，錯誤發(fā)現(xiàn)率FDR＜0.05為篩選條件，使用R3.6.0軟件中的edgeR包篩選出在甲狀腺癌中差異表達的miRNA。使用ggplot2和heatmap軟件程序包繪制差異基因的火山圖和熱圖。

1.3 Cox/Lasso回歸分析

結合患者的生存信息，首先對上述得到的差異miRNA進行單因素Cox回歸分析，計算每個miRNA與甲狀腺癌患者生存的風險比（hazard ratio,HR）和P值，以P＜0.05的標準篩選出與甲狀腺癌患者預后顯著相關的miRNA。將這些miRNA進一步行Lasso回歸分析，目的減少基因之間共線性的影響，防止后續(xù)構建的風險模型變量過度擬合。Lasso回歸使用交叉驗證以確定參數(shù)，得到合適的模型。再將Lasso回歸得到的miRNA進行多因素Cox回歸分析，計算每個miRNA的多因素回歸系數(shù)，構建風險評分方程。

1.4 風險預后模型的建立與分析

根據(jù)上述多因素Cox回歸分析的結果，構建基于miRNA表達的風險評分方程risk score。按照文獻報道使用公式：Risk score=β1×miRNA1EXP+β2×miRNA2EXP+......+βn×miRNAnEXP[15-16]。式中β為相應miRNA的多因素回歸系數(shù)，miRNAEXP為相應miRNA的表達量。根據(jù)risk score數(shù)值的中位值，將PTC患者分為高風險評分組和低風險評分組。利用R3.6.0軟件繪制模型預測預后的列線圖，并比較兩組患者之間生存的差異。利用R3.6.0軟件繪制模型的ROC曲線和校準曲線以評價模型的敏感度和特異性。使用Survival ROC軟件程序包計算受試者工作特征曲線（ROC）下面積（AUC）的數(shù)值。

2 結果

2.1 差異miRNA的篩選

利用|logfoldchange|≥2，F(xiàn)DR＜0.05為篩選條件，在PTC組織中共篩選到差異表達的miRNA75個。其中，上調表達的有70個、下調表達的有5個。根據(jù)FDR數(shù)值的排序前十位差異表達的miRNA，見表1。圖1為差異miRNA相應的火山圖和熱圖。圖中紅色表示與正常甲狀腺組織相比，該基因在PTC組織中表達上調；綠色表示與正常甲狀腺組織相比，該基因在PTC組織中表達下調。

表1 PTC組織中前十位差異表達的miRNATable 1 Top 10 differentially-expressed miRNA in papillary thyroid carcinoma (PTC) tissues

2.2 Cox/Lasso回歸分析

圖1 PTC組織中差異表達的miRNA火山圖(A)和熱圖(B)Figure 1 Volcano plot(A) and heat map(B) of differentially-expressed miRNA in PTC tissues

對上述差異表達的miRNA進行單因素Cox回歸分析，計算相應miRNA與甲狀腺癌患者的HR和P值，結果共有9個miRNA與甲狀腺癌患者生存相關（P＜0.05）：hsa-mir-6730、hsamir-4709、hsa-mir-196a-2、hsa-mir-146b、hsamir-6860、hsa-mir-509-3、hsa-mir-513c、hsamir-515-1、hsa-mir-551b。進一步使用lasso回歸對這些miRNA進行篩選，見圖2A，并使用交叉驗證建立模型，見圖2B。結果共有8個miRNA（hsamir-6730,hsa-mir-4709,hsa-mir-196a-2,hsa-mir-146b,hsa-mir-6860,hsa-mir-509-3,hsa-mir-513c,hsa-mir-515-1）納入分析模型。

2.3 風險評分方程的建立

對上述8個miRNA進行多因素Cox回歸分析，根據(jù)相應的回歸系數(shù)，建立風險評分方程。

Risk score=0.41×hsa-mir-196a-2EXP-0.14×hsa-mir-146bEXP-0.22×hsa-mir-4709EXP+0.83×hsa-mir-509-3EXP-0.03×hsa-mir-513cEXP+0.36×hsa-mir-515-1EXP-0.33×hsa-mir-6730EXP-0.63×hsa-mir-6860EXP。根據(jù)此方程，計算每位甲狀腺癌患者的risk score數(shù)值，并根據(jù)risk score數(shù)值的中位值，將甲狀腺癌患者分為高風險評分組和低風險評分組。

2.4 風險預后模型分析

根據(jù)上述結果，利用R3.6.0軟件繪制基于8個miRNA的組合預測甲狀腺癌患者生存的列線圖，見圖3。

2.5 風險預后模型的評價

圖2 Lasso回歸分析篩選miRNA(A)和交叉驗證結果(B)Figure 2 miRNA screened by Lasso regression analysis(A)and cross validation results(B)

利用R3.6.0軟件繪制模型的ROC曲線以評價模型的敏感度和特異性，見圖4。我們所構建的miRNA模型預測患者3年生存率和5年生存率的ROC曲線下面積數(shù)值AUC分別為0.860和0.896，這表明模型具有較好的敏感度和特異性。同時，校準曲線的結果也顯示構建的模型可靠，這些結果表明我們構建的miRNA預后模型可以較準確的預測甲狀腺癌患者的生存，見圖5。

2.6 生存分析

利用R3.6.0軟件中的survival程序包分析高低風險評分組患者的生存差異，繪制生存點圖，圖中綠色點表示未死亡患者，紅色點表示死亡患者。隨著risk score數(shù)值的升高，高風險評分組患者的死亡人數(shù)顯著多于低風險評分組，說明高風險評分組的患者生存率較差，見圖6。

3 討論

圖3 基于8個miRNA表達預測甲狀腺癌患者生存的列線圖Figure 3 Nomogram of thyroid cancer patients’ survival predicted based on eight miRNA expression

圖4 風險預后模型的ROC曲線Figure 4 ROC curve of risk prognosis model

圖5 預測模型的校準曲線Figure 5 Calibration curve of prediction model

圖6 高低風險評分組患者的生存點圖Figure 6 Survival points diagram of patients in high and low risk groups

甲狀腺癌是最常見的內分泌惡性腫瘤。大多數(shù)甲狀腺癌起源于甲狀腺濾泡細胞。甲狀腺癌可分為乳頭狀癌（75%～80%）、濾泡狀癌（10%～15%）、甲狀腺未分化癌（0.2%～2%）和由濾泡旁C細胞來源的髓樣甲狀腺癌（5%～10%）[17]。許多研究表明，miRNA與甲狀腺癌的發(fā)生和進展有關。Yue等[18]研究表明，miR-7可能通過p21活化激酶1（PAK1）調節(jié)甲狀腺癌細胞的生長、遷移和侵襲；Wang等[19]研究表明，碘可以通過下調miR-422a、上調MAPK1的表達，促進甲狀腺癌的發(fā)生；Dong等[20]研究發(fā)現(xiàn)miR-141可以通過調控胰島素受體底物2（IRS2）調控甲狀腺癌細胞生長和轉移，提示miR-141可以作為治療甲狀腺癌患者的潛在分子靶點；此外，miR-497也被認為是一種通過抑制腦源性神經營養(yǎng)因子發(fā)揮作用的甲狀腺癌腫瘤抑制因子[21]。因此，miRNA可以作為甲狀腺癌診斷、治療或者預后的標志物。

本研究中，我們通過TCGA數(shù)據(jù)庫構建了基于8個miRNA表達的組合預測甲狀腺癌患者預后的風險模型。在這8個miRNA中，有7個miRNA被報道與甲狀腺癌或者其他腫瘤相關。MiR-146b被證明在乳頭狀甲狀腺癌中表達上調，并且與腫瘤的侵襲有關[22-23]；Fu等[24]利用生物信息學的方法挖掘TCGA數(shù)據(jù)庫發(fā)現(xiàn)hsa-mir-196a-2在甲狀腺癌中表達顯著上調，其與腫瘤進展、淋巴結轉移和局部浸潤有關，并且可以作為甲狀腺癌的獨立預后危險因素；miR-509的過表達可以抑制三陰性乳腺癌Hs578T細胞增殖、誘導細胞凋亡和抑制細胞侵襲[25]；Tang等[26]利用生物信息學分析發(fā)現(xiàn)，miR-6860和miR-509-3與甲狀腺癌患者的總體生存率相關，可作為甲狀腺癌的預后因子；miR-513c被證明與神經母細胞瘤[27]、肝細胞癌[28]的細胞遷移、入侵和增殖等有關；has-mir-4709被證明與結腸癌[29]、胰腺導管腺癌[30]患者的生存相關。

綜上所述，我們構建了基于8個miRNA預測甲狀腺癌患者生存的風險模型，該模型具較好的敏感度和特異性，模型中高風險評分組患者的死亡風險較大。研究指出[31]，建立多基因預后模型比患者的臨床病理指標更能提供精確的預后指導，對個體化治療方案的選擇具有重要的參考價值。我們所構建的8個miRNA模型可有效預測甲狀腺癌患者的預后，但仍需要進一步使用大規(guī)模的多中心的循證醫(yī)學證據(jù)加以驗證。