基于生物信息学探索肺腺癌预后相关基因及免疫浸润分析①

黄宏玉 巫艳彬 覃寿明（广西医科大学第一附属医院呼吸与危重症医学科，南宁530021）

肺癌是目前全球范围内第二常见癌症，也是全球癌症相关死亡的主要原因。非小细胞肺癌（nonsmall cell lung cancer，NSCLC）约占肺癌的85%［1-2］，而肺腺癌（lung adenocarcinoma，LUAD）在NSCLC中的占比为40%［3-4］。目前，肺癌主要以手术、放化疗、分子靶向治疗、免疫治疗［5］、中医药治疗等多学科联合的综合治疗模式为主，这极大地改善了肺癌患者的预后。但由于个体差异和耐药性，其临床疗效有限。因此，急切需要寻找肺癌的致病机制以及更高效、精准的诊治方法［6］。

越来越多的证据表明肿瘤微环境（tumor immune microenvironment，TME）在肿瘤发生发展中的重要性。肿瘤的进展主要归因于一个独特而适宜的TME。TME主要由多种基质细胞维持，包括肿瘤干细胞、癌内皮细胞、癌相关成纤维细胞和浸润免疫细胞［7-8］。癌细胞与其支持细胞之间的协同作用导致了肿瘤的恶性表型出现，如永生增殖、抗凋亡、逃避免疫监视等［9］。因此，肿瘤浸润性免疫细胞（TIC）在TME中可作为治疗效果的指标［10］。

本研究从TCGA［11］（https：∕∕cancergenome.nih.gov∕）数据库中下载594例转录组RNA序列数据，用CIBERSORT和ESTIMATE计 算LUAD样本 中TICs的比例以及免疫成分评分（Immune score）和基质成分 评 分（Stromal score）［12］。利 用Immune score和Stromal score共享的DEGs构建PPI网络并进行单变量Cox回归分析，然后利用PPI网络中的核心节点和单变量Cox回归分析得到的最高显著性因子进行交叉 分 析，获 得CCR2、CD79A、BTK、CD19、CD22、CD28、PTPRC共7个预后关键基因。我们选择其中一个基因（CD22）进行生存率和临床病理特征相关性分析、Cox回归分析、基因集富集分析（GSEA）以及与TICs的相关性分析。CIBERSORT算法对LUAD患者进行免疫细胞浸润分析，初步探究免疫浸润细胞在LUAD中的作用。本研究结果为LUAD提供预后相关因子，并初步评估LUAD中免疫细胞浸润情况，为进一步研究LUAD的治疗和分子机制提供了新思路。

1 资料与方法

1.1 资料594例LUAD患者（正常组59例、肿瘤组535例）的转录组RNA-seq数据和相应的临床数据均从TCGA数据库下载。

1.2 方法

1.2.1 免疫评分、基质评分和综合评分的生成使用R语言（4.0.0版本）装载的ESTIMATE包，估计每个样本的免疫、基质成分及两者之和在TME中的比例，以3种分数形式表现：Immune Score、Stromal Score和Estimate Score，分别表示免疫组分、基质组分和两者组分之和的评分，即各自评分越高，对应组分在TME中的比例越大。

1.2.2 生存分析 生存分析采用R语言装载的survival和survminer的软件包。其中535例肿瘤样本中的513例有详细的生存时间记录可用于生存分析，时间跨度从0～19.8年。Kaplan-Meier法绘制生存曲线，log-rank作为统计学显著性检验；P＜0.05则认为差异有统计学意义。

1.2.3 Immune Score和Stromal Score的高分和低分组之间的DEGs生成 将535个肿瘤样本的Immune Score和Stromal Score分别与各自的中位分数比较，分为高分组和低分组。使用R语言的limma包对基因表达进行分析，并通过高分组与低分组比较产生DEGs。将log2（高分组∕低分组）转化后，折叠变化（fold change）＞1且错误发现率（FDR）＜0.05的DEGs视为差异有统计学意义。

1.2.4 功能和通路富集分析 利用R语言，借助于cluster Profiler、enrich plot和ggplot2包，对363个DEGs进行GO和KEGG富集分析。只有P值和q值均小于0.05才可以认为富集显著。

1.2.5 热图 用R语言编制DEGs的热图，并用pheatmap包进行分析。

1.2.6 不同临床分期的评分差异分析 从TCGA下载与LUAD标本相对应的临床病理特征数据。采用R语言对不同临床分期进行统计分析，以Wilcoxon秩和检验或Kruskal-Wallis秩和检验为显著性检验。

