掃碼聯系客服
代謝重編程是腫瘤細胞的主要特征。由腫瘤微環境惡化驅動的葡萄糖、脂質和氨基酸代謝重編程與腫瘤的發生發展存在密切聯系，近年來得到了越來越廣泛的關注。今天小編為大家解析的這篇文章就是以肺腺癌（lung adenocarcinoma, LUAD）轉移相關的代謝重編程為側重點，基于生信信息學方法對公共數據庫中的數據進行挖掘，進一步揭示其潛在的細胞間的變化和可能的致病機理。這篇文章于今年3月發表于Communications Biology(IF = 6.268)。
一、摘要
本研究的主要目的是探究與肺腺癌（lung adenocarcinoma, LUAD）轉移相關的代謝重編程，揭示其潛在的細胞間的變化。基于TCGA數據庫中394例LUAD患者的基因表達譜數據，作者鑒定出11個與轉移相關的代謝基因，這些基因涉及在糖酵解和脂質代謝途徑中。并且，研究通過無監督聚類算法定義了3個代謝重編程表型：MP-I、MP-II、MP-III。在TCGA-LUAD隊列和6個GEO表達譜共計1235個樣本中，MP-III型樣本的糖解酶水平最高、脂質代謝水平最低，并表現出較高的轉移能力和較低的生存率。分析結果表明，TP53和KEAP1基因的突變以及SETD2和PBRM1的缺失可能是導致MP-III型患者代謝重編程的主要因素。LUAD樣本單細胞轉錄組數據進一步驗證從正常到MP-II、MP-III，再到MP-I的代謝軌跡。細胞間通訊分析結果發現，在ANGPTL通路中MP-III型細胞群與內皮細胞和成纖維細胞具有獨特的相互作用，而在VEGF通路中其與內皮細胞的相互作用更強。鑒于此，研究認為，糖酵解-脂質失衡模式能夠導致與轉移相關的代謝重編程表型。該研究為進一步了解致癌驅動因子和腫瘤微環境之間的相互作用提供參考，對于LUAD轉移的治療提供策略。
二、材料與方法
1. 數據來源
本研究中，作者從TCGA、GEO和CPTAC數據庫共獲取10個LUAD隊列。兩個轉移性LUAD數據集（TCGA-LUAD和GSE11969）的患者類型包括遠端轉移、淋巴結轉移和非轉移性。其中，TCGA-LUAD隊列共包括394例LUAD患者樣本（25例遠端轉移、160例淋巴結轉移以及209例非轉移性樣本）。此外，對于TCGA-LUAD和GSE11969這兩個數據集，作者分別下載了390和385個突變和拷貝數變異（CNV）數據。GSE11969數據集共包括18個遠端轉移、15個淋巴結轉移和35個非轉移性LUAD樣本的mRNA表達譜數據。用于生存分析的5個LUAD表達譜，納入的標準為：手術切除后未接受輔助治療的I期患者樣本數要超過50例。
此外，作者選擇了2個LUAD單細胞表達譜數據：GSE131907和GSE123902。在GSE131907數據集中，作者從11位tLung患者、10位腦轉移患者和11個正常肺組織中共提取了208506個細胞。在GSE123902數據集中，作者從7位tLung患者的原發性組織中提取了18511個細胞。
2. 數據預處理
TCGA-LUAD數據集，測序平臺為Illumina HiSeq 2000 platform，fpkm表達矩陣經log2(x+1)轉換；對于單細胞表達譜（scRNA-seq）數據，測序平臺為Illumina HiSeq 2500。將基因的Ensembl ID轉換成symbol。芯片數據GSE31210、GSE50081和GSE68465，測序平臺為Affymetrix；GSE42127數據集測序平臺為Illumina、GSE13213測序平臺為Agilent。體細胞突變數據從Illumina Genome Analyzer DNA Sequencing GAIIx平臺獲取，包含17821個非同義突變。拷貝數變異（CNV）數據使用GISTIC算法進行處理。
3. 差異分析、富集分析
差異基因的鑒定使用student’s t-test，基于超幾何分布模型對差異基因進行富集分析，數據集來源于KEGG數據庫（https://www.genome.jp/kegg/pathway.html）。
4. 代謝重編程表型的鑒定
首先，作者基于代謝相關基因的表達相關性使用層次聚類的方法將患者分成不同的簇（或組）。在聚類之前，去除掉表達量均為0的樣本或細胞，并對表達矩陣進行歸一化處理。接下來，作者通過公式計算MP-score以評估不同簇的患者之間糖酵解和脂代謝基因之間的表達差異。計算公式如下：
5. 相關性及生存分析
使用Wilcoxon秩和檢驗來評估不同簇患者之間缺氧評分、干性評分、增殖評分、免疫評分和腫瘤突變負荷（TMB）的差異。基因間的表達相似性、突變情況、拷貝數變異（CNV）和蛋白表達情況使用Pearson相關系數進行評估。生存時間被定義為從最初的手術切除到死亡或最近一次隨訪的時間。使用Kaplan–Meier法繪制生存曲線。使用單因素cox模型分析患者的生存狀態與臨床特性之間的相關性。
6. 單細胞轉錄組分析
基于前20個主成分，使用R package Seurat中的RunPCA和JackStraw函數來對上皮細胞進行無監督聚類。使用RunUMAP函數對聚類結果進行可視化。使用Featureplot函數對每個細胞中mRNA的平均表達水平進行可視化。使用SingleR對單細胞類型進行注釋。使用R package CellChat分析細胞間通訊。配體-受體互作及其相關的信號通路從CellChat數據庫（https://www.cellchat.org）下載。該方法是根據兩個細胞群間的基因表達、信號因子和細胞百分比來推斷配體-受體間的潛在相關作用強度。
三、結果
3.1 鑒定LUAD轉移相關的代謝重編程表型
首先，基于TCGA數據庫中的表達譜數據，使用Student’s t-test方法進行差異分析（篩選標準為：FDR < 0.05），在轉移和非轉移性LUAD樣本間共鑒定出431個差異基因，其中，顯著下調表達基因151個，顯著上調表達基因280個。分析結果發現，上調表達基因顯著富集到HIF-1信號通路，該通路是糖酵解途徑中的重要通路。先前的研究發現，該通路能夠通過激活糖酵解基因的轉錄來增強糖酵解作用，在癌癥的發生發展過程中起關鍵作用。相比之下，下調表達基因顯著富集到脂代謝相關的生物途徑，包括甘油磷脂代謝（Glycerophospholipid metabolism）和醚脂代謝（Ether lipid metabolism）等。由此作者推測，糖代謝和脂代謝可能在誘導LUAD轉移過程中扮演不同的角色（圖1 a-c）。
因此，作者將顯著富集到糖酵解/糖異生（Glycolysis/Gluconeogenesis）途徑的8個下調基因（ALDOA, ENO1, GAPDH, GPI, LDHA, PGAM1, PGM2, TPI1）以及顯著富集到甘油磷脂代謝（Glycerophospholipid metabolism）和醚脂代謝（Ether lipid metabolism）這兩個脂代謝通路的3個上調基因（PLPP1, GPD1L, PLD3）提取出來，進行更進一步的研究。
首先，作者基于上述11個基因在全部LUAD患者樣本中的表達模式進行無監督聚類。分析結果將全部LUAD患者劃分成3個簇（圖1 d），由于這些基因的表達模式與糖、脂代謝途徑顯著相關，因此作者將這些簇定義為代謝重編程表型（metabolic reprogramming phenotypes, MPs）。接下來，作者計算了每個簇的MP得分，并根據得分的由低到高將簇分別命名為：MP-I、MP-II、MP-III。通過構建混淆矩陣，作者發現，MP-I型樣本的患者類型主要為轉移性，MP-II型樣本的患者類型主要為淋巴結轉移性，而MP-III型樣本的患者類型主要為遠端轉移性（圖1 f）。此外，生存分析結果表明，MP-III患者的生存率顯著低于MP-I和MP-II患者（圖1 g）。
3.2 不同代謝表型的分子特性和臨床性狀差異
分析結果表明，MP-III型患者處于stage IV和stage III的比例顯著高于MP-I和MP-II。此外研究發現，從MP-I ~ III，患者的缺氧評分、干性評分和增殖評分逐漸增加。MP-II型患者的免疫評分最高，而MP-III型患者的免疫評分最低，兩者之間的免疫評分存在顯著差異（p < 0.0001）,見圖2。
3.3 LUAD代謝表型的驗證
上述結果表明，MP-III型大多為未發生任何轉移的早期患者。由此推斷，這些患者可能會有較高的轉移風險，進而導致較差的預后。因此，作者首先基于11個代謝基因的表達譜數據，對來自5個公共數據集的stage I患者進行層次聚類。分析結果發現，stage I患者也可聚類為3種表型（圖3 a-e）。生存分析結果表明，在GSE31210、GSE50081、GSE13213、GSE42127數據集中，MP-III型的stage I患者生存率均顯著低于MP-I和MP-II型。除此之外，作者發現，在全部五個數據集中，MP-I和MP-II型的生存差異并不顯著，但是可以看出在GSE13213、GSE42127和GSE68465這3個數據集中，相較于MP-I型，MP-II型患者的生存率是偏低的，但差異并不顯著。
此外，研究從CPTAC數據庫中獲取了110例LUAD患者的表達譜數據，從蛋白水平上評估11個代謝基因的表達模式。從分析結果中可以看到，基于11個代謝基因的蛋白表達譜，同樣將110例LUAD患者聚成了3個表型。此外，研究觀察到每個基因的mRNA表達水平和蛋白表達水平均呈顯著正相關（圖4）。
3.4 在單細胞水平上描述代謝表型的轉錄軌跡
基于單細胞數據集GSE131907，26436個上皮細胞被分成了14個簇，其中包括正常細胞簇（Cluster 1、11、13）、原代腫瘤細胞簇（Cluster 4、8、9、10）、腦轉移瘤細胞簇（Cluster 0、2、3、5、6、7、12），見圖5 a。圖5 b表現為每個細胞中11個代謝基因的平均表達水平。接下來，作者基于11個基因的表達矩陣對樣本進行聚類以確定上皮腫瘤細胞的代謝表型，分析結果表明，MP-III型患者的腦轉移瘤細胞的比例顯著高于MP-I和MP-II型（圖5 e）。
接下來，作者根據11個代謝基因的mRNA表達變化對上皮細胞進行排序，以構建代謝軌跡。軌跡中的7個轉錄狀態分別代表不同的細胞代謝軌跡。每個細胞狀態中各表型所占比例見圖5 f，可以很明顯的看出，正常上皮細胞顯著富集在state 1和state 4這兩個狀態，然而，處于state 4階段的細胞11個代謝基因的表達水平普遍較低，這一點在TCGA數據集中并沒有觀察到。糖-脂質代謝失衡的代謝軌跡似乎主要是從部分正常上皮細胞開始的，主要在State 1末期分化成MP-I型腫瘤細胞，然后形成一個具有兩個主要細胞命運的分支結構：Cell fate 1和Cell fate 2。通過追蹤Cell fate 1的代謝軌跡，作者發現MP-I型腫瘤細胞逐漸分化成MP-III型（state 3）。此外，對Cell fate 2代謝軌跡的追蹤結果表明，MP-II和MP-III型腫瘤細胞主要富集在State 6和State 7這兩個獨立的分支。
3.5 揭示代謝表型下潛在的細胞間通訊
在該部分分析中，作者從患者原發性腫瘤組織中提取6372個腫瘤上皮細胞、2373個干細胞和35506個免疫細胞。如圖5 a所示，作者在10個細胞群兩兩間共鑒定出41個具顯著意義的配體-受體對。正如所推斷的那樣，作者發現MP-III型細胞群表現出強烈的信號輸出，其次為MP-II和MP-I型細胞群（圖6 b）。
接下來，與3種代謝表型相關的重要配體-受體對被進一步劃分為15個信號通路（圖6 c）。各配體-受體對在信號通路中的相對貢獻率以及在各細胞群中的作用如圖6 b、c所示。此外，分析結果發現，在3個MP細胞群中，MP-II是EGF信號通路的主要受體，可能單獨發揮作用，也可能是與髓細胞、MAST、NK和內皮細胞共同作用（圖7 g-i）。
四、總結
在本研究中，分析結果表明在轉移性LUAD樣本中顯著上調和下調表達的差異基因分別富集在糖酵解和脂代謝途徑中，由此推測糖脂代謝的失衡可能會促進LUAD的轉移。不足之處在于，（1）首先，本研究是基于11個代謝基因的mRNA表達譜來定義代謝表型的，因此，糖酵解和脂質代謝活性需要通過估計其代謝物來驗證；（2）其次，在本研究中發現的細胞間相互作用應使用更多的方法來進行驗證，因為它們是通過配體及其受體的mRNA表達預測的；（3）盡管MP-III群體與其他群體的相互作用在另一個scRNA seq數據集中得到了驗證，但想要獲得更可靠的結果仍需進一步的實驗驗證；（4）此外，其他改變腫瘤微環境的代謝重編程模式，如代謝物等，仍待進一步探究。總的來說，研究鑒定了可能與LUAD轉移相關的三種代謝表型，發現這一代謝重編程背后的四種致癌事件及其與其他細胞間的相互作用（尤其是ANGPTL信號通路）為防止LUAD轉移提供了候選治療策略。
