基于特异性增强子的低级别胶质瘤分型研究

严光灿 田 伟 刘美娜

哈尔滨医科大学公共卫生学院卫生统计教研室，150081 黑龙江 哈尔滨

脑胶质瘤是神经系统最常见的恶性肿瘤，按恶性程度可将其分为I～IV级[1]。由于患者异质性高，预后差异大，2016年新的分类标准增加了形态学和分子特征；2021年进一步增加了甲基化等特征[2-3]。临床上II～III期的胶质瘤患者被统称为低级别胶质瘤(lower grade glioma，LGG)，相对于I期和IV期胶质瘤，LGG患者的肿瘤异质性更高，如肿瘤微环境中的细胞成分及含量不同、预后差异大等[4-5]；即使病理类型相同的患者在接受相同治疗后，治疗结局仍有很大差异[6]，精准的分型有助于临床医生选择更合理的治疗方案。免疫相关基因(immune related genes, IRGs)异常与免疫浸润、肿瘤发生发展、治疗抵抗等过程息息相关，增强子对IRGs的表达有重要的调控作用。增强子是DNA序列上重要的远端调控元件，其转录本称为增强子RNA(enhancr RNA，eRNA)，eRNA含量可代表增强子的功能活性。增强子在各类肿瘤/正常组织中广泛表达，而不同组织、癌症间，其表达有明显的特异性[7]。因此，研究特异性增强子调控的靶基因，有助于LGG患者精细、合理地分型；分析各型患者间的特征差异，对促进精准医疗实施、改善患者预后及相关机制研究具有重要意义。

1 资料与方法

1.1 数据获取与样本筛选

从癌症基因组图谱数据库(the cancer genome atlas, TCGA；www.xena.ucsc.edu/)获取LGG患者的基因表达数据、临床表型数据及预后生存数据作为训练数据集；从中国胶质瘤基因组图谱数据库(chinese glioma genome atlas, CGGA；www.cgga.org.cn/)获取相同类型数据作为验证集。由于CGGA数据包括2个批次的数据，因此需要进行批次校正。增强子与靶基因的调控关系数据从eRNA数据库(eRic；www.hanlab.uth.edu/eRic/)中获取。并从ImmPort 数据库(www.immport.org)和InnateDB数据库(www.innatedb.com)中获取IRGs，从MSigDB数据库中获取免疫相关通路基因集。样本纳入标准：(1)原发性肿瘤；(2)肿瘤分期为II～III期；(3)样本含有基因表达数据、临床表型数据及预后生存数据。

1.2 LGG分型分析

获取LGG特异性eRNA，确定其调控的IRGs；利用IRGs对LGG样本进行K-means聚类分析；聚类数设置为2～8，分析聚类热图的相关性，联合手肘法，选取拐点处的值为最佳聚类数，确定LGG分型；利用CGGA数据进行相同分析，验证聚类结果的稳健性。聚类分析通过R包“ConsensusClusterPlus”实现。

1.3 免疫浸润与通路富集分析

肿瘤微环境中免疫细胞浸润对肿瘤进展、药物敏感性及患者预后有重要影响。利用“CIBERSORT”算法对做LGG患者进行免疫浸润分析，比较各亚型间免疫细胞浸润模式的差异[8]。富集分析探索特异性增强子调控的IRGs所参与的生物学过程，包括KEGG富集分析、GO富集分析。分析结果利用R包“clusterProfiler”实现。

1.4 临床特征分析

生存分析比较各亚型间患者生存时间差异， log-rank检验对差异进行统计推断，并绘制生存曲线。检验分析各型间肿瘤分期、TERTp突变、MGMTp甲基化、IDH突变、染色体1p/19q共缺失的构成差异。秩和检验分析各亚型间基因突变数、TERT基因表达量的差异。

1.5 基因集差异分析

基于免疫相关通路基因集的基因集差异分析(gene set variation analysis，GSVA)：计算每个样本的通路富集得分(得分大于0表示通路活性上调，小于0表示通路活性下调)，比较各型间的通路活性差异，定义|logFC| > 0.5且校正P<0.05的通路为亚型间的差异通路。GSVA利用R包“GSVA”完成。

2 结果

2.1 特异性增强子调控的靶基因筛选

纳入排除标准选择得到训练集495个样本、验证集426个样本。筛选出227个LGG特异性增强子，50个特异性增强子调控的免疫相关性靶基因，包括AKT2、HIF1AN、LRSAM1、PTCH1等。

2.2 LGG分型

TCGA样本聚类发现聚类数为5时效果最好，组内相关性高，组间相关性低(图1A)；手肘法结果显示,拐点处值的为5(图1B)，确定最佳聚类数为5类。验证集CGGA样本聚类分析后得到相同结果，表明聚类结果稳定。LGG分为5型，分别记为A1～A5亚型，其构成比分别为 22.6%、21.6%、13.7%、26.3%和15.8%。

A:K-means聚类热图；B:手肘法结果

2.3 免疫浸润分析

免疫浸润分析结果显示：不同亚型患者间免疫细胞浸润模式不同；22种免疫细胞中有18种在各亚型间存在明显差异(P<0.05)。见图2。其中，A1亚型患者主要富集M2型巨噬细胞、活化型肥大细胞及单核细胞；A2型富集M2型巨噬细胞；A3型富集嗜酸性粒细胞、活化型NK细胞，M2型巨噬细胞含量最低；A4型富集单核细胞和M2型巨噬细胞；A5型富集M1、M2型巨噬细胞，静息型CD4记忆T细胞和CD8+T细胞，单核细胞和NK细胞含量最低。

