基于TCGA 数据库的肺腺癌自噬基因预后模型建立与应用

陈彬彬 叶威 黄约诺 方涛 林雪娇 董晓亭

肺腺癌（lung adenocarcinoma，LUAD）是肺癌中常见的病理类型，建立准确预测LUAD 患者预后的工具有利于指导临床诊治[1]。自噬（autophagy）是溶酶体将细胞内受损、变性或衰老的蛋白质以及细胞器进行消化降解的过程[2]。自噬在肿瘤发生、发展中发挥着双刃剑的作用：在正常情况下，细胞自噬防止致癌的蛋白质和细胞器等积累，维持细胞的稳定；在肿瘤细胞中，细胞自噬可以降解非功能细胞器提供营养给癌细胞促进生长[2]。ATGs 的下调能够通过抑制Hedgehog 信号通路增强抗LUAD 的疗效[3]。本研究利用出公共数据库筛选出生存相关的自噬基因（autophagy gene，ATG），建立LUAD-ATG 预后模型并加以验证其应用价值。

1 材料与方法

1.1 数据下载整理 从癌症基因数据库（The cancer Genome Atlas，TCGA）中选择LUAD 患者的RNA 表达数据，提取与患者相关的临床数据。从人类自噬数据库中（Human Autophagy Database http：//www.autophagy.lu）查找自噬相关基因。本文的所有下载数据截至2021 年1 月23 日。

1.2 自噬差异基因表达 从TCGA 数据库中筛选出LUAD 组织和正常组织样本中ATGs 的表达数据。通过R（v4.0.3）软件运用Wilcox 检验筛选LUAD 癌组织样品和正常样品中差异表达的ATGs。筛选标准为：｜logFC｜＞1，FDR＜0.05。为了进一步探索潜在生物学功能，调用中“org.Hs.eg.db”“enrichplot”和“gplot2”包进行GO 功能注释分析，P＜0.05 被认为是有意义的生物学通路。

1.3 预后相关的自噬基因 采用Survival 包先后运用单、多因素COX 回归筛选出影响患者生存期的ATGs（P 值＜0.05）。此外，运用多因素Cox 回归计算出预后相关ATGs 的相关系数（coefficient，coef）并计算每个样本的风险值（risk score，RS），从而建立预后模型。

RS=a1*coef1+a2*coef2......an*coefn

a 代表预后相关TAGs 的表达量，coef 代表预后相关TAGs 的相关系数。

然后根据RS 的中位数将所有患样本分为高风险组和低风险组，运用Kaplan-Meier 生存曲线评价高、低风险组的生存差异，P＜0.05 被认为生存差异有统计学意义。进一步判断独立预后因素，通过Survival 包对年龄、性别、TNM 分期、临床分期和RS 多个因素进行单因素COX 分析。TMN 分期能够对肿瘤进行准确的分期，指导临床肿瘤治疗，其中T 是指照肿瘤大小，N 是指淋巴状态，M 是指远处转移状态。为了排除混杂因素的干扰，对上述多个影响因素进行多因素COX 分析。最后采用受试者工作特征（receiver operator characteristic，ROC）曲线评估预测模型的准确性。其中曲线下面积（area under the curve，AUC），AUC＞0.06 被认为是可以接受的预测。

1.4 临床相关性分析 运用R 软件中beeswarm 包绘采用Wilcox 检验来比较预后相关的ATGs 与肿瘤大小、淋巴结状态、远处转移、临床分期、性别、年龄之间的相关性。P＜0.05 被认为差异有统计学意义。

2 结果

2.1 自噬相关差异基因的筛选 TCGA 数据库中共下载了535 例LUAD 组织样本和59 例癌旁组织样本的mRNA 表达数据及对应的临床信息。采用Wilcox 检验所得P 值的对数-lg（pvalue）为纵坐标，以差异倍数的对数值log2（fold change，FC）为横坐标绘制火山图。以｜logFC｜＞1，FDR＜0.05 为标准筛选出30 个差异表达ATGs，其中HSPA5、IFNG、P4HB、ATIC、PTK6、ERO1A 等18 个基因上调，DLC1、CCL2、DAPK2、PRKCQ、FOS 等12 个基因下调，进一步通过箱线图展示ATG 在肿瘤组织和癌旁组织的差异表达情况（P＜0.05），见图1。

