基于5hmC-Seal测序技术检测糖尿病肾病患者血浆游离DNA中5-羟甲基胞嘧啶的变化

储金林，马睿瑶，李琳琳

(新疆医科大学药学院，乌鲁木齐 830017)

表观遗传学是一种非DNA序列改变而进行的可遗传修饰调控[1]。DNA 羟甲基化作为重要的表观遗传修饰之一，在个体发育，细胞成熟，分化和自我更新换代过程中起着关键的作用[2]。5-羟甲基胞嘧啶(5-hydroxymethylcytosine, 5hmC)是5-甲基胞嘧啶的羟甲基化形式，可通过TET(Ten-eleven translocation)酶氧化形成，是活性DNA去甲基化过程中相对稳定的中间体，被认为是重要的表观遗传特征[2]。近年来，血浆游离DNA(Cell-free DNA, cfDNA)的5hmC谱已被证明是癌症及糖尿病等相关疾病早期检测和进展监测的有力工具[3-5]。糖尿病肾病(Diabetic kidney disease, DKD)是糖尿病患者主要的微血管并发症，随着人口老龄化的加重，糖尿病的发病率呈逐年升高趋势，研究发现约30%～40%的糖尿病患者发展为DKD[6]。由于患者一直处于高血糖状态，一旦发病，进展较快，病情很难得到有效控制，并且目前临床上对于DKD的早期诊断较差，会出现部分患者漏诊错诊，确诊时已难以治愈疾病。另外，由于DKD患者体内会出现复杂的一系列代谢紊乱，在进展到终末期肾病时，相对于其他类型的肾脏疾病，对于DKD的治疗更加困难并且预后不好[7]。因此，及时防治DKD对于缓解患病压力有着重大的意义。此外，目前已有部分对DKD的基因表达谱的研究[8-10]，但羟甲基化基因的具体通路调控是否参与DKD形成过程中的分子机制还需进一步研究。因此，本文对糖尿病和DKD患者的外周血样本进行5hmC-Seal测序，拟寻找与DKD相关的通路中的羟甲基化生物标志物，现报道如下。

1 资料与方法

1.1 一般资料收集北京大学第三医院的糖尿病患者血浆(14例，年龄47.9±6.9岁，男性9例，女性5例)，以及糖尿病肾病患者血浆(17例，年龄51.1±16.7岁，男性9例，女性8例)。纳入标准如下：糖尿病的诊断基于美国糖尿病协会于2020年公布的糖尿病医疗护理标准[11]，在至少8 h不摄入热量的情况下，空腹血糖值≥126 mg/dL(7.0 mmol/L)和/或口服75 g溶解在水中的无水葡萄糖的葡萄糖负荷实验，2 h后血糖值≥200 mg/dL(11.1 mmol/L)。DKD的诊断标准依据KDIGO 2020临床实践指南和2014年中国糖尿病学会糖尿病肾病防治共识，主要包括尿白蛋白/肌酐比值(UACR)≥30 mg/g和/或估算的肾小球滤过率(eGFR)<60 mL·min-1·(1.73 m2)-1，并持续三个多月[12]。排除标准：有证据证明还有其他疾病；同时接受其他临床研究的试验性治疗。本研究严格按照赫尔辛基宣言进行，并得到了北京大学第三医院伦理委员会的批准(批准号：LA2020390)。所有受试人员均知情同意。

1.2 方法

1.2.1 cfDNA的提取及测序 取患者外周血2 mL，在4 ℃离心机以1 350×g离心15 min，取上清；4 ℃，13 500×g离心5 min，上清即血浆。血浆中cfDNA的提取主要采用QIAamp循环核酸试剂盒(QIAGEN)；并使用KAPA Hyper Prep Kit试剂盒以及DNA Clean&Concentrator Kit试剂盒进行建库；最后将患者的文库混合使用Nextseq500/550 High Output Kit v2(75 cycles)测序试剂盒在Illumina NextSeq500测序平台上进行双端38-bp进行高通量测序，可通过文库构建时加入的Index来区分每个患者。用DEseq2包找出差异羟甲基化基因(Differentially hydroxymethylated genes, DhMGs)，筛选条件为|log2FC|>0.5且FDR<0.05(FC为差异倍数，FDR为校正后的P值)。

