范以東，秦剛，劉金富，蘇國威，肖世富，劉俊良，李威材，吳廣濤

1廣西中醫(yī)藥大學研究生院，南寧 530222

2廣西中醫(yī)藥大學第一附屬醫(yī)院骨病創(chuàng)傷骨科，南寧 530022

骨肉瘤（osteosarcoma，OS）是一種主要由間充質(zhì)細胞引起的惡性腫瘤，主要發(fā)生在長管狀骨中[1]。OS 是發(fā)病率僅次于多發(fā)性骨髓瘤的第二常見的原發(fā)性骨惡性腫瘤，占所有骨腫瘤的15%[2]。隨著醫(yī)療水平的提升，OS 患者的5 年生存率可達到70%，但治療效果不佳且轉(zhuǎn)移的患者，其5 年生存率低于30%[3]。雖然新輔助化療的出現(xiàn)提高了OS 患者的5 年生存率，但OS 患者的生存情況仍然較差。二代測序（next-generation sequencing，NGS）技術的出現(xiàn)徹底改變了臨床研究、基礎醫(yī)學和應用醫(yī)學，通過NGS 技術有可能發(fā)現(xiàn)新的預后生物標志物，為OS 患者的個體化治療提供依據(jù)，有助于預測患者預后和提高治療效果[4]。因此，確定新的靶點至關重要。鐵死亡是一種鐵依賴性，區(qū)別于壞死和凋亡的細胞死亡模式，其特征在于脂質(zhì)活性氧（reactive oxygen species，ROS）的積累[5]，由Stockwell 等于2012 年首次提出[6]。鐵死亡在細胞形態(tài)學和生物學方面各有特點，形態(tài)學方面表現(xiàn)為線粒體體積變小、雙層膜密度增加、線粒體嵴減少或消失，生物學方面主要表現(xiàn)為谷胱甘肽過氧化物酶4（glutathione peroxidase 4，GPX4）失活，導致脂質(zhì)過氧化物積累，F(xiàn)e2+氧化脂質(zhì)產(chǎn)生ROS，從而發(fā)生鐵死亡[7]。研究表明，鐵死亡可以影響多種腫瘤的發(fā)生發(fā)展，包括結(jié)直腸癌[8]、膠質(zhì)母細胞瘤[9]、乳腺癌[10]、胃癌[11]等。本研究應用生物信息學方法通過公共數(shù)據(jù)庫探討鐵死亡相關基因（ferroptosis related gene，F(xiàn)RG）在OS 中的表達及與患者預后的關系，并分析FRG 在OS 中的作用機制，構(gòu)建FRG相關的OS 預后模型，旨在為OS 的診斷、治療及預后評估提供參考依據(jù)，現(xiàn)報道如下。

1 材料與方法

1.1 數(shù)據(jù)來源與整理

從UCSC Xena 數(shù)據(jù)庫（http://xena.ucsc.edu/）中獲取88例OS患者的轉(zhuǎn)錄組測序數(shù)據(jù)和其中85例患者的臨床資料（2 例患者生存狀態(tài)不明，1 例患者總生存時間不明，故予以刪除）。在基因型-組織表達（Genotype-Tissue Expression，GTEx）數(shù)據(jù)庫（https://www.gtexportal.org/home/）中獲取396例正常骨組織樣本，對所有GTEx 數(shù)據(jù)進行l(wèi)og2（x+1）轉(zhuǎn)換。使用“l(fā)imma”包對數(shù)據(jù)進行合并和矯正。從FerrDb數(shù)據(jù)庫（http://www.zhounan.org/ferrdb）中下載FRG，所有數(shù)據(jù)下載截止日期為2022年12月12日。

1.2 獲取差異表達的FRG

以錯誤發(fā)現(xiàn)率（false discovery rate，F(xiàn)DR）＜0.05，|log2FC|≥2 為過濾條件篩選正常骨組織與OS組織中的差異表達基因。從篩選的差異基因列表中進一步提取差異表達的FRG。

1.3 基因富集分析

利用R 語言“org.Hs.eg.db”“enrichplot”包對差異表達的FRG 進行基因本體（Gene Ontology，GO）和京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genomes，KEGG）富集分析，探索OS 中FRG 的生物學功能。

