CN111564177A - 基于dna甲基化的早期非小细胞肺癌复发模型构建方法 - Google Patents

Abstract

本发明涉及生物医学领域，公开了一种基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，用以解决目前非小细胞肺癌复发模型对非小细胞肺癌复发预测不合理的问题。本发明针对甲基化位点，搜索并合并其中的相关性高的甲基化位点，形成相应的甲基化区间；再根据甲基化芯片的类型，挑选出能够被芯片探针检测的甲基化位点/区间作为训练数据建立包含1个甲基化区间和8个甲基化位点的甲基化风险得分模型，即为非小细胞肺癌复发模型；所述甲基化风险得分模型具有多个得分阈值，这些得分阈值可将甲基化风险得分模型计算出来的得分对应到高得分组、中得分组、低得分组中的一组。本发明适用于非小细胞肺癌复发预测。

Description

基于DNA甲基化的早期非小细胞肺癌复发模型构建方法

技术领域

本发明涉及生物医学领域，特别涉及基于DNA甲基化的早期非小细胞肺癌复发模型构建方法。

背景技术

非小细胞肺癌是世界范围内的头号致死癌症。低剂量CT能够在人群中筛选出早期肺癌患者，对这些患者进行手术切除和术后辅助治疗，可以显著提高其整体生存率，但仍然有超过25％的早期(T2N0M0)患者会经历术后复发，而术后复发是预后较差的主要原因。考虑到肺癌存在广泛的异质性，在早期肺癌患者中精确识别其中的高复发人群是重要的科学问题。最新的美国国家综合癌症网络(National Comprehensive Cancer Network,NCCN)指南推荐6个高风险因素，包括肿瘤分化状态、楔形切除术和胸膜浸润等。然而，这些高风险因素的预测效能不足，亟须新的有效模型用于预测早期肺癌复发。

DNA甲基化是一种表观遗传学修饰，其信号改变与癌症的发生发展密切相关。Sandoval等人开发了一个基于HumanMethylation450芯片(后简称为450K)的非小细胞肺癌复发模型，用于预测I期非小细胞肺癌肺(non-small-cell lung cancer，NSCLC)患者的复发风险。该模型涉及5个基因，分别为HIST1H4F,PCDHGB6,NPBWR1,ALX1和HOXA9，这5个基因的甲基化程度用10个探针的beta值表示，分别为：cg10723962，cg22723502，cg12260798，cg16104915，cg12600174，cg18507379，cg18617005，cg26205771，cg07770968和cg14996220。对于一个基因对应多个探针的情况，用这些探针beta值的均值代表。在检测样本中，对于每一个基因，若其甲基化的beta值不小于0.4,认为其发生甲基化，否则认为该基因未发生甲基化。若该样本中这五个基因均未发生甲基化，或仅一个基因发生甲基化，认为该样本为低风险样本，否则为高风险样本。该非小细胞肺癌复发模型存在以下问题：

1.该模型在外部数据的效能不足。当把模型应用于目前常见癌症数据库(如TheCancer Genome Atlas)的早期肺腺癌样本中，基于该模型不具有区分生存的效能(P＝0.28)。

2.该模型的部分探针不能被简化亚硫酸氢盐测序(Reduced RepresentationBisulfite Sequencing,RRBS)所覆盖，即难以应用于RRBS检测平台。

3.该模型没有对阈值(0.4)附近的分类不清情况进行有效处理，如特定基因在第一个样本中的beta值是0.399，在第二个样本beta值是0.401，二者beta值差异不大(仅0.002)，但前者被认为是甲基化，而后者被认为是非甲基化，这显然是不合理的

发明内容

本发明要解决的技术问题是：提供一种基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，用以解决目前非小细胞肺癌复发模型对非小细胞肺癌复发预测不合理的问题。

为解决上述问题，本发明采用的技术方案是：基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，包括如下步骤：

步骤1：收集早期的非小细胞肺癌病例，并对病例患者的肿瘤以及癌旁正常组织提取DNA，并进行甲基化测序；

步骤2：对步骤1的测序数据，计算其甲基化位点的覆盖度和beta值；

步骤3：针对甲基化位点，搜索并合并其中的相关性高的甲基化位点，形成相应的甲基化区间；

