DTC预后标志物及其应用、DTC预后评估模型的构建方法与流程

dtc预后标志物及其应用、dtc预后评估模型的构建方法
技术领域
1.本发明属于癌症预后评估技术领域，特别涉及dtc预后标志物及其应用、dtc预后评估模型的构建方法。


背景技术：

2.甲状腺癌是内分泌系统中最常见的恶性肿瘤之一，在我国发病率呈逐年上升趋势，每年新发病例占全部恶性肿瘤的5％～7％。分化型甲状腺癌是甲状腺癌中恶性程度较低的一种类型，但由于目前缺少敏感性和特异性较高的血清学指标，导致该病早期易出现漏诊和误诊。手术切除联合131放射粒子治疗是目前治疗甲状腺癌的主要方案；但相关研究显示，术后3年内20％～30％的分化型甲状腺癌患者可出现淋巴结转移或远处转移，其1年生存率小于30％。远低于早期被发现复发高风险，且接收规范治疗患者的生存率。
3.分化型甲状腺癌(dtc)有许多预后分类系统，以美国癌症联合委员会/国际癌症控制联盟/肿瘤、淋巴结和转移(ajcc/uicc tmn)分期使用最广泛。但随着对dtc疾病特点的认识不断加深，我们意识到对于dtc这类相对缓慢进展、患者存活期长的恶性肿瘤，仅根据初始治疗时的临床病理特征划分风险对预测患者最终临床结局存在局限性，可能的原因包括：
⑴
未考虑dtc的肿瘤组织学亚型导致分层不足[如乳头状甲状腺癌(ptc)的高细胞亚型、柱状细胞弥漫硬化型、实体亚型和滤泡状甲状腺癌(ftc)的广泛浸润型，是生物学侵袭性较高的dtc亚型]；
⑵
部分分化较差、病情进展较快的dtc，在初始阶段可能难以识别，不能做出准确的生物学行为预判；
⑶
未考虑治疗的影响，而手术范围、是否行淋巴结清扫、是否行131i治疗及131i亲和力、tsh抑制治疗等都可能使疾病的预后发生变化。这种局限性会导致根据初始风险评估而制定的某种治疗或随访策略，面对不同的病情发展走势，变得不再适用，引发诊疗不足或过度。现有风险分层的实用性受到了挑战。因此需要进一步研究以更好地定义预测复杂临床情况患者的预后风险。
[0004]
warburg效应，也称为有氧糖酵解，是一种现象，其中各种类型的癌细胞的特征是葡萄糖过度转化为乳酸作为其能量底物，而不管氧气水平如何。越来越多的证据表明，在癌症中加速糖酵解影响治疗结果，当与邻近的正常组织(相比大多数癌症表现出的葡萄糖摄取)增加显著。此外，据报道，糖酵解增加可促进血管生成和侵袭性癌症的生长。由糖酵解肿瘤细胞产生的乳酸在抑制抗癌免疫细胞中起关键作用，然后在抗癌治疗后促进肿瘤复发。乳酸在肿瘤微环境(tme)中的大量积累，将细胞外ph值降低至6.0-6.5，阻碍了t细胞的功能和增殖速度。tme中高浓度的乳酸会影响抗肿瘤治疗，这导致抑制糖酵解途径的建议，因此乳酸的产生可能提供一种有效和潜在的策略来增强抗癌药物。
[0005]
癌症基因组图谱(the cancer genome atlas，癌症基因图集)通过大规模高通量基因组测序和芯片技术绘制人类肿瘤的基因组图谱，旨在探索肿瘤发生发展的潜在分子机制，提高人类诊断、治疗和预防肿瘤的能力。目前，数据库包括详细的临床数据、基因组突变数据、转录组数据、mirnas测序数据等。geo(gene expression omnibus)数据库主要由芯片数据组成，如snp芯片数据、比较基因组数据、microrna芯片数据等，此外还包含一些非芯片
数据，如基因表达序列分析数据、核糖体序列标签连续分析数据、质谱数据、蛋白质组学数据和新一代高通量测序数据等。随着人类基因组计划的完成和生物信息学的迅速发展，对这些数据库资源中相关的肿瘤数据进行搜索、处理和分析，就可以找出治疗癌症的关键靶点。
[0006]
申请号为202011004518.4的中国专利公开了“骨肉瘤预后标志物及预后评估模型”，是通过识别骨肉瘤预后相关的糖酵解基因，构建骨肉瘤预后模型。该技术在mrna层面对骨肉瘤的特征进行整合，所以采用组学技术中最常见的mrna差异表达对甲状腺癌患者进行预后预测是比较可靠的。但该模型设计没有纳入临床信息，需要找到更合理的、便于推动临床分期改变的mrna构建模型。