1.4 預后模型的構(gòu)建

采用單因素Cox 回歸模型分析差異表達的FRG 對OS 患者總生存期的影響，初步篩選與預后相關的基因。然后使用R 語言“glmnet”包進行最小絕對收縮和選擇算子（least absolute shrinkage and selection operator，Lasso），最終依據(jù)最佳λ值構(gòu)建一個基于FRG 的預后模型。每例患者的風險值公式：風險值=，其中xi為每個FRG 的表達量，Coefi為相應的基因系數(shù)，將其求和后即為該患者的風險評分。根據(jù)風險評分的中位值將患者分為高風險組和低風險組。

1.5 預后模型的分析

應用R 語言“survival”包繪制高風險組和低風險組患者的生存曲線，比較兩組患者的生存情況，采用受試者工作特征（receiver operating characteristic，ROC）曲線分析總體風險評分對預后的預測價值。應用R 語言“pheatmap”包繪制兩組患者的風險評分、生存狀態(tài)分布圖和模型基因表達熱圖。對兩組患者進行主成分分析（principal component analysis，PCA），觀察兩組患者的分布是否具有差異。同時分析風險評分、轉(zhuǎn)移情況和性別是否為OS 患者預后的獨立影響因素。

1.6 單基因富集分析

為探究預后相關基因在OS 進展過程中的潛在機制，采用單基因基因集富集分析（gene set enrichment analysis，GSEA）預測9 個預后相關基因的KEGG 下游途徑。通過R 語言“clusterProiler”包對9 個預后相關基因進行單基因GSEA 分析，將基因表達矩陣中的正常樣本刪除，保留腫瘤樣本作為GSEA 分析的輸入文件。同時在MSigDB 數(shù)據(jù)庫（https://www.gsea-msigdb.org/gsea/msigdb）中下載“c2.cp.kegg.v7.5.1.symbols.gmt”基因集，用以評估可能的相關途徑。|NES|＞1、p-val＜0.05、q-val＜0.25 的通路被認為是顯著富集的，展示富集最顯著的前5 條通路。

1.7 預后相關基因的驗證

1.7.1 預后相關基因的選擇 通過查閱文獻發(fā)現(xiàn)，酰基輔酶A 合成酶家族成員2（acyl-CoA synthetase family member 2，ACSF2）、脂肪酸去飽和酶2（fatty acid desaturase 2，F(xiàn)ADS2）在OS 的相關研究中仍是空白，而參與預后模型構(gòu)建的其他基因與OS 的關系已有了一定的研究。因此，本研究選擇ACSF2、FADS2作為驗證的目標基因。

1.7.2 儀器與設備 U2OS、MG63 細胞系購自武漢普諾賽生命科技有限公司，hFOB1.19（人SV40 轉(zhuǎn)染成骨細胞）及hFOB1.19 細胞專用培養(yǎng)基均購自上海滬震生物科技有限公司，去基因組與逆轉(zhuǎn)錄一管化三代預混液（MR05101M）購自莫納生物科技有限公司，熒光定量聚合酶鏈反應（quantitative polymerase chain reaction，qPCR）儀Light-Cycler?96 購自瑞士Roche 公司，梯度PCR 儀Biometra Tone 96 購自德國耶拿公司，超微量分光光度計NANODROP ONE 購自美國Thermo Fisher Scientific 公司。

1.7.3 細胞培養(yǎng) 以hFOB1.19 成骨細胞系作為對照組，OS 細胞系U2OS 和MG63 作為實驗組。應用hFOB1.19 細胞專用培養(yǎng)基培養(yǎng)hFOB1.19 細胞，培養(yǎng)基包含DMEM/F12、0.3 g/L G418 及10%胎牛血清，將細胞置于34 ℃、5%CO2培養(yǎng)箱中孵育，隔2天更換1 次新的培養(yǎng)基。當瓶壁細胞密度達90%左右時，按1∶2 傳代。U2OS、MG63 細胞在37 ℃、5%CO2培養(yǎng)箱中孵育。