步骤4：根据甲基化芯片的类型，挑选出能够被芯片探针检测的甲基化位点/区间作为训练数据建立甲基化风险得分模型，即为非小细胞肺癌复发模型；所述甲基化风险得分模型具有多个得分阈值，这些得分阈值可将甲基化风险得分模型计算出来的得分对应到高得分组、中得分组、低得分组中的一组；

甲基化风险得分模型的得分计算公式如下：

ModelScore＝1/(1+exp(-0.174+β_{chr10:101238903-101238917}×1.372+β_{chr18:22172331}×0.031+

β_{chr6:10421227}×0.503+β_{chr6:158516326}×0.532+β_{chr9:23820981}×1.131-β_{chr1:58576441}×1.237-

β_{chr11:9759297}×1.137-β_{chr17:72590998}×0.359-β_chr7:1587474×0.321))

其中，β_{chr10:101238903-101238917}、β_{chr18:22172331}、β_{chr6:10421227}、β_{chr6:158516326}、β_{chr9:23820981}、β_{chr1:58576441}、β_{chr11:9759297}、β_{chr17:72590998}、β_chr7:1587474分别是甲基化区间chr10:101238903-101238917以及甲基化位点chr18:22172331、chr6:10421227、chr6:158516326、chr9:23820981、chr1:58576441、chr11:9759297、chr17:72590998、chr7:1587474的beta值，ModelScore是模型的得分。

进一步的，步骤1可采用简化亚硫酸氢盐进行甲基化测序。如果为了降低成本，针采用步骤1采用靶向测序也是可以的。

进一步的，步骤2计算测序数据的甲基化位点的覆盖度和beta值之后，为了排除覆盖度不足的位点，可排除至少1/4样本中覆盖度小于10X的甲基化位点。

步骤3在合并高度相关的甲基化位点时，理论上很多其他方法也是可行的，比如滑窗法或者寻找断点的方法(circular binary segment等)。本发明可通过贪婪算法合并高度相关的甲基化位点，对于每一个甲基化位点，以该位点为种子进行贪婪搜索，分别将种子位点左侧和右侧最邻近的甲基化位点与种子位点合并，形成候选甲基化区间，分别计算两个候选甲基化区间的平均beta值和种子位点beta值在样本中的皮尔森相关系数，若二者中至少一个大于相关性阈值，则选取相关性更高的甲基化区间作为下一次搜索的种子区间，以此类推，反复迭代，直至左侧和右侧的相关系数均小于相关性阈值，搜索停止，最后形成最终的甲基化区间。

进一步的，通过随机矩阵的方式，申请人评估过0.3，0.5和0.7三个相关性阈值，发现0.3和0.5的假阳性较高，所以以上的相关性阈值优选为0.7。另外，在RRBS中相似性大于0.7的区域，在450k芯片的相关性也非常高。

具体的，所述的甲基化芯片类型包括：RRBS，WGBS，450K或850K芯片。

具体的，步骤4中所述显著差异指：Beta值差异不小于0.1，显著性P<0.05。

具体的，所述甲基化风险得分模型的多个得分阈值可利用高斯混合模型拟合甲基化风险得分模型的得分分布，并利用期望最大化求解模型参数，通过网格搜索优化得到。由实施例可知，所述甲基化风险得分模型的得分阈值可包括0.3和0.6，高得分组≥0.6，中得分组＝0.3～0.6，低得分组≤0.3。这里给出了最优模型阈值(≥0.6为高风险组，≤0.3为低风险组)，但不排除其他阈值的预测效能，即相同模型的不同阈值能够达到相同的目的。

本发明的有益效果如下：

1.本发明建立的模型针对早期(I/II期)非小细胞肺癌的患者，而现有大多数模型用所有肺癌样本进行建模，其模型不一定适用于早期患者。

2.本发明建立模型的训练数据用的是8版TNM分期中的IB期和IIA期，样本的均一性更好。

3.模型构建时同时考虑单个甲基化位点和多甲基化位点构成的甲基化区域，对于检测单位点甲基化常见的低敏感性和高噪音的限制，有更强的耐受性。

4.由表1可知，建立的模型既包含基因启动子区域，还包含增强子和转录因子结合位点；而其他基于甲基化的模型大多仅包含启动子区域。

5.建立的模型可应用于多种甲基化检测平台，比如RRBS、WGBS、450K和850K芯片平台，而大多数现有模型仅能应用于一种甲基化检测平台。

