基于TCGA數據庫探究可溶性因子對乳腺癌患者的預后作用

劉孟晨 孫貽安 李靜蔚

乳腺癌是以乳房組織中乳腺上皮細胞的增殖失控為主要表現的惡行異質性疾病。有統計表明，僅2018年在全球范圍內便有60多萬人死于乳腺癌[1]，且存活的人也因焦慮[2]、放療[3]、乳房切除[4]等因素極易誘發抑郁癥等精神類疾病。因此，乳腺癌對女性群體的生命及身心健康產生了極大的威脅。值得注意的是，乳腺癌患者的死亡原因大多可歸咎于轉移及相關并發癥[5]。在轉移的過程，腫瘤細胞從原發腫瘤脫離并進入血流[6]。這些循環腫瘤細胞最終滯留在遠處器官的毛細血管床中，導致繼發部位產生轉移性集落[5]。因此，可溶性因子在轉移等過程中發揮重要作用，它們的變化往往能加劇乳腺癌患者病情惡化甚至死亡的風險[7]。

可溶性因子是機體內可溶于血清的各類細胞因子，可通過增加炎性損傷、干預信號傳導、調節免疫等方式參與各類疾病的發展[8]。研究表明，骨橋蛋白、表皮生長因子和IL-6等可溶性因子會參與乳腺癌細胞和巨噬細胞之間的細胞相互作用[9]。成纖維細胞分泌的可溶性因子可募集單核細胞進而促進三陰性乳腺癌患者反應性基質的激活使其更具備侵襲性[10]。乳腺癌微環境下所分泌的可溶性因子甚至可觸發異常信號傳導，刺激骨微環境中的破骨細胞分化，進而加大骨轉移的風險[11]。這些研究表明，可溶性因子極大程度參與了乳腺癌的發病進展及預后，但大多研究并未基于整體性及系統性對可溶性因子的作用進行篩選和分析。基于此，本研究擬基于癌癥基因圖譜（the cancer genome Atlas,TCGA）數據庫中乳腺癌患者的基因表達等數據，探究可溶性因子在乳腺癌中的變化與影響。

1 資料和方法

1.1 數據的采集 通過TCGA數據庫（https://portal.gdc.cancer.gov/projects）下載乳腺癌患者的轉錄組數據及臨床數據。刪除生存數據缺失的樣本后，通過R語言中的limma包對腫瘤組織與正常組織進行差異分析，以FDR值＜0.05為條件篩選差異基因。通過基因集富集分析（gene set enrichment analysis,GSEA）數據庫（http://www.gsea-msigdb.org/gsea/msigdb/cards/NABA_SECRETED_FACTORS.html）下載可溶性因子數據集。并通過survival包對樣本的生存時間和生存狀態與可溶性因子基因相關表達數據進行生存分析得到可溶性因子預后相關基因。將預后相關可溶性因子數據集與TCGA中差異基因取交集后得到可溶性因子預后相關基因表達數據。以熱圖和森林圖的方式對可溶性因子預后相關基因的差異表達及對乳腺癌患者預后的影響進行展示。

1.2 蛋白質-蛋白質相互作用（protein-protein inter‐action,PPI）、基因本體論（gene ontology,GO）與京都基因和基因組百科全書（Kyoto Encyclopedia of Genes and Genomes,KEGG）分析 通過基因搜索工具（search tool for recurring instances of neighbouring genes,String）數據庫（https://string-db.org/cgi/input.pl）對交集基因進行PPI分析，并根據TCGA數據中乳腺癌患者轉錄組數據的表達計算可溶性因子預后相關基因的相關性，最后通過R語言對可溶性因子預后相關基因進行GO、KEGG分析。