技术实现要素：

[0007]
本发明通过生物信息学分析，对dtc的rna-deq数据和疾病关联数据进行大规模系统分析探索糖酵解相关基因和dtc复发的规律，阐明与dtc复发相关的糖酵解基因的生物学特性。建立dtc复发预后模型，筛选关键糖酵解相关基因，得出与dtc预后有关的生物靶标，为未来研究阐明dtc的致病机理提供了理论依据，也为今后dtc的诊断治疗及病人预后提供新的思路。
[0008]
本发明通过limma算法筛选出在dtc和正常甲状腺组织中差异表达的糖酵解相关基因。limma是一种r/bioconductor软件包，为分析基因表达实验数据提供了集成解决方案。limma涵盖了在基因表达分析的每个主要步骤的功能，从数据导入，预处理，质量评估和标准化，到线性建模，差异表达分析和基因特征分析等。limma包的核心组件是能够将基因方式的线性模型拟合到基因表达数据，以评估差异表达分析。
[0009]
本发明中，ssgsea算法被用于计算dtc患者的糖酵解活性评分。ssgsea算法计算基因集评分：在分析来自微阵列数据的全基因组转录谱时，一个典型的目标是找到与由特定表型(例如，肿瘤与正常细胞)定义的不同样本类别显着差异相关的基因。这些发现可用于深入了解潜在的生物学机制或对新样本进行分类(预测表型)。基因集富集分析(gsea)通过评估与特定生物过程、染色体位置或实验结果相关的先验定义的基因集是否在差异列表的顶部或底部富集，从而解决了这个问题。表达的基因按样本类别中基因表达差异的某种度量进行排序。排名指标的例子是分类表型的倍数变化(例如，肿瘤与正常细胞)和pearson相关性的连续表型(例如，年龄)。富集为基因组成员的协调上调或下调以及某些相应生物过程的激活或抑制提供了证据。gsea根据数据集中样本集合的表型差异生成基因集的富集分数，而ssgsea为每个样本和基因集配对计算单独的富集分数，与表型标记无关。通过这种方式，ssgsea将单个样本的基因表达谱转换为基因集富集谱。基因集的富集分数代表生物过程的活动水平，其中基因集的成员被协调上调或下调。这种转变使研究人员能够根据生物过程和途径的活动水平而不是通过单个基因的表达水平来表征细胞状态。在处理转换后的数据时，目标是找到在感兴趣的表型中具有不同活性的生物过程，并使用这些过程活动的测量来表征表型。因此，这里的好处是ssgsea投影将数据转换到更高级别(路径而不是基因)空间，表示可以应用分析方法的一组更具生物学可解释性的特征。
[0010]
在本发明中，lasso算法被用于筛选出最稳健的基因，并基于cox回归分析，构建了dtc复发预测模型。lasso-cox算法构建基因预测复发模型：cox比例风险模型是医学中最流
行的将协变量与生存数据联系起来的方法之一。已有研究表明cox模型至少与神经网络和其他机器学习模型一样好，甚至更好。带有套索惩罚变量选择的cox回归模型通常用于从数千个基因分析中识别一些预后生物标志物，并为更简单和更便宜的临床应用获得简约模型。在考虑协变量的数量p(通常可以是20,000个基因产物)与数据库中的患者数量相关时，n(通常只有几百个)，由于高维数会出现各种问题，包括所选基因缺乏稳定性和过度拟合。这个p＞＞n问题被称为“维度灾难”。已经为广义线性模型(例如cox回归)提出了套索泛化，以提高性能和稳定性。特别是，弹性网和自适应弹性网是正则化程序，可以在存在高度相关变量的情况下克服套索的一些稳定性问题。这尤其发生在基因特征中。自适应弹性网还确保了额外的理论特性，以恢复真正的生物标志物。作为类似套索的方法，岭惩罚允许控制估计量的方差。虽然没有选择，但岭回归已被证明有望通过高维微阵列数据进行可靠的生存预测。
[0011]
rpart算法被用于本发明的决策树构建。基于rpart算法构建临床决策树：决策树是以历史数据为基础的归纳学习算法。它从一组无次序、无规则的数据中推理出决策树表现形式的分类规则，并使用该规则预测未来情况，被广泛地应用于数据挖掘领域。决策树由节点、分支和叶节点三个基本部分组成，树中的节点是一个测试条件，这个测试条件将决策树分为多个分支，每个分支代表了该测试条件每个可能的答案，位于决策树最顶端的节点为根节点。每个分支会连接另一个节点，或者连接决策树的末端(即叶节点)，决策树建好后，需要预测的数据对象从根节点出发，根据经过节点的测试条件选择相应分支，并最终到达叶节点。从根到叶节点的一条路径对应着一条合取规则，整个决策树对应一组析取表达式规则。rpart需要从cran上下载安装。使用时以rpart()和prune()两个函数为主，前者是用来拟合一个树模型，后者用来根据“成本复杂性”对生成的树进行剪枝。之所以要剪枝，是因为若不加任何限制，最后生成的树必然能完全拟合原始数据，这样的树在实际应用中毫无意义，因为树的枝节太多，而不能反映数据内在大规律；而从另一个极端情况来看，若树的枝节太少，那么必然也会带来很大的预测误差。综合看来，要兼顾树的规模和误差的大小，因此通常采用一个叫“成本复杂性”的标准来对树进行限制，最后达到的目的是使误差和数的规模都尽可能小。误差的计算通常基于交叉验证等方法，即用一部分训练样本建立模型，而剩下的样本用来作验证看模型的预测误差大小。
[0012]
本发明的优势在于：
[0013]
1)本发明基于糖酵解相关基因构建dtc预后预测模型，实现了分化型甲状腺癌患者预后的风险分层，显著地将高低风险的患者分开，进而可以预测分化型甲状腺癌的临床结果，指导个体化治疗，具有较高的临床应用价值。
[0014]
2)本发明所构建模型经过了在训练数据集和验证数据集中均表现为极高的预测灵敏性和准确性，因此具有较高的临床推广价值。
[0015]
3)本发明找到与分化型甲状腺癌复发相关的免疫细胞亚型，并探索了复发相关免疫细胞标志物。这为之后开发药物提供了理论依据和潜在靶点。
[0016]
4)本发明建立的模型是从开源的公共数据库中下载甲状腺癌患者的基因表达数据和临床数据，解决了样品收集难，测序费用高，以及对病人随访的问题。
附图说明
[0017]
图1.糖酵解评分与甲状腺癌的复发和各种临床病理特征有关。(a)火山图显示tcga队列中正常甲状腺和癌症样本之间的degs。(b)维恩图展示了175个候选糖酵解相关degs的筛选方法。(c)术后复发的患者中，糖酵解评分明显升高。(d)kaplan-meier分析显示，糖酵解分数较高的患者表现出较差的rfs。(e)热图描述了糖酵解分数和临床病理特征之间的相关性。
[0018]
图2.建立一个糖酵解相关的基因签名。(a)51个grg在5个geo队列(gse33630、gse60542、gse58545、gse35570和gse27155)中共同差异化表达。(b)27个差异化的糖酵解基因具有预后价值。(c)使用cox算法构建4个基因的糖酵解风险评分(grs)。(d)在合并或各自的grg中，按糖酵解相关基因对dtc复发的roc曲线，有无复发的患者。(e-g)训练、验证和测试集中4个糖酵解相关基因的kaplan-meier曲线。
[0019]
图3.基于grs和传统特征，建立了dtc患者复发的综合风险因素数据库。(a,b)临床病理特征、grs和rfs之间关联的单变量cox回归分析。(c)通过roc曲线分析比较grs与传统临床指标的诊断效果。(d)多变量cox回归分析表明，n期和grs是dtc患者rfs的独立风险因素。
[0020]
图4.基于grs的提名图和临床决策树的构建。(a)基于grs和n期的提名图的构建。1年、3年和5年复发的风险用箭头表示。例如，使用n0和高grs的患者。(b)残差分析以测试提名图是否成立。(c)校准分析以测试提名图的预测准确性。(d)roc分析，比较提名图和单独的n阶段对grs的预测敏感性。(e)dca分析，比较提名图、grs和n-stage的净生存收益。(f)基于grs和n期的临床决策树的构建。(g)k-m分析，比较不同风险组的生存差异。
[0021]
图5.高危患者表现出免疫途径的激活。(a)火山图显示高危和低危患者之间的degs。(b)高危组的gsea结果。(c)上调的degs的go术语。(d)不同风险组的免疫状态评分。(e)糖酵解评分和免疫评分的相关性分析。(f)ssgsea算法来评估不同风险组的免疫细胞浸润程度。(g)timer算法评估grs和免疫细胞浸润程度之间的相关性。
具体实施方式
[0022]
为使本发明要解决的技术问题、技术方案和优点更加清楚，下面将结合具体实施例进行详细描述。
[0023]
dtc复发风险评估系统的构建具体包括以下步骤：
[0024]
1)检索数据并处理rna序列数据：从tcga公共数据库(the cancer genome atlas，https://portal.gdc.cancer.gov/)中确认并下载183例thca队列患者的组织样本rna-seq数据；从geo公共数据库(gene expression omnibus，https://www.ncbi.nlm.nih.gov/geo/)中确认并下载包含甲状腺癌相关组织与正常组织信息的gse28735的表达谱，其中包含了45个肿瘤样本和45个正常甲状腺组织；
[0025]
2)定义致癌糖酵解差异基因集：在genecard网站(http://www.genecards.org/)以“glycolysis”为关键词进行检索，下载糖酵解相关基因共2538个。在r语言环境下，用edger和limma包对tcga-thca队列的表达矩阵进行标准化、数据拟合及差异分析，并以|log2fc|》1，fdr《0.05为标准筛选差异显著的基因(图1a)(结果以火山图展示)。将糖酵解相关基因与log2fc》1，fdr《0.05的上调基因取交集，得到致癌的糖酵解差异基因集，包括175
个基因(图1b)；所述175个致癌的糖酵解差异基因为：ereg、actbl2、tnfsf11、gckr、styk1、cdkn2b-as1、aldh3b2、retn、grp、mmp13、hnf1b、kiss1、h1-3、syt5、slc5a12、styxl2、myo1a、linc01554、shh、atp6v0a4、kng1、igfbp1、mmp3、chat、gabra1、agrp、dusp9、serpinb5、uca1、serpinb2、trpc5、adcy8、klk6、hp、vgf、mmp1、dhrs2、linc00242、mtfr2、diaph3、nek2、pimreg、tpsab1、ren、gsdma、alk、nudt11、gpt、gpr162、cpt1c、birc5、aurkb、gls2、nos2、nox1、ppp1r1b、cxcl8、pmaip1、cadn2a、hk3、adam12、tp63、s100b、nmu、trim63、tubb3、ccl20、klk10、cxcl14、chi3l1、slpi、muc1、lamb3、apoc1、igf2bp2、abcc3、cdh3、pnp、erbb3、nrip1、nr1d1、p3h3、ethe1、ptp4a3、pdlim7、pc、dysf、ecm1、ca12、cdh2、adm、isg15、pxdn、ak1、ndufa4l2、lpl、dram1、ddb2、mthfd1l、me3、sphk1、slc1a5、clip3、runx1、aldh3b1、adora1、ryr1、thrsp、plxna4、ca2、hdac9、mecom、nos3、gpr4、mki67、cenpf、kl、nox4、rin1、mlxipl、trem2、tyms、e2f1、tk1、cytor、irf5、bbc3、diras3、macc1、runx2、spp1、myo1g、bcat1、scd、fcer1g、ccl18、sfn、aldh1a3、fn1、anxa1、serpina1、igfbp3、mcam、insr、aqp1、lgals1、inf2、srebf1、bcl2l1、tgm2、itgb4、nt5e、tgfbi、p4ha2、cd276、bax、tnfrsf10b、epha2、tgfb1、sting1、cdkn1a、myh10、efhd2、gpx1、cd44、ccnd2、ccnd1、met、bhlhe40、lrrk2、krt19、cldn1、timp1、ctsh、lgals3。
[0026]
3)探索糖酵解对dtc患者预后的影响：基于上一步定义的致癌糖酵解差异基因集，用gsva包计算每一个患者的糖酵解评分。wilcoxon检验被用于比较rfs( )和rfs(-)患者之间的糖酵解评分高低(图1c)；用survival包选取糖酵解评分的最佳截断值，将患者分为高糖酵解评分组和低糖酵解评分组；并用kaplan-meier方法统计糖酵解评分对dtc患者预后的影响(图1d)。进行糖酵解分数和临床病理特征之间的相关性分析，结果如图1e所示。
[0027]
4)lasso算法选取稳健的糖酵解基因：首先在5个跨平台geo数据集中验证175个致癌糖酵解基因的表达水平。得到51个差异共表达基因后(图2a)，用kaplan-meier方法进一步筛选出27个有复发意义的基因(图2b)。最后用lasso算法进行十倍交叉验证，挑选出最稳健的4个基因，即adm，cd44，mki67，tyms(图2c-图2d)。
[0028]
5)构建糖酵解相关风险评分(grs)：用多因素cox比例风险回归分析degs与患者生存时间的关联性。根据多因素cox分析结果，推断出糖酵解相关风险评分。即grs＝0.176
×
exp(adm) 0.73
×
exp(cd44) 0.452
×
exp(mki67) 0.833
×
exp(tyms)，计算风险分数。再根据最佳风险评分将病人分为高风险组和低风险组，并进行生存分析绘制患者的生存曲线和生存状态图(图2e-f-g)，由此建立起联合grg的预后模型。绘制grs的roc曲线判断其预测能力，auc值》0.7则说明该模型具有一定的预测能力。
[0029]
6)构建grs-临床综合模型和临床决策树：用多因素cox例风险回归分析grs和临床病理特征对患者复发的影响(图3)。并基于多因素cox回归的结果，用regplot包构建grs-临床综合模型。然后将诺莫图展示综合模型以量化风险。再通过残差检测法检验诺莫图是否成立。绘制诺莫图的roc曲线、校正曲线和dca曲线来验证综合模型预测效能的灵敏性、准确性和净临床收益。最后，开发出grs-临床综合模型对应的临床决策树以帮助临床决策(图4)。
[0030]
7)高风险人群机制探索：临床决策树将患者分为高风险组和低风险组。用limma算法统计并筛选出两组之间的差异表达基因(|log2fc|》＝1，fdr《0.05)。并将差异基因上传至metascape网站，进行go和kegg富集分析。并使用gsea算法探索不同风险人群中的通路富
集情况(图5)。
[0031]
以上所述是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明所述原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。