1.2.7 PPI网络建设 用STRING数据库构建PPI网络，然后用3.8.0版的Cytoscape进行重建。采用交互关系置信度大于0.95的节点构建网络。

1.2.8 Cox回归分析 用R语言进行单变量Cox回归分析，并列出符合P＜0.05的基因。

1.2.9 基因集富集分析 从分子签名数据库（Molecular Signatures Database）下载c2.cp.kegg.v7.1.symbols.gmt作为GSEA的目标集合，使用从Broad Institute下载的软件GSEA-4.0.3进行GSEA。用所有肿瘤标本的全部转录组进行GSEA，只有NOMP＜0.05和FDRq＜0.05的基因集被认为是有意义的。

1.2.10 TICs丰度分布 应用CIBERSORT法估计所有肿瘤样本的TICs丰度分布，然后进行质量过滤，仅选择P＜0.05的535个肿瘤样本进行后续分析。

2 结果

2.1 评分与LUAD患者的生存相关 为了建立免疫和基质与存活率的相关性，分别对LUAD患者的Immune Score、Stromal Score和Estimate Score与相应中位数的比较确定高低分组后进行Kaplan-Meier生存分析。Immune score或Stromal Score的评分越高，说明TME中的免疫或基质成分越多。Estimate score是免疫得分和基质得分的总和，表示两种成分在TME中的综合比例。如图1所示，免疫组分、基质组分和综合得分的比例都与总生存率呈正相关（图1A～C），提示TME中的免疫组分和基质组分均可作为判断LUAD患者预后的指标。

图1 评分与LUAD患者生存率的相关性Fig.1 Correlation between score and survival rate of pa⁃tients with LUAD

2.2 LUAD患者临床病理分期与评分的相关性研究 为了确定免疫组分和基质组分的比例与临床病理特征的关系，从TCGA数据库中分析了LUAD患者的相应临床信息。免疫组分评分与临床病理特征的关系如图2A～D所示,Stromal Score与临床分期和TMN分期的M分期呈负相关（图2A,P=0.015;图2C,P＜0.01）;基质组分评分与临床病理特征的关系如图2E～H所示,Immune Score与TNM分期的T分期呈负相关(图2F,P＜0.01);图2I～L则显示了综合评分与临床病理特征的关系,可见Estimate Score与临床分期和TMN中T和M分类的进展呈负相关(图2I,P=0.022；图2J,P=0.027;图2K,P=0.022)。这些结果提示免疫组分和基质组分的比例与LUAD的进展有关，如肿瘤浸润和转移。

图2 免疫评分和基质评分与临床病理分期特征的相关性Fig.2 Correlation between clinicopathological staging characteristics with immune score and matrix score