注:*P<0.05；**P <0.01；***P<0.001。

2.4 富集分析

GO富集分析发现IRGs主要参与c-Jun-氮末端激酶(JNK)活性调节、多种细胞分化(如脂肪细胞、白细胞)、造血功能调节以及多种酶活性调节(如内切酶、蛋白酶)、受体配体活性、蛋白/组蛋白脱乙酰基酶活性等多条信号通路(图3A)。KEGG结果显示IRGs主要参与酒精性肝病、多种病原体感染(如沙门氏菌、志贺杆菌、人免疫缺陷病毒)以及cAMP等多种生物学过程(图3B)。

A:GO富集分析结果；B:KEGG富集分析结果

2.5 临床特征分析

生存分析结果显示A5亚型预后最差，中位生存期约为2年，其余亚型间预后生存无统计学差异。见图4A。A3和A4组疗效最好，超过60%的患者治疗有效；A5组疗效最差，超过30%的患者出现病情进展。见图4B。A5组III期患者比例最高，其次是A2，A3组最低。见图4C。IDH基因突变、MGMTp甲基化均表现为在A5组比例最低，其余各组较高；TERTp突变在A5组最高，A2组次之，最低是A1组；染色体1p/19q共缺失在A2最高，A5最低。见图4D～4G。总突变数在A5组最高，平均突变数约为40，A3组最低，平均突变数不到20；TERT基因表达在A2和A5组最高，A1和A4组最低。见图4H、图4I。

A：生存分析结果；B：疗效构成比；C：分期构成；D：1p/19q共缺失构成；E：IDH突变构成；F： MGMTp甲基化构成；

2.6 基因集差异分析

GSVA结果发现115条差异通路，选择|logFC|最大的前20条通路进行可视化后发现，这些通路在A5组患者中均过度激活，在A3组患者中活性则均受到抑制，其余3组活性无明显差别。见图5。

3 讨论

随着医疗技术发展，许多疾病的疗效均显著提升，但由于LGG患者肿瘤异质性高、个体差异大，疗效提升有限。对LGG患者进行精准的分型，针对不同亚型患者选择不同的治疗方案，有利于提升治疗效果、改善患者预后。以往的LGG分型研究由于缺乏特异性的标志物，分型效果不理想。增强子有高特异性，且可直接影响肿瘤的发生发展，对LGG患者精准分型有重要意义。本研究通过筛选特异性增强子调控的IRGs获得5个LGG亚型，各型间生存、临床特征、基因通路活性及免疫浸润模式等存在显著差异。但增强子对基因的调控作用能多大程度影响患者预后仍需进一步探索，后续研究需在本研究的基础上深入探索，并利用实验验证机制假设。

在肿瘤微环境中，除肿瘤细胞外，还存在大量的免疫细胞、基质等非肿瘤细胞成分。免疫细胞浸润模式的不同对肿瘤发生进展、药物敏感性及患者预后有重要影响。M2型巨噬细胞可促进肿瘤细胞免疫逃逸、加快肿瘤进展和转移[9]；本研究显示A5组M2型巨噬细胞浸润程度最高，这表明A5的不良预后可能与M2型巨噬细胞引起的免疫逃逸有关。单核细胞可分化为巨噬细胞，主要包括M0、M1、M2 3个亚型，其中M0对免疫逃逸有抑制作用，M1和M2则可促进肿瘤细胞的免疫逃逸[10]。预后最差的A5组单核细胞含量最低，说明单核细胞可能更多地分化为抑制免疫逃逸的M0型巨噬细胞；免疫细胞的分化可能对LGG患者的生存结局有重要影响。本研究发现免疫细胞成分和含量不同，会导致患者间预后差异很大，这与各免疫细胞的功能不同有关，但具体的机制有待后续的实验加以佐证。

图5 GSVA差异活性通路

GSVA结果显示各型间多条信号通路活性有显著差异，表明各组间预后及疗效差异可能与这些通路活性不同有关。c-Jun-氮末端激酶(JNKs)是一类蛋白激酶家族，在细胞增殖分化、衰老凋亡等生物学过程中发挥重要作用[11]。本研究筛选的靶基因富集在多条与JNKs活性调节相关的通路，这提示JNKs相关通路可能是影响LGG患者疗效与生存的重要通路；A5亚型预后差可能与JNKs相关通路失调导致细胞增殖分化紊乱、癌细胞异常增殖有关。此外，这些靶基因还参与多种细胞分化、造血功能调节等生物学过程，均是与癌细胞增殖、癌组织生长相关的关键通路。

各型间临床分期存在显著差异，A5亚型III期患者比例最高，与其预后最差相一致；其余亚型间分期也存在明显差异，生存结局却无显著差异，说明除肿瘤分期外，还有其他因素影响LGG患者生存，也说明进行LGG分型十分必要。结果显示：TERTp突变比例高，IDH突变、1p/19q联合缺失、MGMTp甲基化比例低的患者预后差，这与以往的研究结果相一致[2,12-13]。A2组III期患者占比及TERTp突变比例都仅次于A5，但生存却与其余各组无差异，很可能与A2组患者1p/19q共缺失比例最高有关。

综上，本研究获得多个特异性增强子和其参与调控IRGs；共获得5个LGG亚型且分型结果稳定，可为临床治疗提供参考；不同亚型患者间预后生存、免疫浸润模式、通路活性、临床及分子特征存在统计学差异，可为相关机制研究提供参考。