6.模型构建时，本发明利用了高斯混合模型(Gaussian mixture model，GMM)拟合模型的得分分布，在考虑样本覆盖度的情况下，允许存在得分中间组，即不确定组，能够解决得分在阈值附近分类不清的问题，从而降低错误分组的可能性。

7.建立的模型在多套外部数据中能得到良好的验证。

附图说明

图1是实施例的流程图。

图2是贪婪搜索的原理图。

图3是ENCODE 4个细胞系的H3K27ac信号图。

图4是模型得分的分布图。

图5是训练数据的KM曲线图。

图6是模型在外部数据GSE39279中的KM曲线图。

图7是模型在外部数据TCGA中的KM曲线图。

图8是实施例设计的列线图。

具体实施方式

为了解决目前非小细胞肺癌复发模型对非小细胞肺癌复发预测不合理的问题，本发明提供了一种基于DNA甲基化的早期非小细胞肺癌复发预测方法，主要包括建立非小细胞肺癌复发模型和使用非小细胞肺癌复发模型预测。

一、建立非小细胞肺癌复发模型

1、收集早期的非小细胞肺癌病例，并对病例患者的肿瘤以及癌旁正常组织提取DNA，并进行甲基化测序；

2、对步骤1的测序数据，计算其甲基化位点的覆盖度和beta值；

3、针对甲基化位点，搜索并合并其中的相关性高的甲基化位点，形成相应的甲基化区间；

4、根据甲基化芯片的类型，挑选出能够被芯片探针检测的甲基化位点/区间作为训练数据建立甲基化风险得分模型，即为非小细胞肺癌复发模型；所述甲基化风险得分模型具有多个得分阈值，这些得分阈值可将甲基化风险得分模型计算出来的得分对应到高得分组、中得分组、低得分组中的一组；

甲基化风险得分模型的得分计算公式如下：

ModelScore＝1/(1+exp(-0.174+β_{chr10:101238903-101238917}×1.372+β_{chr18:22172331}×0.031+

β_{chr6:10421227}×0.503+β_{chr6:158516326}×0.532+β_{chr9:23820981}×1.131-β_{chr1:58576441}×1.237-

β_{chr11:9759297}×1.137-β_{chr17:72590998}×0.359-β_chr7:1587474×0.321))

其中，β_{chr10:101238903-101238917}、β_{chr18:22172331}、β_{chr6:10421227}、β_{chr6:158516326}、β_{chr9:23820981}、β_{chr1:58576441}、β_{chr11:9759297}、β_{chr17:72590998}、β_chr7:1587474分别是甲基化区间chr10:101238903-101238917以及甲基化位点chr18:22172331、chr6:10421227、chr6:158516326、chr9:23820981、chr1:58576441、chr11:9759297、chr17:72590998、chr7:1587474的beta值，ModelScore是模型的得分。

本发明在确定所述甲基化风险得分模型的多个得分阈值时，可利用高斯混合模型拟合甲基化风险得分模型的得分分布，并利用期望最大化求解模型参数，通过网格搜索优化得到的。

二、使用非小细胞肺癌复发模型预测

非小细胞肺癌复发模型建立之后，就可以利用得到非小细胞肺癌复发模型进行非小细胞肺癌复发预测。

以下再结合附图和实施例对本发明做近一步说明。

实施例一种基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，如图1所示，流程如下：

1.收集57例IB/IIA期的非小细胞肺癌病例，对患者的肿瘤以及癌旁正常组织提取DNA，并进行简化亚硫酸氢盐测序。

2.数据清洗，对原始测序数据，利用TrimGalore软件去除接头，得到clean data，随后使用bismark软件计算甲基化位点的覆盖度和beta值，在后续分析中去除低覆盖度的甲基化位点，即排除在至少1/4样本中覆盖度小于10X的甲基化位点。

3.搜索并合并相关性高的甲基化位点，形成的甲基化区间。这里可通过贪婪算法合并高度相关的甲基化位点，如图2所示。对于每一个甲基化位点，以该位点为种子进行贪婪搜索，分别将种子位点左侧和右侧最邻近的甲基化位点与种子位点合并，形成候选甲基化区间，分别计算两个候选甲基化区间的平均beta值和种子位点beta值在样本中的皮尔森相关系数，若二者中至少一个大于0.7，则选取相关性更高的甲基化区间作为下一次搜索的种子区间，以此类推，反复迭代，直至左侧和右侧的相关系数均小于0.7，搜索停止，最后形成最终的甲基化区间。显然，只有部分甲基化位点能够通过贪婪算法扩展成甲基化区间，这使得甲基化区间和单碱基的甲基化位点同时进入模型训练。