1.7.4 qPCR 檢測ACSF2、FADS2基因表達 當細胞培養(yǎng)至密度達90%左右時，應用Trizol 法提取樣本總RNA，并采用MR05101M 合成互補DNA（complementary DNA，cDNA），然后進行qPCR。反應條件：95 ℃預變性10 min，95 ℃變性10 s，60 ℃退火30 s，共40 個循環(huán)。以β-肌動蛋白（β-actin）作為內(nèi)參，采用2-△△Ct法計算ACSF2、FADS2mRNA 相對表達量。設置3 個復孔并計算均數(shù)和標準差，引物序列見表1。

表1 β-actin、FADS2、ACSF2 引物序列

1.8 統(tǒng)計學方法

采用R 語言4.1.1 版對生物信息學數(shù)據(jù)進行統(tǒng)計分析，采用Lasso 回歸分析篩選預后相關基因并構(gòu)建預后模型；采用Kaplan-Meier 法繪制生存曲線，組間比較采用Log-rank 檢驗；采用ROC 曲線評估預后模型對OS 患者1、3、5 年生存率的預測價值；采用單因素和多因素Cox 回歸模型分析OS 患者總生存期的獨立影響因素。采用GraphPad Prism 8.4.0 軟件對qPCR 實驗數(shù)據(jù)進行統(tǒng)計分析，符合正態(tài)分布的計量資料以均數(shù)±標準差（±s）表示，多組間比較采用單因素方差分析，多組間兩兩比較采用LSD-t檢驗。以P＜0.05 為差異有統(tǒng)計學意義。

2 結(jié)果

2.1 差異表達的FRG

通過差異分析共獲得表達上調(diào)基因1404 個，下調(diào)基因1597 個。從FerrDb 數(shù)據(jù)庫中共下載得到258 個FRG，將其與正常骨組織和OS 組織中的差異表達基因取交集，共得到57 個差異表達的FRG。（圖1）

圖1 差異表達基因的箱線圖

2.2 GO 及KEGG 富集分析

GO 富集分析結(jié)果顯示，在生物過程（biological process，BP）方面，F(xiàn)RG 主要與對營養(yǎng)水平的反應、缺氧反應、細胞對化學應激的反應等有關；在細胞成分（cellular component，CC）方面，F(xiàn)RG 主要與線粒體外膜、細胞器外膜、黑素體等有關；在分子功能（molecular function，MF）方面，F(xiàn)RG 主要與鐵離子結(jié)合、泛素樣蛋白連接酶結(jié)合、轉(zhuǎn)氨酶活性等有關。KEGG 通路富集分析結(jié)果顯示，F(xiàn)RG 主要參與多種疾病的神經(jīng)退行性變、線粒體自噬、化學致癌-活性氧以及鐵死亡等途徑。（圖2）

圖2 OS中FRG基因的GO和KEGG富集分析

2.3 預后模型的構(gòu)建

采用單因素Cox 回歸模型分析差異表達的FRG對患者總生存期的影響，共得到12個預后相關基因，繪制森林圖（圖3A）。隨后將上述12 個基因納入Lasso 回歸模型進一步篩選，回歸模型的類型為Cox，采用10 倍交叉驗證，最大迭代次數(shù)設為1000，篩選最佳lambda 值來確定最優(yōu)模型，最終得出預后最相關的9 個FRG（圖3B、3C）。其中，B 細胞淋巴瘤/白血病-2 相互作用蛋白3（B-cell lymphoma/leukemia-2 interacting protein 3，BNIP3）、FADS2、血管內(nèi)皮生長因子A（vascular endothelial growth factor A，VEGFA）為高風險基因，ACSF2、芳香烴受體核轉(zhuǎn)運因子樣蛋白（aryl hydrocarbon receptor nuclear translocator like，ARNTL）、葡萄糖-6-磷酸脫氫酶（glucose-6-phosphate dehydrogenase，G6PD）、細胞因子信號轉(zhuǎn)導抑制因子1（suppressor of cytokine signaling 1，SOCS1）、轉(zhuǎn)化生長因子β受體1（transforming growth factor beta receptor 1，TGFBR1）、磷酸葡萄糖酸脫氫酶（phosphogluconate dehydrogenase，PGD）為低風險基因，分別計算每例患者的風險評分，按照風險評分的中位值將患者分為高風險組（42 例）和低風險組（43 例）。

圖3 預后模型基因的篩選

2.4 預后模型分析