1.2.2 基因集变异分析(Gene set variation analysis，GSVA) 将筛选得到的DhMGs进行GSVA富集中的KEGG通路分析，主要分析方法通过GSVA R包来实现，比对文件为：c2.cp.kegg.v6.2基因集。GSVA可通过基因表达数据对基因富集结果进行量化，我们使用FDR<0.05和|log2FC|>0.2为条件找到一些有生物学意义的差异信号通路进行后续分析。

1.2.3 DhMGs调控网络的构建与分析 利用STRING(https://string-db.org/)数据库对感兴趣的差异信号通路里的DhMGs进行相互作用分析，并建立蛋白互作网络(Protein-protein interaction, PPI)。通过Cytoscape3.9.1(https://cytoscape.org/)软件构建基因调控网络，再使用MCODE插件筛选出关键模块里的关键枢纽基因。

1.2.4 受试者工作特征曲线(Receiver operating characteristic,ROC) 对最终的关键枢纽基因进行ROC分析，判断用于诊断糖尿病肾病的能力，若曲线下面积(Area under the curve,AUC)≥0.8，则认为该基因具有较高的临床诊断价值[13]。

2 结果

2.1 DKD患者5hmC整体分布特征全基因组5hmC图谱显示(图1A)，与糖尿病患者(14例)比较，DKD患者(17例)的整体5hmC分布水平显著降低(P<0.001)，5hmC分布水平在启动子，内含子，外显子，5’UTR，3’UTR等基因功能区域内都有显著差异(P<0.001)。根据主成分分析可得，所有5hmC位点可以将糖尿病患者和糖尿病肾病患者区分开，见图1B。

注：A，5hmC分布分布；B，主成分分析图，D为糖尿病组(黄色)，DKD为糖尿病肾病组(红色)。

2.2 DKD和糖尿病患者DhMGs的差别以FDR<0.05和|log2FC|>0.5的条件下筛选DhMGs，在DKD组中上调基因(红色)有572个，下调基因(绿色)有1 072个，见图2A。这些差异位点主要分布在内含子、基因间和启动子区域，见图2B，与图1A显示的结果相一致。

注：A，火山图；B，差异基因分布区域图。

2.3 DhMGs的GSVA信号通路富集分析GSVA无监督聚类通路富集(c2.cp.kegg.v6.2基因集)分析结果显示，以|log2FC|>0.2，FDR<0.05为筛选条件，发现神经活性配体受体相互作用通路在DKD组中显著上调，肾脏上皮通路，趋化因子信号通路，T细胞受体信号通路，黏附连接通路和Wnt信号通路在DKD组中显著下调，差异有统计学意义，见表1，同时也发现这些通路能将糖尿病患者和糖尿病肾病患者区分，见图3。有多项研究认为Wnt信号通路与糖尿病肾病的发生机制密切相关[14-16]，因此选定此通路为下一步分析的目标通路，分析此通路中的羟甲基化基因能否作为区分糖尿病和糖尿病肾病的生物标志物。

表1 GSVA富集结果

注：D为糖尿病组(蓝色)，DKD为糖尿病肾病组(红色)。

2.4 关键基因筛选结果利用STRING数据库构建了富集在Wnt信号通路中的18个DhMGs的PPI网络。以Combined scores>0.4为筛选条件，结果共显示18个节点，50条边，并使用Cytoscape软件进行可视化，结果见图4。接着使用MCODE插件依据MCODE_score>6.8识别出关键模块(图4左)，包括11个关键基因，分别为CTNNB1、GSK3B、CREBBP、CUL1、NFATC1、CCND3、WNT11、PRKACB、DKK1、SMAD3、CCND1。颜色越深，degree值越大，说明这个基因在疾病过程中越关键。

图4 PPI网络图

2.5 DKD和糖尿病患者的5hmC标志物将通过MCODE算法筛选得到的11个关键基因进行ROC曲线分析，结果表明11个基因组合能很好地诊断糖尿病肾病，AUC=1，说明诊断准确性能达到100%(图5A)。接着分析每个基因的诊断准确性，发现CCND1、WNT11、CREBBP、CUL1、CTNNB1、GSK3B、PRKACB、NFATC1、DKK1的诊断准确性较高，满足AUC≥0.8，具有较好的诊断价值(图5B-J)。另外，SMAD3和CCND3的羟甲基化水平用于诊断糖尿病肾病的准确性不高，AUC<0.8(图5K-L)。

图5 ROC分析结果

3 讨论

早期发现DKD对于防止进展为肾功能衰竭很重要，因此，找到一些有意义的生物标志物可能对DKD的诊断和治疗产生重大的临床影响。本项研究的重点是从患者血浆cfDNA 中发现DKD特异性羟甲基化生物标志物，从而对DKD进行分子诊断。

本研究的数据表明，DKD患者5hmC分布水平显著低于糖尿病患者，并且5hmC信号主要富集在启动子、外显子和内含子等基因功能区。有研究显示，基因功能区发生羟甲基化可以影响基因的转录水平，表明这些区域通过5hmC在调节基因表达中发挥关键作用[17]。本研究对糖尿病和DKD患者血浆cfDNA中5hmC修饰基因进行差异分析，找到了1 644个DhMGs，这些5hmC差异性位点在基因体上主要分布于启动子、外显子和内含子等区域，基因功能区发生羟甲基化可以影响基因的转录。所以，进一步通过GSVA富集分析探究这些羟甲基化修饰基因的生物学功能，发现糖尿病和DKD两组之间存在显著差异的疾病通路，如：神经活性配体受体相互作用通路，肾脏上皮通路，趋化因子信号通路，T细胞受体信号通路，黏附连接通路和Wnt信号通路。其中，Wnt信号通路已被多项研究证明与DKD发生的分子机制密切相关，同时后续数据分析也得出CTNNB1为关键枢纽基因。CTNNB1基因编码的蛋白质β-catenin在细胞黏附中起关键作用并 调节 Wnt 介导的转录。先前的研究[18]揭示了DKD 患者和小鼠模型的肾小球和足细胞中 Wnt/β-catenin 通路转录物和蛋白质的表达增加。DKD中Wnt/β-catenin 通路的上调不仅可以促进足细胞存活，还可以导致细胞脱离和足细胞丢失。同时，这篇报道也称足细胞中Wnt/β-catenin通路的下调可能对终末分化很重要。并且，它又会增强足细胞凋亡的易感性。因此，CTNNB1的平衡表达对于DKD中肾小球滤过屏障的维持至关重要。也有文献[19]报道六味地黄丸可通过下调Wnt/β-catenin 信号通路，对糖尿病肾病相关性骨质疏松症有一定的治疗效果。

Wnt信号通路中的糖原合成酶激酶3β(GSK3β)已报道[20]在DKD发展初期其表达水平已经升高，并与糖尿病肾损伤的症状进展相关，并认为尿脱落细胞中的GSK3β可能作为一种新的生物标志物用于预测DKD进展。而本文也得出血浆cfDNA中的羟甲基化GSK3β可以作为一种诊断DKD的生物标志物。本研究发现Wnt信号通路中的WNT11基因目前还没有报道与DKD相关的文献，事实上这个基因在肾脏的发育过程中发挥重要作用[21]，WNT11功能损害导致肾小管异常和继发性肾小球囊肿。以上可以揭示Wnt信号通路的相关基因参与了糖尿病患者的肾损伤进程，尽管有些基因还未报道与DKD的关系，但因为有对临床患者的血浆或者病灶直接来源的相关测序手段，能够发现新的、足够准确的生物标志物，或者发现可能的、关键的发病机制，当然，这需要后续大量实验的深入验证。

本研究ROC分析结果显示，Wnt信号通路富集到的9个关键基因诊断准确性较高(AUC≥0.8)，分别有CCND1、WNT11、CREBBP、CUL1、CTNNB1、GSK3B、PRKACB、NFATC1、DKK1。但因样本数量有限，所以认为这些基因具有一定的潜在诊断价值。

综上所述，本研究通过分析血浆样本基因组的高通量数据，初步探讨了参与糖尿病出现肾损伤的羟甲基化修饰基因以及这些基因参与的信号通路，发现了CCND1、WNT11、CREBBP、CUL1、CTNNB1、GSK3B、PRKACB、NFATC1、DKK1具有潜在的诊断DKD的价值，为DKD的后续机制研究及诊疗策略奠定了基础。