4.以450K芯片平台为例，为了使模型可以应用于450K芯片平台，我们挑选出能够被450K芯片探针检测的甲基化位点/区间作为训练数据进行后续非小细胞肺癌复发模型建立。

构建甲基化时，对于训练数据，我们将复发患者的癌症样本作为高风险组，未复发患者的癌症样本和所有正常样本作为低风险组，利用双侧Wilcoxon秩和检验计算出两组中显著差异(Beta值差异不小于0.1，显著性P<0.05)的甲基化位点/区域。

6.对于上述差异显著的甲基化区域/位点，使用弹性网络回归进行特征选择和模型训练，以防止过拟合。其中，alpha＝0.2,罚分参数lambda通过十倍交叉验证获得，最终得到包含1个甲基化区间和8个甲基化位点的甲基化风险得分模型，即为非小细胞肺癌复发模型，该模型涉及的甲基化区间/位点以及其模型系数如表1所示，模型使用的参考基因组版本是hg38。我们结合UCSC基因组浏览器，对表1中甲基化区间/位点进行注释，发现模型中除了常见的基因启动子区域，还包括增强子和转录因子结合位点。其中，chr11:9759297在ENCODE的4个细胞系中均有强的H3K27ac的信号覆盖，如图3，证明该位点是一个广泛的增强子。

表1.甲基化区间/位点及其系数

甲基化风险得分模型的得分计算公式如下：

ModelScore＝1/(1+exp(-0.174+β_{chr10:101238903-101238917}×1.372+β_{chr18:22172331}×0.031+

β_{chr6:10421227}×0.503+β_{chr6:158516326}×0.532+β_{chr9:23820981}×1.131-β_{chr1:58576441}×1.237-

β_{chr11:9759297}×1.137-β_{chr17:72590998}×0.359-β_chr7:1587474×0.321))

其中，β_{chr10:101238903-101238917}、β_{chr18:22172331}、β_{chr6:10421227}、β_{chr6:158516326}、β_{chr9:23820981}、β_{chr1:58576441}、β_{chr11:9759297}、β_{chr17:72590998}、β_chr7:1587474分别是甲基化区间chr10:101238903-101238917以及甲基化位点chr18:22172331、chr6:10421227、chr6:158516326、chr9:23820981、chr1:58576441、chr11:9759297、chr17:72590998、chr7:1587474的beta值，ModelScore是模型的得分。

为了方便计算，实施例还针对以上的得分计算公式设计了如图8所示的列线图，列线图使用方法为：每个位点的甲基化值可以对应到第一行的Points，所有Points相加，对应倒数第二行的Total Points，倒数第二行的total Points与最后一行的Risk得分对应，即得到模型得分。

为了增强模型的实用性，需要确定模型得分阈值。在训练集中，模型得分在样本中呈现双峰分布，如图4所示，说明样本存在应该分成2-3个组，即两个峰和/或两峰之间的中间状态。为此，我们还利用高斯混合模型拟合模型得分分布，并利用期望最大化求解模型参数，得μ₁＝0.190,σ₁＝0.092,μ₂＝0.701和σ₂＝0.167。对于每个正态分布，我们在其均值周围一倍标准差中，通过网格搜索确定阈值。我们发现，中间组比例较高时，虽然模型的准确性较高，但是模型能够解释的样本量会较少。最终，在模型样本覆盖度大于80％的情况下，中间组在得分在0.3～0.6内，可以使模型具有最好的区分效能，即得分≥0.6为高得分组，0.3～0.6为中得分组，得分≤0.3为低得分组。若该样本模型得分≥0.6，则被认为是高得分组(高风险组)；若模型得分≤0.3,则被认为是低得分组(低风险组)；若在0.3～0.6之间，则被认为是低得分组(未确定组)。相对于低得分组，高得分组的患者可以建议密切随访或者提早干预治疗。

经测验，实施例的模型应用效果如下：

1、在训练数据中，高得分组比低得分组具有更高的复发风险，如图5所示。相对于低得分组，高得分组的风险比为9.14(P<0.001)，且该模型独立于常见的临床因素(包括年龄、性别，TNM分期、胸膜浸润、分化程度、吸烟状态和术后化疗)，如表2所示。