將高風險組和低風險組患者進行PCA，從圖中點的分布可以看出兩組患者的分布具有顯著差異，以PC1 為分割線，低風險組患者主要分布在圖中的左側(cè)，高風險組患者主要分布在右側(cè)（圖4A）。生存分析顯示，低風險組患者的生存情況明顯優(yōu)于高風險組，差異有統(tǒng)計學意義（P＜0.01）（圖4B）。使用ROC 曲線驗證風險評分對OS 患者1、3、5 年生存率的預測價值，結(jié)果顯示，風險評分預測OS 患者1、3、5 年生存率的曲線下面積分別為0.815、0.845、0.838（圖4C）。在風險曲線中，隨著風險評分的增大，死亡患者逐漸增加，生存時間逐漸縮短（圖5A、5B）。風險熱圖顯示，隨著風險評分的增大，BNIP3、FADS2、VEGFA的表達量增加，ACSF2、ARNTL、G6PD、SOCS1、TGFBR1、PGD的表達量降低（圖5C）。

圖4 高風險組（n=42）和低風險組（n=42）OS患者的生存分析與預后預測

圖5 預后模型中OS患者生存風險狀態(tài)圖

2.5 預后影響因素分析

單因素和多因素分析結(jié)果顯示，轉(zhuǎn)移和風險評分均是OS 患者預后的獨立影響因素（P＜0.01），性別不是OS患者預后的獨立影響因素（P＞0.05）。（圖6）

圖6 OS患者預后的影響因素分析

2.6 單基因GSEA 分析

GSEA 結(jié)果顯示，這些預后相關基因主要與原發(fā)性免疫缺陷、趨化因子信號通路、嗅覺轉(zhuǎn)導、自噬的調(diào)控、維甲酸誘導基因-Ⅰ（retinoic acid inducible gene-Ⅰ，RIG-Ⅰ）樣受體信號通路、剪接體、細胞質(zhì)DNA 感應通路等密切相關。其中，PGD與VEGFA未發(fā)現(xiàn)符合條件的顯著富集通路。（圖7）

圖7 OS患者預后相關基因的單基因GSEA分析

2.7 ACSF2、FADS2 mRNA 相對表達量的比較

U2OS 和MG63 細胞中FADS2mRNA 相對表達量均高于hFOB1.19細胞，差異均有統(tǒng)計學意義（P＜0.05）；MG63 細胞中ACSF2mRNA 相對表達量高于hFOB1.19 細胞，差異有統(tǒng)計學意義（P＜0.05）。hFOB1.19 和U2OS 細胞中ACSF2mRNA 相對表達量比較，差異無統(tǒng)計學意義（P＞0.05）。（表2）

表2 不同細胞中ACSF2、FADS2 mRNA相對表達量的比較（±s）

表2 不同細胞中ACSF2、FADS2 mRNA相對表達量的比較（±s）

注：*與hFOB1.19細胞比較，P＜0.05

細胞hFOB1.19 U2OS MG63 F值P值ACSF2 mRNA相對表達量1.000±0.000 2.349±1.285 10.416±0.972*89.79＜0.01 FADS2 mRNA相對表達量1.000±0.000 5.242±1.229*2.836±0.491*23.25＜0.01

3 討論

雖然新輔助化療的出現(xiàn)提高了OS 患者的無復發(fā)生存率和總生存率，但隨后的幾十年間，OS 的治療效果并沒有得到太大提升，似乎到了瓶頸期[12]。鐵死亡自2012 年被提出后，10 余年間相關研究不斷增多，逐漸成了醫(yī)學領域的研究熱點，在多種類型的腫瘤中均有相應研究。雖然有研究表明鐵死亡與OS 的發(fā)生發(fā)展關系密切，但其作用機制尚不清楚。通過生物信息學探索與預后有關的FRG 有望成為一種突破OS 治療瓶頸的新方法。

本研究篩選出57 個差異表達的FRG，通過GO與KEGG 富集分析發(fā)現(xiàn)，這些FRG 與缺氧反應、線粒體外膜、細胞器外膜、鐵離子結(jié)合等有關，F(xiàn)RG主要參與多種疾病的神經(jīng)退行性變、線粒體自噬、化學致癌-活性氧以及鐵死亡等途徑，這些均與腫瘤的發(fā)生發(fā)展密切相關。缺氧反應在體內(nèi)腫瘤微環(huán)境中具有重要地位。有研究表明，缺氧條件促進了OS 細胞的侵襲和遷移[13]。