图1 肺腺癌和癌旁样本的差异表达的自噬基因

为了揭示差异TAGs 的潜在功能，用clusterProfiler 包对差异表达的ATGs 进行功能分析，分别展示了分子功能（molecular function）、生物过程（biological process）、细胞组成（cellular component）三个部分GO分类条目（P＜0.05）。图2 显示，主要与内源性凋亡、未折叠蛋白反应、ErbB-2 信号通路等相关。

图2 差异表达的自噬基因的GO 富集分析

2.2 ATG 预后模型的建立 通过单因素Cox 回归分析以P ＜0.01 为阈值筛选出ERN1、ATG16L2、VEGFA、CCR2、ATG2A、DDIT3、CFLAR、IRGM 等共有28 个可能影响患者总体生存期的ATGs，其中HR＞1表示高风险基因，HR＜1 表示低风险因素，P＜0.05 表示差异有统计学意义（见图3）。进一步采用多因素Cox 分析优化，筛选出RAC1、SQSTM1、CD46、NRG3、IKBKB、VMP1、WIPI1、FKBP1B、IKBKE 在内的9 个预后相关ATGs。其中CD46、NRG3、IKBKB、VMP1、IKBKE、WIPI1 被认为是危险因素，高表达与较差生存期相关；高表达的RAC1、SQSTM1、FKBP1B 被认为是保护因素，过度表达能够延长生存期。见表1。

图3 肺腺癌差异表达的自噬基因单因素COX 回归分析注：Hazard ration（HR）为风险比；pvalue 为P值

表1 9 个自噬基因的多因素COX 回归分析

根据风险值公式计算的各个样本RS 值并将LUAD 患者分为高风险组、低风险组。为了进一步验证结果准确性，图4A 运用Kaplan-Meier 生存曲线进行分析发现低RS 组具有更长的生存期（P＜0.05）。此外，通过将所有的患者跟进RS 进行排名，以观察其生存率。随着RS 越高，死亡患者数量越多（见图4B-C）。图4D 显示，LUAO 自噬相关模型中ATGs 在高风险组、低风险组的表达量。

图4 风险值与生存时间/生存状态、基因表达量之间的关系

2.3 生存分析及ROC 曲线 纳入Cox 回归分析的自变量中的连续变量为RS，年龄、TMN 分期、病理分期等为分类变量。单因素Cox 回归模型分析结果显示，RS、肿瘤原发灶情况是影响癌总体生存期的危险因素（P＜0.05），见图5A。多因素Cox 回归模型分析结果显示，RS、肿瘤原发灶情况是LUAD 患者预后的独立影响因素（P＜0.05），见图5B。

图5 单因素（A）和多因素（B）分析患者预后影响因素

ROC 曲线分析影响因素与5 年的生存率之间的关系，RS 的ROC 曲线的AUC 值最大，为0.659，认为RS 可作为LUAD 独立预后指标。见图6。

图6 LUAD 预后预测模型的ROC 曲线分析

2.4 不同临床特征下自噬相关基因表达水平情况进一步分析9 个预后相关ATGs 与临床特征的相关性（见图7-8），IKBKB 与临床分期、肿瘤原发灶情况、肿瘤淋巴结情况显著相关（P＜0.05）；NRG3 在肿瘤淋巴结情况显著相关（P＜0.05）；FKBPIB 与性别具有相关性（P＜0.05），在女性中FKBPIB 的表达量更多。

图7 CD46、FKBP1B、IKBKB、IKBKE、NRG3 与临床特征之间的相关性