1.3 風險模型的構建及主成分分析、隨機臨近嵌入分析 Glmnet包是構建Lasso回歸模型的重要工具，通過樣本中所有可溶性因子預后相關基因的表達及危險系數的乘積作為風險值為基礎建模。以所有樣本風險值的中位數為界限，劃分高低風險組。通過survival、survminer、timeROC包比較高低風險組患者的生存差異并獲取受試者工作特征曲線。使用Rtsne包對高低風險組進行主成分分析和T分布隨機臨近嵌入分析展現高低風險組基因的分布情況并進行降維可視化。

1.4 高低風險組免疫浸潤模式分析 基于單樣本基因富集分析（single sample GSEA,ssGSEA）的方法對每個樣本進行免疫打分，并根據模型中高低風險組中可溶性因子集內核心基因的表達進行免疫細胞及免疫功能的差異分析。

2 結果

2.1 數據的分析及處理 通過對TCGA數據集中的數據進行整理，發現數據庫中共有1 222個乳腺癌相關樣本，其中正常樣本有113個，乳腺癌樣本1 109個。GSEA數據庫中可溶性因子基因共有343個，其中有42個基因與乳腺癌患者的預后相關（圖1a）。將TCGA數據庫中乳腺癌患者體內的差異基因與預后相關的可溶性因子取交集后共得到27個預后相關的可溶性因子的差異基因（圖1b），并通過森林圖（圖1c）及熱圖（圖1d）的形式展現交集基因的風險因素及表達。從圖1c的結果中可以看出，除FGF8及WNT11為高風險致病基因外，其他可溶性因子均對患者預后影響較小，這可能是由于眾多可溶性基因綜合作用下單基因產生的生物學效應相互之間存在一定的抵消。從圖1d的結果來看，FGF8、WNT11、CCL23、FLT3LG、TNFSF12、IL-7、CXCL1、FGFBP1、IL-17B、S-100B、SFRP1、NRG1、IFNE、IL-4、CXCL14、IL-16、XCL1在正常組中表達明顯，但在腫瘤組中表達較低，故這些基因可能有抑癌的作用。CCL25、CCL5、CXCL13、CXCL9、IL-12B、IL-18、IL-24、LTA、LTB、TNFSF12、WNT7B等基因在正常組中表達較低但在腫瘤組中表達較高，這表明這些基因可能有抑癌的作用。

圖1 數據的分析與處理（a：42個可溶性因子預后相關基因；b：可溶性因子預后相關基因與差異基因的韋恩圖；c：預后相關可溶性因子森林圖；d：預后相關可溶性因子熱圖）

2.2 PPI網絡圖及GO、KEGG分析 PPI網絡圖顯示，在可溶性因子預后相關基因相互間存在較強的相互作用（圖2a），其中相關度≥10的靶點為IL-7、IL-4、IL-16、CCL5、CXCL13、CXCL9、CXCL1、LTA，這表明這些靶點與其他基因之間具有較強的相互作用。相關性分析發 現，FLT3LG、CXCL9、IL-12B、IL-18、IL-24、LTA、LTB、TNFSF12的表達之間存在明顯的正相關（圖2b）。GO分析表明,可溶性因子預后相關基因在生物過程中可以參與補體的激活，循環免疫球蛋白介導的體液免疫反應，體液免疫反應等（圖2c）。在細胞組分中涉及免疫球蛋白復合物、免疫球蛋白復合物，循環的漿膜外側血液微粒子等。在分子功能中多涉及抗原結合、免疫球蛋白受體結合、趨化因子活性等。KEGG分析（圖2d）表明，可溶性因子預后相關基因與病毒蛋白與細胞因子和細胞因子受體的相互作用、細胞因子-細胞因子受體的相互作用、趨化因子信號傳導途徑、造血細胞譜系、細胞粘附分子、NF-κB信號傳導途徑、Th1和Th2細胞分化、TNF信號傳導途徑、IL-17信號傳導途徑、癌癥中的PD-L1表達和PD-1檢查點通路、JAK-STATd等信號傳導途徑密切相關。

圖2 PPI、GO、KEGG分析（a：PPI分析圖；b：相關性分析圖；c：GO分析圖；d：KEGG分析圖）