表2.cox回归分析结果

2、在外部数据GSE39279中，高得分组比低得分组具有更高的复发风险，如图6所示。相对于低得分组，高得分组的风险比为4.64(P<0.001)，且该模型独立于常见的临床因素(包括年龄、性别、TNM分期、组织学类型和吸烟状态)，如表2所示。

3、在外部数据TCGA中，高得分组比低得分组具有更高的复发风险，如图7所示。相对于低得分组，高得分组的风险比为2.32(P<0.05)，且该模型独立于常见的临床因素(包括年龄、性别、TNM分期、吸烟状态和术后化疗)，如表2所示。

与已有的临床常见复发指标和另一个非小细胞肺癌复发的模型相比，我本发明的非小细胞肺癌复发模型具有更好的预测效能，见表3。

表3.模型效能与比较

Claims (8)

1.基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，包括如下步骤：

步骤1：收集早期的非小细胞肺癌病例，并对病例患者的肿瘤以及癌旁正常组织提取DNA，并进行甲基化测序；

步骤2：对步骤1的测序数据，计算其甲基化位点的覆盖度和beta值；

步骤3：针对甲基化位点，搜索并合并其中的相关性高的甲基化位点，形成相应的甲基化区间；

步骤4：根据甲基化芯片的类型，挑选出能够被芯片探针检测的甲基化位点/区间作为训练数据建立甲基化风险得分模型，即为非小细胞肺癌复发模型；所述甲基化风险得分模型具有多个得分阈值，这些得分阈值可将甲基化风险得分模型计算出来的得分对应到高得分组、中得分组、低得分组中的一组；

甲基化风险得分模型的得分计算公式如下：

ModelScore＝1/(1+exp(-0.174+β_{chr10:101238903-101238917}×1.372+β_{chr18:22172331}×0.031+β_{chr6:10421227}×0.503+β_{chr6:158516326}×0.532+β_{chr9:23820981}×1.131-β_{chr1:58576441}×1.237-β_{chr11:9759297}×1.137-β_{chr17:72590998}×0.359-β_chr7:1587474×0.321))

其中，β_{chr10:101238903-101238917}、β_{chr18:22172331}、β_{chr6:10421227}、β_{chr6:158516326}、β_{chr9:23820981}、β_{chr1:58576441}、β_{chr11:9759297}、β_{chr17:72590998}、β_chr7:1587474分别是甲基化区间chr10:101238903-101238917以及甲基化位点chr18:22172331、chr6:10421227、chr6:158516326、chr9:23820981、chr1:58576441、chr11:9759297、chr17:72590998、chr7:1587474的beta值，ModelScore是模型的得分。

2.如权利要求1所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，步骤1采用简化亚硫酸氢盐进行甲基化测序。

3.如权利要求1所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，步骤2计算测序数据的甲基化位点的覆盖度和beta值之后，需排除至少1/4样本中覆盖度小于10X的甲基化位点。

4.如权利要求1所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，步骤3通过贪婪算法合并高度相关的甲基化位点，对于每一个甲基化位点，以该位点为种子进行贪婪搜索，分别将种子位点左侧和右侧最邻近的甲基化位点与种子位点合并，形成候选甲基化区间，分别计算两个候选甲基化区间的平均beta值和种子位点beta值在样本中的皮尔森相关系数，若二者中至少一个大于相关性阈值，则选取相关性更高的甲基化区间作为下一次搜索的种子区间，以此类推，反复迭代，直至左侧和右侧的相关系数均小于相关性阈值，搜索停止，最后形成最终的甲基化区间。

5.如权利要求4所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，所述相关性阈值为0.7。

6.如权利要求1所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，步骤4所述的甲基化芯片类型包括：RRBS，WGBS，450K或850K芯片。

7.如权利要求1所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，所述甲基化风险得分模型的多个得分阈值是利用高斯混合模型拟合甲基化风险得分模型的得分分布，并利用期望最大化求解模型参数，通过网格搜索优化得到的。

8.如权利要求7所述的基于DNA甲基化的早期非小细胞肺癌复发模型构建方法，其特征在于，所述甲基化风险得分模型的得分阈值包括0.3和0.6，高得分组≥0.6，中得分组＝0.3～0.6，低得分组≤0.3。

Images