2.3 Immune Score和Stromal Score共 享 的DEGs相关基因的富集 为确定TME中免疫组分和基质组分基因谱的变化，对高、低分样本进行了比较分析。与来自Immune Score（得分高的样本与得分低的样本）的776个DEGs相比，其中613个基因上调，163个基因下调（图3B～D）。同样，从Stromal Score中获得792个DEGs，由678个上调基因和114个下调基因组成（图3A、C、D）。Venn图的交叉分析显示，Immune Score和Stromal Score共有297个高表达基因共享，66个低表达基因共享。这些DEGs（共363个基因）可能是TME临床结局的决定因素。基因本体（GO）富集分析结果表明，DEGs与白细胞增殖、淋巴细胞增殖、单个核细胞增殖、淋巴细胞活化的调节和T细胞活化等免疫相关GO项对应（图3E）；京都基因和基因组百科全书（KEGG）的富集分析显示DEGs显著富集在B细胞受体信号通路、趋化因子信号通路、细胞因子-细胞因子受体相互作用通路、造血细胞谱系通路以及病毒蛋白与细胞因子和细胞因子受体的相互作用抗病毒相关通路（Viral protein interaction with cytokine and cytokine receptor)上(图3F）。因此，DEGs的整体功能与免疫相关活动有关，提示免疫组分的参与是LUAD中TME的主要特征。

图3 热图、韦恩图和GO和KEGG的富集分析（适用于DEGs）Fig.3 Heat map，Wayne map and enrichment analysis of GO and KEGG（applicable to DEGs）

2.4 PPI网络与单变量COX回归 利用美国国立普通医学科学研究所（NIGMS）的Cytoscape软件构建了基于STRING数据库的PPI网络。363个基因之间的相互作用如图4A所示，条形图表示按节点数排列的前30个基因（图4B）。对LUAD患者的生存率进行单变量COX回归分析，以确定363个DEGs中的与预后相关的关键基因（图4C）。然后，对PPI网络中的领先节点与按单变量COX回归P值排序的前30个因子进行交叉分析，得到CCR2、CD79A、BTK、CD19、CD22、CD28、PTPRC共7个因子在上述各分析结果中重叠。（图4D）。

图4 蛋白质-蛋白质相互作用网络和单变量COXFig.4 Protein-protein interaction network and univariate COX

2.5 CD22表达与LUAD患者生存期及TNM分期的关系CD22高表达组患者比CD22低表达的患者存活时间更长（图5C）。结合临床特点进行CD22分析，Wilcoxon秩和检验显示正常标本中CD22的表达明显高于肿瘤标本（图5A）。在来自同一患者的正常组织和肿瘤组织的配对分析中也观察到相同的结果（图5B）。以上结果表明，CD22在TME中的表达与LUAD患者的预后呈正相关。另外CD22表达随着临床病理分期和TNM分期中T分期的进展而下降（图5D～G）。

图5 样本中CD22的分化表达及其与LUAD患者生存率和临床病理分期特征的相关性Fig.5 Differentiation and expression of CD22 in sample and its correlation with the survival rate and clini⁃copathological staging characteristics of LUAD pa⁃tients

2.6 CD22在TMN中的作用 分别在CD22高表达组和低表达组使用GSEA软件，并与表达水平的中位数进行比较。如图6A所示，CD22高表达组中主要富集有免疫相关活性通路，如B细胞受体、T细胞受体、趋化因子、细胞黏附分子、JAK-STAT信号通路；对于CD22低表达组，代谢基因通路富集，包括蛋白酶体、碱基切除修复、RNA聚合酶、柠檬酸循环、氧化磷酸化和磷酸戊糖途径。

图6 CD22高表达和低表达样品的GSEAFig.6 GSEA of samples with high expression of CD22 and low expression of GSEA

2.7 CD22与TICs比例的相关性研究 为了进一步证实CD22表达与免疫微环境的相关性，采用CIBERSORT算法分析了肿瘤浸润免疫亚群的比例，构建了21种LUAD样本的免疫细胞图谱（图7A、B）。差异及相关分析结果显示，共有12种TICs与CD22的表达相关（图8A）。其中7种TICs与CD22表达呈正相关，分别为记忆B细胞、幼稚B细胞、静息肥大细胞、静息突状细胞、静息的CD4+记忆T细胞、CD8+T细胞、调节T细胞；5种TICs与CD22表达呈负相关，分别为活化的树突状细胞、巨噬细胞M0、巨噬细胞M2、激活的肥大细胞、中性粒细胞、活化的NK细胞。这些结果进一步证实CD22水平影响TME的免疫活性，CD22可能是TME的免疫状态的潜在标志物。

图7 肿瘤标本中的TICs丰度及相关分析Fig.7 TIC profile and correlation analysis in tumor speci⁃mens

图8 TICs比例与CD22表达的相关性Fig.8 Correlation between TICs proportion and CD22 expression

3 讨论

本研究的目标是从TCGA数据库中找出影响LUAD患者生存率和TNM分期的TME相关基因。结果证实CD22与肿瘤免疫微环境中的免疫活性有关，一系列生物信息学分析表明CD22可能是LUAD患者TME状态的指标。大量研究揭示了免疫微环境在肿瘤发生中的重要性，探索TME重塑的潜在治疗靶点具有重要意义。本研究对TCGA数据库中LUAD数据的转录组分析表明，TME中的免疫成分和基质成分影响患者的预后，TME中免疫成分和基质成分的比例与LUAD的进展（临床分期和浸润）显著相关。这些结果提示探讨肿瘤细胞与免疫细胞的相互作用，为寻找更有效的LUAD治疗方案提供了新观点。

近几年，免疫治疗在肿瘤临床治疗上取得了很大的进展，免疫检查点抑制剂（ICIS）已被批准用于晚期NSCLC患者的一线药物［13］。虽然ICIS对治疗非小细胞肺癌显示出良好的疗效，但ICIS与免疫相关不良事件也不容忽视［14-15］。有研究报道程序性死亡配体1（PD-L1）不是决定NSCLC患者是否适合接受免疫治疗的最佳指标，因此有必要找出一些新的非小细胞肺癌免疫治疗候选方案［16］。此外，TICs的丰度与NSCLC的五年生存率显著相关［17］。本研究对TCGA数据库中LUAD的转录组分析，发现CD22的表达减少与晚期临床分期、肿瘤浸润和预后不良密切相关。因此，CD22可能是LUAD患者TME潜在的预后因素和治疗靶点。

CD22为Ⅰ型跨膜糖蛋白，相对分子质量为140 kD，是唾液酸结合免疫球蛋白样凝集素家族成员［18-19］。CD22的细胞外部分包含7个Ig结构域，其末端的V状结构域能与包含α-2，6连接唾液酸的糖蛋白配体特异性结合［20］。CD22是近年来发现的一种新兴抗原，是B细胞表面抑制性辅助受体之一，其与配体的不同结合形式以及在细胞表面的不同分布状态均影响B细胞活性，与B细胞的发展、分化和功能密切相关［21］。CD22在细胞内通过下游通路减弱BCR信号发挥负向调节作用，与骨髓内皮细胞表面配体作用促进B细胞归巢。CD22对B细胞系的限制以及在抗体介导的交联后迅速内化，使其成为干预B细胞免疫及治疗某些B细胞肿瘤的良好靶标［22-25］。研究发现，CD22在大多数B细胞恶性肿瘤细胞表面都有表达，包括急性淋巴细胞白血病、非霍奇金淋巴瘤和慢性淋巴细胞白血病。而CD22在非小细胞肺癌中的作用仍存在争议。有研究表明CD22抗原广泛表达于肺癌细胞，是抗体治疗的靶点，检测CD22的作用可能有助于更好地了解肺癌的发病机制和肿瘤侵袭性［18］。而另外一项研究也说明了CD22在肺癌细胞表面的表达水平不高，而且这些肺癌细胞不能被抗CD22免疫毒素杀死［26］。因此，CD22在非小细胞肺癌中的表达、作用机制及靶点治疗方面仍需探索。

CD22已成为免疫疗法的一个新兴靶点。目前，以CD22为靶点的治疗方法主要有单克隆抗体药物、抗体偶联物（ADC）、CAR-T疗法三大类。CD22已成为除CD19之外CAR-T疗法研究领域的又一火热靶点。CAR-T的概念自提出以来，在恶性血液系统肿瘤中获得较好效果［27］。临床上CAR-T疗法的3个靶分子（CD19、CD20、CD22）均是针对B细胞［28］，而目前所知绝大多数实体瘤如非小细胞肺癌、黑色素瘤中不具备理想的靶向抗原，仍需要更长时间的探究。因此，我们进一步分析了CD22表达与LUAD中TME的关系。GSEA结果显示，CD22高表达组显著富集于免疫相关信号通路，如B细胞受体、T细胞受体、趋化因子、细胞黏附分子、JAK-STAT信号通路。CD22低表达组的代谢途径包括蛋白酶体、碱基切除修复、RNA聚合酶、柠檬酸循环、氧化磷酸化和磷酸戊糖途径。这些结果提示CD22可能是肿瘤微环境中免疫和代谢状态的指标。因此，CD22表达随LUAD的进展而下调，TME由免疫状态转为向代谢状态，表明CD22可能在LUAD中发挥抗肿瘤作用。

CD22在B淋巴细胞的功能至关重要。一些文章表明肿瘤中的B细胞浸润可能是NSCLC患者生存的一个良好预后生物标志物。B细胞可能参与了LUAD患者的抗肿瘤免疫反应［17，29］。肿瘤内B细胞被认为在肿瘤微环境中发挥重要作用，特异性B细胞受体（BCR）和抗体在抗肿瘤和潜在的促肿瘤反应中的作用越来越明显。血清中肿瘤相关抗原自身抗体已被证实为早期发现癌症的生物标志物［30-32］，也可作为疾病晚期有用的预后标志物［33］。肿瘤中B细胞数量也与多种实体瘤的良好预后相关［34-35］。研究显示肿瘤浸润性B细胞的数量与具有特定突变驱动因子的LUAD患者的生存率之间存在相关性。本研究中CIBERSORT分析显示，LUAD患者的记忆B细胞与CD22表达呈正相关，提示CD22可能是TME维持免疫活性状态的重要因素。

总之，利用ESTIMATE算法，通过对TCGA数据库中LUAD样本的功能富集分析，确定LUAD中TME预后相关基因有CCR2、CD79A、BTK、CD19、CD22、CD28、PTPRC。挑选CD22进行单基因分析，发现CD22上下调可能会影响TME的免疫和代谢状态。这些结果为制定LUAD的预后生物标志物及进一步研究其分子机制提供了理论依据。未来我们会在分子、细胞及组织水平上验证本次研究结果的准确性，并进一步探讨CD22与LUAD肿瘤微环境中B细胞之间的调控关系。