2.3 風險模型的構建及主成分分析、隨機臨近嵌入分析 使用Lasso回歸分析構建風險模型，最終篩選出14個基因參與構建風險模型，風險模型計算公式為：風險值=CCL25×（-0.267 761 762 190 859）+CXCL1×（-0.0343324083608361）+CXCL13×（-0.0427683660862547）+CXCL14×（-0.0498671701605376）+FGF8×1.2580095296421+FGFBP1×（-0.0758992553296801）+FLT3LG×（-0.159 083 190 747 296）+IL-12B×（-0.088 335 107 649 380 9）+IL-24×（-0.0514706044485947）+IL-4×（-1.22386864779545）+NRG1× （-0.252 446 906 132 768 ）+S-100B×（-0.0277742304519763）+WNT11×（-0.0736641709388857）+WNT7B×0.123 111 236 484 031。

以風險值的中位數為界，將患者劃分為高、低風險兩組，依據患者的分組及風險值繪制患者的生存狀態（圖3a）及風險值熱圖（圖3b）。使用主成分分析、隨機臨近嵌入分析將患者的分組信息進行可視化（圖3c、3d）。通過比較兩組間患者的總生存率，發現高風險患者的生存時間顯著低于低風險組（圖3e，P＜0.01）。通過R語言繪制的ROC曲線發現，患者1～3年的AUC均＞0.65，這表明本模型較為可靠（圖3f）。最后通過單因素及多因素分析，發現可溶性因子預后相關基因的致病性與乳腺癌患者的年齡、性別、分期均顯著相關（圖3g、3h）。

圖3 風險模型的構建及主成分分析、隨機臨近嵌入分析（a：高、低風險組生存狀態圖；b：高、低風險組風險值熱圖；c：高、低風險組主成分分析圖；d：高、低風險組隨機臨近嵌入分析；e：高、低風險組生存分析；f：高、低風險組ROC曲線下面積圖；g：高、低風險組單因素Cox回歸；h：高、低風險組多因素Cox回歸）

2.4 高低風險組免疫浸潤模式分析 通過ssGSEA對高低風險組中免疫細胞及免疫功能進行了差異分析，研究表明可溶性因子可以很大程度上改變免疫細胞在高低風險組的評分。在高風險組中所有的免疫細胞及免疫功能評分均發生了降低。這表明aDCs、B細胞、CD8+T細胞、DCs、iDCs、巨噬細胞、肥大細胞、中性粒細胞、NK細胞、pDCs、輔助T細胞等免疫細胞在高風險患者中自身的功能有明顯的抑制，見圖4。

圖4 高低風險組免疫浸潤模式分析（a：高、低風險組免疫細胞差異分析；b：高、低風險組免疫功能分析；**P＜0.01）

3 討論

可溶性因子是人體內各類免疫細胞、內分泌細胞、神經細胞分泌的具有高度誘導性細胞因子，深度參與了乳腺癌發展、轉移、惡化等過程[12]。

