子痫前期(preeclampsia)是妊娠期特有的一种多系统进展性疾病，其特点是妊娠20周以后出现新发高血压和蛋白尿，或出现新发高血压和终末器官功能障碍伴或不伴蛋白尿，占所有妊娠的2%～8%[1]。在临床实践中，子痫前期通常分为早发型(妊娠34周内)和晚发型(妊娠34周后)、轻度和重度(基于血压、临床表现和蛋白尿程度)[2]。研究证实从妊娠早期开始使用阿司匹林可以降低子痫前期的患病率[3-4]。然而，目前对阿司匹林用药的适宜人群、开始和结束用药的时机及剂量等仍在不断探索中，阿司匹林在临床研究中的应用并未如预期可明显降低子痫前期的发病率，其用于预防子痫前期仍存在一定的局限性[5]，当前唯一有效治疗子痫前期的方法是终止妊娠。因此对子痫前期的早期预测和诊断极其重要，了解子痫前期发生和发展的分子机制可能会改善治疗现状。

由于患有子痫前期的孕妇症状通常在分娩后缓解，且分娩后可以检测到胎盘的组织病理学变化，因此胎盘功能不全一直被认为是导致子痫前期的根本原因，胎盘标志物对预测子痫前期可能具有特异性和灵敏性。由于样本量小或数据分析不充分，许多芯片研究未能确定独特的胎盘分子标志物。转录组测序(RNA-seq)可对包括编码和非编码转录本在内的转录组进行全测序，对探索疾病机制和生物标志物有参考价值[6]。Kaartokallio等[7]针对子痫前期胎盘的RNA-seq数据进行了分析，观察到子痫前期患者的胎盘存在血管功能和免疫平衡紊乱，并鉴定出了一些或许可以预测和诊断子痫前期的差异表达基因。Ren等[8]通过RNA-seq发现早发型和晚发型重度子痫前期的分子机制不同，晚发型轻度子痫前期可能没有胎盘特异性致病因素。这些研究表明分析子痫前期胎盘组织的转录组数据将有助于发现子痫前期发展的分子机制，以及筛选出预测、诊断子痫前期的标志物。

随着各种生物信息学工具和公共数据库的出现，研究人员利用生物信息学方法能够高效且经济地从高通量测序数据中挖掘出疾病的致病基因及有潜力的诊断靶标。一般的生物信息学方法在处理高维数据时往往存在困难，因为高维数据之间的相关性和交互作用复杂且难以解释。加权基因共表达网络分析(weighted gene co-expression network analysis，WGCNA)是一种描述大样本基因间强关联性的高级生物信息学方法[9]。WGCNA的独特优势在于能将基因表达数据转化为共表达模块，从而深入了解可能导致相关表型特征的信号网络[10]。WGCNA可用于寻找与疾病高度相关的基因模块，并能够将模块的特征基因或关键基因与样本的表型特征联系起来。该方法被广泛应用于各种疾病研究，对鉴定候选生物标志物或治疗靶点有很大帮助[11-13]。另外机器学习算法能够从新视角分析大群体基因测序或微阵列数据，随机森林算法的优势在于其分类表现较好，抗过拟合能力较强，具有较好的鲁棒性，对噪声和异常值有较好的容忍性[14]。来自美国国家生物技术信息中心(National Center for Biotechnology Information，NCBI)的基因表达综合数据库(Gene Expression Omnibus，GEO)是当今最大、最全面的公共测序数据资源，包含了多种疾病的多组学高通量测序的原始数据，提供了大量的数据资源。

本研究整合了来自NCBI GEO的4个RNA-seq数据集(每个数据集都包含子痫前期孕妇的胎盘组织样本和正常孕妇的胎盘组织样本)，对比鉴定了4个数据集共有的差异表达基因并对这些差异表达基因进行了功能富集，然后通过WGCNA构建了子痫前期的共表达网络，鉴定出与子痫前期高度相关的关键基因模块，而且这些关键基因在随机森林预测模型中得到验证。

1、材料和方法

1.1 数据集选择与预处理