本研究篩選了9 個基因進行預后模型的構(gòu)建。ACSF2 主要位于線粒體內(nèi)，可以通過與輔酶A（coenzyme A，CoA）形成硫酯來催化脂肪酸代謝中的初始反應[14]。研究發(fā)現(xiàn)，ACSF2 在潰瘍性結(jié)腸炎中表達顯著下調(diào)，且與免疫相關通路顯著相關，鐵死亡抑制劑可逆轉(zhuǎn)ACSF2的表達，證明ACSF2在介導鐵死亡過程中具有重要作用[15]。FADS2是鐵死亡的限速因子，具有催化脂肪酸去飽和的作用，可誘導宿主細胞中的鐵死亡，并阻礙丙型肝炎病毒的復制[16]。ARNTL是哺乳動物生物鐘的核心組成部分，可通過抑制egl-9家族缺氧誘導因子2（egl-9 family hypoxia inducible factor 2，EGLN2）的轉(zhuǎn)錄來抑制鐵死亡，而阻斷ARNTL降解或抑制EGLN2活化可以誘導鐵死亡發(fā)生[17]。G6PD 被認為是肝細胞癌患者預后的危險因素，可通過抑制細胞色素p450 氧化還原酶（cytochrome p450 oxidoreductase，POR）的表達促進腫瘤細胞生長、侵襲和轉(zhuǎn)移，同時減少鐵死亡[18]。PGD可以減少細胞中的糖酵解，增加戊糖磷酸途徑（pentose phosphate pathway，PPP），從而保護細胞免受ROS 的侵害[19]。鐵死亡與ROS 關系密切，未來需要更多研究探索PGD 與鐵死亡之間的關系。SOCS1是細胞因子信號轉(zhuǎn)導抑制因子家族成員，在細胞生長和分化過程中發(fā)揮重要作用，可通過與p53 相互作用抑制腫瘤發(fā)生發(fā)展[20]。溶質(zhì)載體家族7成員11（solute carrier family 7 member 11，SLC7A11）和亞精胺/精胺N1-乙酰轉(zhuǎn)移酶1（spermidine/spermine N1-acetyltransferase 1，SAT1）在鐵死亡過程中具有關鍵作用，且SOCS1 能夠抑制谷胱甘肽（glutathione，GSH）的表達，證明SOCS1 能夠提高細胞對鐵死亡的敏感性[21]。TGFBR1又稱激活素受體樣激酶（activin receptor-like kinase，ALK）4/5，在細胞中主要參與氧化應激反應[22]。在腎近端腎小管上皮細胞中，ALK4/5 信號通路已被證明與鐵死亡相關，阻斷ALK4/5 信號通路可抑制鐵死亡[23]。VEGFA 被認為是一種重要的血管生長因子，研究發(fā)現(xiàn)，當子宮內(nèi)膜基質(zhì)細胞發(fā)生鐵死亡時，p38-促分裂原活化的蛋白激酶（mitogen-activated protein kinase，MAPK）/信號轉(zhuǎn)導及轉(zhuǎn)錄激活因子6（signal transduction and activator of transcription 6，STAT6）信號轉(zhuǎn)導會抑制VEGFA 和白細胞介素-8（interleukin-8，IL-8）的表達[24]。因此，鐵死亡可能通過旁分泌VEGFA和IL-8對鄰近病變產(chǎn)生血管生成作用，從而在子宮內(nèi)膜異位癥的進展過程中發(fā)揮關鍵作用。

本研究成功構(gòu)建了基于FRG 的骨肉瘤風險預后模型，發(fā)現(xiàn)ACSF2、ARNTL、G6PD、PGD、FADS2、SOCS1、BNIP3、TGFBR1、VEGFA等9 個FRG 能夠作為OS 患者的預后生物標志物，并采用qPCR 技術驗證了ACSF2和FADS2在OS 細胞與正常成骨細胞中的表達情況，為OS 患者的臨床治療和預后評估提供了參考。然而本研究也存在一定的局限性，需要更多的樣本及基礎研究進一步驗證。