本研究表明，可溶性因子能通過多種途徑影響乳腺癌的發生、發展。如CXCL1、CXCL9、CXCL13、CX‐CL14等可通過影響機體的免疫反應在乳腺癌中發揮作用。其中CXCL1是腫瘤相關巨噬細胞分泌的最豐富的趨化因子，可以將各種基質細胞募集到腫瘤環境中，構建促進腫瘤細胞生長、血管生成和轉移的環境[13]。CXCL9可通過調節乳腺癌中免疫細胞浸潤的豐度，進而干預乳腺癌的發展[14]。CXCL13是重要的趨化因子之一，B淋巴細胞聚集的標志物，在B淋巴細胞的歸巢、遷移和積累中起關鍵作用[15]。中性粒細胞被認為是提供抵御入侵病原體的主要免疫細胞，它們的募集在很大程度上被CXCL1和CXCL2等可溶性因子所誘導[16-17]。CCL5可募集T細胞、B細胞、嗜酸性粒細胞、嗜堿性粒細胞、中性粒細胞、巨噬細胞和成纖維細胞等各類免疫細胞[18-19]，這與本研究后續免疫浸潤分析的結果高度重合。CCL5還可通過PI3K/Akt通路激活αvβ3整合素以促進細胞遷移，而αvβ3整合素又能激活IKKα/β和NF-κB通路形成聯機反應[20]。當然，其他白介素家族所導致的炎癥損傷也不能被忽視。有研究表明，IL-4由巨噬細胞和T淋巴細胞產生，可控制T輔助細胞的成熟，IL-4誘導的谷氨酰胺代謝對乳腺癌的生長有較大的刺激作用[21]。IL-7最近參與了產生IL-17A的先天性樣T細胞的選擇性擴增和功能，包括自然殺傷T細胞和Th17細胞[22-23]。此外，IL-7水平與前列腺癌患者的不良預后相關[24-25]，與乳腺癌患者的腫瘤侵襲性正相關[26]。還有研究表明，FGF8作為一種可調節乳腺癌生長的一種可溶性生長因子，可通過調控類固醇激素的分泌影響乳腺癌的預后[27]。這些研究均顯示可溶性因子自身具有較高的致病性，而本研究基于這些關鍵基因作為預后模型在后續的驗證中也表現出極強的特異性。

更重要的是，可溶性因子不僅自身具有一定致病性，還可以與其他因子的相互作用下形成聯級反應。KEGG結果中，可溶性因子所涉及的NF-κB、Th1/Th2、PD-L1/PD-1等通路均與乳腺癌的發生發展密切相關。據報道，乳腺癌中NF-κB的激活可上調細胞周期蛋白D1、細胞周期蛋白依賴性激酶2（CDK2）和c-Myc[28-30]的表達，從而加速細胞周期進程并導致細胞增殖失控。NF-κB還調節IL-1β、TNFα、EGFR和HER2等可溶性因子的分泌以促進腫瘤細胞的生長[31]。NF-κB 活性的下降能改變Bcl-2[32]、IAP蛋白（XIAP、cIAP-1/2）[32]和 TNF 受體相關因子（TRAF）1/2[33]等細胞死亡調節基因的表達，進而導致抗凋亡和促存活基因的上調并抑制對化療藥物的細胞凋亡反應。Th1/Th2的平衡也影響著乳腺癌的發生、發展。人體內幼稚的CD4+T細胞在不同的細胞因子刺激下分化為Th1和Th2細胞。Th1細胞產生具有抗癌活性的IFN-γ和IL-12。Th2細胞介導體液免疫并產生IL-4和IL-10。與Th1反應相反，Th2反應實際上可以促進癌癥進程。Th2細胞分泌的IL-4抑制IFN-γ的分泌[34]。IL-10可下調Th1細胞的增殖和IFN-γ的產生[35]。因此Th1/Th2的平衡會極大影響著乳腺癌患者腫瘤微環境的狀況。PD-1/PD-L1是一個免疫檢查點，被認為是乳腺癌治療的重要靶點。PD-1與其配體PD-L1結合，限制T細胞活性，并使乳腺癌細胞逃避免疫系統[36-37]。在正常情況下，PD-1/PD-L1通路通過負性調節免疫反應防止過度刺激并維持對自身抗原的免疫耐受，相反，PD-L1在乳腺癌中過度表達，導致腫瘤微環境中的免疫功能受損[38]。故這些可溶性因子可能是通過影響這些信號傳導干預乳腺癌的進展。

綜上所述，可溶性因子影響著乳腺癌患者的預后。它們自身對乳腺癌組織的損傷及相互間的聯級反應是推動乳腺癌發展的重要因素。此外，可溶性因子還對乳腺癌患者的免疫功能產生了重要影響。但這此過程所涉及的因素較多，故本文并未詳細論述。后續課題組將在嚴格執行倫理審查的情況下廣泛收集臨床病理樣本，通過單細胞測序等手段明確可溶性因子對各類免疫細胞的影響。