图8 RAC1、SQSTM1、VMP1、WIPI1 与临床特征之间的相关性

3 讨论

自噬既是一种促细胞生存机制，也是一种程序性细胞死亡机制。越来越多证据证明自噬不仅能够为癌细胞提供能量，而且能够从各个途径促进癌症侵袭及转移[4]。Wang等[5]对393 例接受了放射治疗的北美非小细胞肺癌患者进行研究，通过对ATG2B、ATG10、ATG12 ATG16L2 在内的4 个ATGs 中的9个潜在功能单核苷酸多态性进行基因分型，发现ATG16L2 rs10898880 CC 变异基因型患者的无局部复发生存率、无局部区域进展生存率和总生存率较好，表明该位点的多态性参与了抗肿瘤治疗和肺部炎症的过程，也提示自噬以某种方式参加肿瘤的病理变化。NRG3 是神经调节蛋白家族之一（其它成员包括NRG1、NRG2 等），是ErbB 家族受体酪氨酸激酶的配体[6]。近年研究显示，NRG-1 与非小细胞肺癌、食管癌、甲状腺癌等多种肿瘤的发生、发展密切相关，其作用过程可能涉及肿瘤相关的信号通路分子改变[7-9]。但是关于NRG3 尚未针对癌症进行过研究。IKBKB 是组成IKK 复合体的重要催化亚基之一，其在IKK 复合体激活NF-κB 过程中发挥重要作用[10]。而NF-κB 异常激活存在肿瘤细胞和肿瘤微环境中，被认为在乳腺癌、肝癌、肺癌、前列腺癌在内的多种肿瘤发生、发展发挥着重要的作用[11-14]。此外，IKBKB的高表达与乳腺癌耐顺铂耐药性相关，可能通过抑制FOXO3 而导致耐药性的产生[15]。Qi等[16]通过A549细胞中转染pcDNA3.1/IKBKB 质粒过表达了IKBKB基因，发现IKBKB 高表达的A549 细胞耐药性明显增加，明显抑制细胞凋亡。这可能与激活NF-κB 信号通路有关，导致耐药性的产生。

对差异表达的ATGs 进行GO 分析发现，较多的基因富集在凋亡信号通路、未折叠蛋白反应，这与以前的研究结论是一致的：肿瘤细胞通过引起的错误折叠的蛋白质和受损的细胞器来消除细胞器的压力[17]。ErbB-2 信号通路在肿瘤发生中发挥着重要的作用，ERBB2 种族突变存在于肺癌患者中[18]。本研究发现，ErbB-2 信号通路上富集到多个ATG，侧面说明ATGs 也可通过这种经典的癌症通路影响肿瘤进展。

进一步使用单、多COX 回归分析后确定了9 个与生存相关的ATGs（RAC1、SQSTM1、CD46、NRG3、IKBKB、VMP1、WIPI1、FKBP1B、IKBKE）并构建预后风险模型。根据风险评分中位值分为高、低风险组，不同的RS 分组与生存相关，拥有不同临床特征（P＜0.05）。通过绘制RS 分布、Kaplan-Meier 分析证实，相较于高RS 组，低RS 组患者的预后更好。ROC曲线具有简单、直观的特点，用来评价模型的准确性。根据ROC 曲线，9 个ATGs 建立的ATGs 预后模型AUC 为0.659，进一步证实预后模型的稳定性和可靠性。此外，对RS 和其他因素进行了单、多因素Cox 回归分析，证实RS 可作为LUAD 患者的独立预后预测因子。

近年出现了许多不同的统计学方法构建的肿瘤预后模型，包括乳腺癌、肺癌、结直肠癌等[19-21]，有利于临床诊治，推进精准医学发展。本研究构建的ATG-LUAD 预后模型仍存下述局限：（1）本模型处于初步建立阶段，仍需要更多临床数据验证模型的可靠性；（2）本模型以TCGA 数据库为基础构建，缺少亚洲人群的转录组数据，在预测国人生存期方面的准确性较差。