从NCBI GEO中检索获得子痫前期患者的胎盘组织RNA-seq数据集，分别为GSE114691、GSE186257、GSE148241和GSE218039。使用SRA-Toolkit软件下载原始数据SRA文件，通过Fastq-dump命令获得FASTQ文件。使用FastQC软件对原始数据进行质量评估，并通过trim＿galore软件进行质量控制。然后利用Hisat2软件将质控后的序列与人类基因组序列(GRCh38)进行比对，使用FeatureCounts软件对基因表达进行定量，最终获得基因表达矩阵。

1.2 统计学处理

使用R 4.2.3软件对测序数据进行统计分析。使用R 4.2.3软件中DESeq2 1.38.3包鉴定子痫前期患者与正常孕妇胎盘组织中的差异表达基因，|log2(FC)|≥1(其中FC为差异倍数)且校准P<0.05被认为是差异表达的基因。使用R 4.2.3软件randomForest 4.7-1.1包建立预测子痫前期的随机森林模型，训练集与测试集分别为每个数据集中70%和30%的数据，通过Predict函数计算模型的预测概率分值，并使用R 4.2.3软件pROC 1.18.2包绘制模型的ROC曲线，计算AUC值以衡量模型的预测性能。

1.3 基因富集分析

利用注释、可视化和集成发现数据库(Database for Annotation，Visualization，andIntegrated Discovery；DAVID对差异表达基因进行基因本体(Gene Ontology，GO)分析，确定基因富集的生物学过程、分子功能和细胞组分[15-16]。在KOBAS 3.0网站进行京都基因与基因组百科全书(Kyoto Encyclopedia of Genes and Genomes，KEGG)通路富集分析[16-17]。通过R 4.2.3软件ggplot2包将P<0.05的通路可视化于气泡图中。

1.4 WGCNA共表达网络分析

基因表达数据通过方差稳定变换(variance stabilizing transformation，VST)及对数转换后，利用R 4.2.3软件WGCNA 1.72-1包构建共表达网络。具有相似表达模式的差异表达基因被归入一个模块，每个模块被赋予一种颜色。通过计算基因的模块隶属度(modulemembership，MM)和基因显著性(gene significance，GS)识别与临床表型相关的关键基因[18]。然后利用modulePreservation函数计算保守性Z-summary得分，以剔除保守性较差的模块。最后使用Cytoscape 3.10.0软件对关键基因网络进行可视化。

2、结果

2.1 研究队列及基因表达分析

通过检索NCBI GEO共收集到4个与子痫前期相关的数据集，其中2个数据集明确了样本的亚型分别为重度子痫前期和早发型重度子痫前期(表1)。4个数据集共包括156例孕妇胎盘组织样本的RNA-seq数据，其中70例为子痫前期患者样本、86例为正常孕妇样本。主成分分析结果显示，子痫前期患者与正常孕妇的胎盘组织样本可分为2个明显的聚类，说明子痫前期患者与正常孕妇胎盘组织间的基因表达差异显著(图1A)；而且在所有数据集中显著上调的差异表达基因数量均多于显著下调的差异表达基因(图1B)。

表14个子痫前期相关数据集的样本信息

2.2 差异表达基因的通路富集分析

有49个差异表达基因是4个数据集共有的，其中44个表达上调、5个表达下调(表2、图2A)。GO功能分析结果显示，在生物学过程中这些差异表达基因显著富集在卵泡刺激素分泌的正向调节通路；在细胞组分中，有10个基因[蛋白C受体(protein C receptor，PROCR)、嗅觉调节素样蛋白3(olfactomedin like 3，OLFML3)、瘦素(leptin，LEP)、谷氨酰胺酰肽环转移酶(glutaminyl-peptide cyclotransferase，QPCT)、黄体生成素亚基β(luteinizing hormone subunit β，LHB)、HtrA丝氨酸肽酶4(HtrA serinepeptidase 4，HTRA4)、抑制素亚基α(inhibin subunit α，INHA)、卵泡抑素样蛋白3(follistatin like 3，FSTL3)、冠毛素2(pappalysin 2，PAPPA2)和糖基磷脂酰肌醇锚定高密度脂蛋白结合蛋白1(glycosylphosphatidylinositol anchored high density lipoprotein binding protein 1，GPIHBP1)]富集在细胞外区域，除OLFML3外其余9个基因均表达上调；在分子功能中，差异表达基因显著富集在激素活性通路(图2B)。KEGG通路富集分析结果显示，最显著的是细胞因子-细胞因子受体相互作用通路，C-X3-C基序趋化因子受体1(C-X3-C motif chemokine receptor 1，CX3CR1)、LEP和INHA均富集在该通路上(图2C)。这些结果说明，在子痫前期患者与正常孕妇胎盘组织之间的差异表达基因可能与激素分泌、免疫应答与免疫调节等有关。

图14个数据集的主成分分析(A)和差异表达基因的火山图(B)

2.3 加权基因共表达网络的建立

为进一步探讨这49个差异表达基因与疾病表型的相关性，选择疾病样本量最大的数据集GSE186257进行WGCNA，挖掘可能与子痫前期发生和发展密切相关的共表达模块。共构建了7个共表达模块，这些模块均独立于其他模块，可见棕色和蓝色模块与胎儿性别相关(均P=0.01，图3A)，并且在这2个模块中INHA、长基因间非蛋白编码RNA 2009(long intergenic non-protein coding RNA 2009，LINC02009)和MIR31宿主基因(MIR31 host gene，MIR31HG)是与胎儿性别相关的重要基因(均GS>0.5，图3B)。然而，这7个共表达模块与患有子痫前期的孕妇是否早产并不相关(图3A)。此外，红色模块和蓝色模块之间具有强相关性(图3C)，并且这2个模块内的8个基因[fms相关受体酪氨酸激酶1(fms related receptor tyrosine kinase 1，FLT1)、INHA、肌球蛋白ⅦB(myosin ⅦB，MYO7B)、蛋白磷酸酶1调节抑制因子亚基1C(protein phosphatase 1 regulatory inhibitor subunit 1C，PPP1R1C)、LINC02009、PAPPA2、FSTL3和LEP]与这2个模块高度相关(均MM>0.8，图3D)。

表24个子痫前期相关数据集共有的49个差异表达基因的表达趋势及差异倍数

图24个子痫前期相关数据集中DEG的维恩图(A)及GO(B)、KEGG(C)富集分析结果

图3子痫前期相关数据集GSE186257加权基因共表达网络的建立

2.4 关键基因模块的验证

为了评估GSE186257网络模块在其他数据集上的表达模式，使用modulePreservation函数计算模块的保守性Z-summary得分(图4A)。蓝色模块在其他3个数据集(GSE114691、GSE148241、GSE218039)均得分较高，表明该基因模块在不同数据集中均表现稳定。然而，Z-summary得分较低的红色和绿色模块在不同数据集中保守性较低，表明这2个基因模块可能能够较好地区分子痫前期的不同亚型。基于以上结果，推测蓝色模块是与子痫前期高度相关的关键模块。蓝色模块内的12个基因中有4个基因(FLT1、PAPPA2、INHA和MYO7B)表达量相对较高(图4B)，并且网络拓扑图显示这4个基因之间存在强相关性(图4C)，提示FLT1、PAPPA2、INHA和MYO7B这4个基因可能是子痫前期网络中的关键基因。

图4数据集GSE186257的共表达模块在数据集GSE114691、GSE148241和GSE218039中的验证

2.5 关键基因预测子痫前期模型的建立

在基因模块中，与模块高度相关的同时与性状也高度相关的基因被定义为该模块的关键基因。通过绘制MM与GS散点图筛选出蓝色模块的6个关键基因，分别为FLT1、PAPPA2、PPP1R1C、MYO7B、LINC02009和INHA(MM>0.8、GS>0.3，图5A)。将这6个基因纳入随机森林模型，ROC曲线分析结果显示，在测试集中，所建立的模型对子痫前期有较好的预测价值(AUC=0.978，图5B)。以上结果表明，基于FLT1、PAPPA2、PPP1R1C、MYO7B、LINC02009和INHA这6个基因所建立的预测模型有助于诊断子痫前期。

3、讨论

子痫前期是导致孕产妇和胎儿发病与死亡的主要原因之一[19]，当前唯一有效治疗子痫前期的方法是终止妊娠，因此开发早期预测及诊断子痫前期的标志物至关重要[20]。胎盘功能不全被认为是导致孕妇患子痫前期的根本原因，然而源自胎盘的分子机制在很大程度上仍不为人所知。RNA-seq技术的发展使许多基因的表达得以测量，通过数据挖掘方法确定正常组织和患者组织之间的差异表达基因有助于了解疾病的发病机制[21]。然而，在目前的许多研究中数据分散及数据量较小的问题普遍存在。因此，本研究整合并深度挖掘了4个RNA-seq公共数据集，获得了更全面的生物学信息，有助于理解子痫前期的发生、发展机制。一般的数据挖掘方法缺乏对大规模的高维数据的系统性分析，WGCNA通过基因之间的相关系数构建分层聚类树，根据聚类树的不同分支将大量基因分为不同的基因模块，从而评估基因模块与临床特征间的关联，该方法在筛选疾病特征标志物和潜在靶点中表现出较其他分析方法更明显的优势[22]。机器学习算法在数据挖掘中的应用也为探究疾病潜在的治疗靶点提供了支撑[23]。本研究通过结合WGCNA与机器学习算法挖掘潜在的子痫前期诊断标志物，并证实这些基因对子痫前期的诊断价值，同时证明本研究所使用挖掘方法的可靠性。

图5子痫前期关键基因的筛选及随机森林模型的预测效能评价

本研究鉴定出的49个共有差异表达基因显著富集在卵泡刺激素分泌的正向调节通路。卵泡刺激素作用于其受体在刺激卵泡发育和成熟中起着关键作用[24]，并且在子痫前期患者的胎盘样本中卵泡刺激素受体mRNA的表达水平显著低于正常孕妇的胎盘样本[25]。细胞外区域是最显著富集且基因数量最多的细胞组分。细胞外区域对细胞维持生理功能发挥着至关重要的作用，包括细胞间信号传递、细胞黏附、细胞外基质形成等，已有研究表明存在于细胞外区域的基因如高迁移率族蛋白B1(high mobility group box 1，HMGB1)具有促炎作用且与子痫前期相关[26]。在本研究中显著富集在细胞外区域的PROCR、HTRA4和LHB已在子痫前期相关研究中被报道。PROCR是一种跨膜糖蛋白，其在胎盘滋养层细胞中表达下调且与子痫前期有关[27]。HTRA4和LHB在子痫前期患者的胎盘绒毛组织中均明显上调，而且已有研究证实相较于健康人群，LHB在子痫前期患者的血浆中显著上调[28]。上述结果提示激素分泌和免疫反应可能在子痫前期发病机制中起着重要作用。KEGG通路富集分析结果显示，差异表达基因显著富集在细胞因子-细胞因子受体相互作用等信号通路，富集在该通路上的基因CX3CR1是一种促血管生成因子，CX3CR1表达可能与子痫前期胎盘血管网发育不全有关[29]。这一结果表明血管生成功能障碍可能是子痫前期的发病机制之一[30]。

为了进一步探索与子痫前期高度相关的关键基因，本研究通过WGCNA将差异表达基因分为7个共表达模块，结果提示蓝色与棕色模块中的INHA、LINC02009和MIR31HG与胎儿性别相关。研究表明孕早期维生素D水平对胎儿性别有一定影响[31]，Murata等[32]发现在成年雌鼠的垂体前叶中INHA和维生素D具有强正相关性，提示INHA与子痫前期患者的胎儿性别有一定关联。本研究结果还提示，蓝色和红色模块中的FSTL3、INHA和LEP与模块高度相关。先前的研究表明，FSTL3和INHA在子痫前期孕妇胎盘样本中明显上调，被定义为枢纽基因，且在预测子痫前期的logistic回归模型中表现优秀[33]。此外，INHA被认为是母体子痫前期的易感基因，可能通过高度失调的免疫和炎症反应促进子痫前期的发展[34]。LEP也被多项研究证明是与子痫前期发病相关的关键基因[35-37]。蓝色模块进一步被鉴定为与子痫前期高度相关的关键模块，FLT1、PAPPA2、PPP1R1C、MYO7B、LINC02009和INHA是与模块高度相关且与性状显著相关的关键基因，其中FLT1、PAPPA2和INHA在所有样本中的表达量相对较高。这些基因之间具有强相互作用，基于此构建的子痫前期的随机森林模型表现优秀(AUC=0.978)。FLT1、PAPPA2和INHA已被研究证实与子痫前期密切相关。FLT1是一种血管内皮生长因子受体(vascular endothelial growth factor receptor，VEGFR)，研究报道VEGFR基因与血管生长和内皮功能障碍有关，这可能在一定程度上解释了子痫前期的发生[38]。多项研究显示，可溶性FLT1在子痫前期早筛和诊断中具有重要价值[39-40]。PAPPA2编码蛋白质分裂形成的胰岛素样生长因子结合蛋白5在子痫前期患者的血浆中显著上调[41-43]。胰岛素样生长因子在刺激绒毛外滋养层侵袭和子宫螺旋动脉重塑过程中发挥着重要作用，PAPPA2在子痫前期患者的血浆中表达增加被认为与子宫胎盘缺血有关[44]。一些血清学研究也证明了PAPPA2早期预测子痫前期的潜在价值[45-46]。上述结果表明，本研究所使用的方法挖掘出的大量靶标与子痫前期的发生、发展高度相关。因此，尽管目前仅在少数研究中涉及的MYO7B[7]，以及尚未有明确研究报道的PPP1R1C和LINC02009同样值得关注，它们也可能参与了子痫前期的发生，并对患者的早期诊断及治疗有潜在价值。

综上所述，本研究利用来自NCBI GEO的子痫前期患者胎盘组织的RNA-seq数据，从4个数据集中初步筛选出49个共有差异表达基因，并明确它们的功能富集途径与激素分泌和免疫反应等相关；进一步通过WGCNA筛选出共表达网络中的6个关键基因(FLT1、PAPPA2、PPP1R1C、MYO7B、LINC02009和INHA)，并在随机森林模型中证实了它们作为子痫前期早期筛查与诊断分子标志物的潜力。这些包括血管生成因子、妊娠相关血浆蛋白、抑制素在内的分子或许能够成为子痫前期早期筛查、诊断的候选标志物，一些新发现的关键基因也可能为子痫前期的治疗提供新的靶点。本研究同时证实WGCNA方法能够显著增强数据挖掘的系统性并提高效率，在筛选疾病诊断标志物和治疗靶点方面的具有较高的优越性和可靠性。

本研究具有一定的局限性，一是本研究的4个数据集之间存在一定异质性，尽管本研究选择了4个数据集的共有差异表达基因，但在后续分析中仍应继续扩大样本量，增加结果的可靠性和可重复性；二是本研究筛选出的候选基因在未来需要更多的工作来验证它们的实际临床应用价值，并探索它们在子痫前期发病相关通路中的调控作用。

参考文献:

[5]赫英东,陈倩.阿司匹林预防子痫前期的局限性和临床应用选择[J].中国实用妇科与产科杂志,2021,37(5):519-522.

[16]宋英娜,杨剑秋,刘俊涛,等.早发型重度子痫前期孕妇胎盘组织中差异表达基因的研究[J].中华妇产科杂志,2014,49 (7):501-505.

基金资助:国家自然科学基金面上项目(81971402)~~;

文章来源:姚瑞倩,喻东,薛赓.基于加权基因共表达网络分析挖掘子痫前期的诊断标志物[J].海军军医大学学报,2024,45(12):1529-1539.