本研究探讨了系统性红斑狼疮与复发性流产之间的关系,并确定IFI27是值得进一步研究的候选生物标志物。
研究文章
* These authors contributed equally
本研究探讨了系统性红斑狼疮与复发性流产之间的关系,并确定IFI27是值得进一步研究的候选生物标志物。
系统性红斑狼疮(SLE)与不良妊娠结局相关,但其与复发性流产(RPL)之间的因果关系及共享的分子特征尚不明确。本研究整合双向双样本孟德尔随机化(MR)与转录组生物信息学分析,以探讨这一关系并鉴定潜在的共享生物标志物。选择FinnGen和英国生物样本库(UK Biobank),因其提供了大规模、无重叠、欧洲血统的全基因组关联研究(GWAS)汇总统计数据。差异表达基因(DEGs)从GSE61635(血液;|log₂ 倍数变化 > 1)和 GSE165004(子宫内膜;|log₂ 倍数变化 > 0.5)使用校正后的 P < 0.05,随后进行功能富集分析、蛋白质-蛋白质相互作用(PPI)分析、枢纽基因筛选、最小绝对收缩与选择算子(LASSO)回归分析、利用GSE50772和GSE198700数据集进行外部验证、受试者工作特征(ROC)分析以及单样本基因集富集分析(ssGSEA)。遗传预测的系统性红斑狼疮(SLE)与自然流产风险的统计学显著增加相关,但效应量较小(逆方差加权法[IVW]比值比[OR] = 1.01,95%置信区间[CI] = 1.00–1.02; P < 0.001)。工具变量强度充足,敏感性分析未发现显著异质性、方向性多效性或具有显著影响的单个变异体。59个共同差异表达基因(DEGs)在抗病毒免疫应答、细胞黏附及凋亡相关过程方面显著富集。IFI27在系统性红斑狼疮(SLE)患者血液中持续高表达,但在复发性流产(RPL)患者的子宫内膜和绒毛组织中表达下调;而CXCL11缺乏一致的外部验证结果。回顾性ROC分析显示,SLE和RPL的曲线下面积(AUC)分别为0.822和0.872。通过计算推断的单样本基因集富集分析(ssGSEA)评分显示,IFI27表达水平与多种免疫细胞特征谱呈相关性,包括辅助性T细胞2(Th2细胞)。这些结果提示IFI27可能是SLE与RPL共有的候选生物标志物;然而,仍需前瞻性临床研究和实验研究进一步验证其生物学与临床意义。
系统性红斑狼疮(SLE)是一种复杂的自身免疫性疾病,其特征是多系统受累和慢性免疫失调1。SLE的病理异常主要归因于适应性免疫反应受损以及抗原-抗体复合物的沉积,从而导致自身免疫介导的组织损伤和器官损害2,3。全球SLE发病率约为每10万人年5.14(1.4–15.13)例,其中女性发病率估计为每10万人年8.82(2.4–25.99)例4。SLE可发生于所有年龄人群,但主要见于育龄期女性5,6。患有SLE的孕妇发生不良妊娠结局的风险增加,包括复发性流产、死产、早产和宫内生长受限7,8。复发性妊娠丢失(RPL)定义为在妊娠20–24周之前发生两次或以上流产9,其报告患病率约为2.6%10,是一种具有临床意义的生殖并发症。约20%的SLE孕妇会发生流产11,SLE被认为是RPL的重要危险因素12。提出的机制包括激素水平改变和免疫失调。抗心磷脂抗体和狼疮抗凝物等生物标志物已被研究作为SLE患者不良妊娠结局的潜在预测因子11。这些自身抗体可结合胎盘滋养层细胞,改变滋养层信号传导、增殖和侵袭能力,调节激素和细胞因子分泌,并增加细胞凋亡,从而导致妊娠结局不良13。此外,β2-糖蛋白I(β2-GPI)作为抗磷脂综合征中的主要抗原,在胎盘组织中表达。抗-β2-GPI抗体与β2-GPI结合会抑制滋养层细胞的生长和分化,导致胎盘发育缺陷。这种相互作用还可促进促炎性环境的形成,表现为破坏性细胞因子的产生和补体激活,从而导致胎盘血栓形成和复发性流产14,15。然而,以往的研究往往缺乏对局部生殖组织(如蜕膜)的全面分析,限制了将全身性生物标志物与局部病理变化相关联的能力。此外,SLE相关妊娠的管理以及不良妊娠结局的预防仍具挑战性。遗传易感性在SLE的发生中起重要作用,遗传变异也被认为参与RPL的发病机制16,17。然而,SLE与RPL之间是否存在因果关系,以及二者共存背后的分子机制和共享基因,目前仍不明确。
孟德尔随机化(Mendelian randomization, MR)是一种成熟的因果推断方法,利用遗传变异作为工具变量来估计暴露因素对疾病结局的因果效应18。通过利用基因型与表型之间的关系,与传统的观察性研究相比,MR 能够减少混杂偏倚和反向因果关系带来的影响。与此同时,基因组微阵列平台和高通量测序技术的进步使得生物信息学分析能够通过转录组谱型识别候选诊断生物标志物和治疗靶点。整合这些互补的方法,可能通过将因果遗传证据与疾病相关基因表达模式相结合,更全面地理解系统性红斑狼疮(SLE)与复发性流产(RPL)之间的关系。因此,本研究旨在探讨 SLE 与 RPL 之间潜在的因果关系,识别共有的候选生物标志物和生物学通路,并优先筛选出值得未来验证的靶点。为实现上述目标,预先设定了以下分析流程:采用双向 MR 评估因果方向;独立进行差异基因表达分析,随后进行转录组整合;通过蛋白质-蛋白质相互作用(PPI)网络分析和最小绝对收缩与选择算子(least absolute shrinkage and selection operator, LASSO)回归实现生物标志物的优先排序;进行外部表达验证和受试者工作特征(receiver operating characteristic, ROC)分析;以及单样本基因集富集分析(single-sample gene set enrichment analysis, ssGSEA),以评估与免疫细胞特征的关联性。该逐步分析流程在图1中进行了总结。

图1。研究设计与分析流程。
上图展示了利用全基因组关联研究(GWAS)汇总统计数据分析系统性红斑狼疮(SLE)与自然流产次数之间关联的双向双样本孟德尔随机化(MR)分析。内容包括工具变量筛选、连锁不平衡剔除、孟德尔随机化(MR)分析及敏感性分析的总结。下图概述了生物信息学分析流程,包括差异表达分析、共同差异表达基因(DEGs)的识别、功能富集分析、蛋白质-蛋白质相互作用(PPI)网络构建、枢纽基因筛选、最小绝对收缩与选择算子(LASSO)回归、外部验证、受试者工作特征(ROC)分析、单样本基因集富集分析(ssGSEA)以及候选生物标志物IFI27的优先排序。IVW,逆方差加权法;KEGG,京都基因与基因组百科全书;GO,基因本体论。 请点击此处查看该图的放大版本。
本研究无需伦理批准,因为仅涉及对公开可用的、去标识化的全基因组关联研究(GWAS)汇总统计和转录组数据集的二次分析。未招募新的参与者,未采集任何生物样本,也未访问任何个体层面的可识别信息。原始的 FinnGen、英国生物样本库(UK Biobank)和基因表达综合数据库(Gene Expression Omnibus)研究均已报告,其伦理批准和知情同意已根据各自机构、国家以及数据库特定的要求获得。本研究使用的所有数据集均按照适用的数据库使用政策、数据访问条件和伦理指南进行访问和分析。作者未尝试对任何参与者进行重新识别。因此,本次二次分析无需额外的书面知情同意。本研究内容包括两部分:孟德尔随机化(MR)分析和生物信息学分析(图1)。本研究完全基于计算分析,使用了公开可用的汇总水平 GWAS 和转录组数据集。未使用任何湿实验试剂或耗材。
MR 分析
GWAS 汇总统计量的数据来源、获取与预处理:
系统性红斑狼疮的 GWAS 汇总统计量数据来源于 FinnGen 第11次发布(finngen_R11_L12_LUPUS;RRID:SCR_022254),该数据来自一个基于芬兰人群的队列。表型定义采用 ICD-10 代码 L93,共包含 423,818 名参与者,其中 777 例病例和 423,041 名对照。FinnGen 的汇总统计量文件以压缩的制表符分隔格式从 FinnGen 公共数据门户下载,并根据数据集作为暴露或结局变量,使用 TwoSampleMR 软件包中的 read_exposure_data() 或 read_outcome_data() 函数导入 R 环境。保留的字段包括:rsID、染色体、基因组位置、效应等位基因、其他等位基因、效应等位基因频率、β 系数、标准误以及关联性 P 值。
自发性流产数量的汇总统计数据通过IEU OpenGWAS资源(ukb-b-419;RRID:SCR_012815)从英国生物样本库(UK Biobank)获取,共包含78,700名参与者。在正向分析中,使用extract_outcome_data(outcomes = "ukb-b-419", proxies = FALSE)提取所选FinnGen单核苷酸多态性(SNPs)与结局之间的关联。在反向分析中,使用extract_instruments(outcomes = "finngen_R11_L12_LUPUS", clump = FALSE)提取与红斑狼疮相关的SNPs,随后从英国生物样本库的汇总统计文件中提取相应SNP的关联数据。用于双向孟德尔随机化(MR)分析的全基因组关联研究(GWAS)数据集特征总结于表1中。
| 性状 | 样本量 | 祖先群体 | 联盟 | 年份 | 全基因组关联研究数据集标识符 |
| 系统性红斑狼疮 | 4,23,818 | 欧洲 | FinnGen (RRID: SCR_022254) | 2024 | finngen_R11_L12_LUPUS |
| 自然流产次数 | 78,700 | 欧洲 | UK Biobank (RRID: SCR_012815) | 2018 | ukb-b-419 |
表1: 用于双向孟德尔随机化分析的全基因组关联研究(GWAS)汇总统计信息。
该表格总结了公开可用的全基因组关联研究数据集,这些数据集被用作正向和反向孟德尔随机化分析中的暴露和结局数据来源,包括样本量、祖先背景、数据来源、数据发布年份以及数据集标识符。
选择 FinnGen 和英国生物样本库(UK Biobank)是因为它们提供了来自非重叠源人群的大规模、公开可获取的、主要为欧洲血统的数据集,并且包含足够的变异位点覆盖范围,适用于两样本孟德尔随机化(MR)分析。据报告,暴露样本与结局样本之间无重叠。由于仅使用了汇总水平的数据,本研究未访问任何个体水平的基因型数据,未进行额外的个体水平标准化处理,也未排除任何参与者。我们依赖于原始全基因组关联研究(GWAS)联盟所实施的样本水平和变异位点水平的质量控制流程。在本分析过程中,进一步在变异位点水平上进行了额外的质量控制,包括显著性筛选、连锁不平衡剔除(clumping)、等位基因一致性校正、工具变量强度评估以及多效性筛查,具体方法如下所述。
FinnGen 第11版根据GRCh38/hg38报告基因组位置,而IEU OpenGWAS的整合数据集及连锁不平衡参考资源则使用与GRCh37兼容的变异注释。因此,暴露和结局相关的变异位点主要通过稳定的rsID进行匹配,而非染色体位置坐标。未进行直接的跨基因组版本位置比对。在孟德尔随机化(MR)分析前,排除了在不同数据集中缺乏明确rsID或等位基因信息不一致的变异位点。FinnGen第11版的汇总统计量采用GRCh38,而OpenGWAS数据则统一至Build 37所使用的参考序列标准。因此,当整合这两个资源时,基于rsID的匹配至关重要。
孟德尔随机化研究设计:
我们严格遵循STROBE-MR指南(补充文件1)19。采用双向双样本孟德尔随机化(MR)设计,评估遗传预测的系统性红斑狼疮与自然流产次数之间的关联。在正向分析中,将系统性红斑狼疮作为暴露因素,自然流产次数作为结局指标;在反向分析中,暴露与结局互换,并重复完整的工具变量筛选、连锁不平衡剔除、数据整合、因果效应估计及敏感性分析流程。使用单核苷酸多态性(SNPs)作为工具变量(IVs)。整个分析流程按以下顺序进行:获取并格式化全基因组关联研究(GWAS)汇总统计数据;筛选与暴露相关的SNPs;剔除重复或注释不完整的变异位点;进行连锁不平衡剔除;提取相应结局的关联数据;排除与结局直接相关的SNPs;整合暴露与结局的等位基因;计算工具变量强度;筛查潜在的混杂表型;估计因果效应;评估异质性与水平多效性;采用孟德尔随机化多效性残差和离群值检验(MR-PRESSO)进行离群值检测;以及进行留一法和单SNP敏感性分析。所有MR分析均使用R 4.4.2版本(RRID:SCR_001905)、TwoSampleMR 0.6.6版本(RRID:SCR_019010)、MRPRESSO 1.0版本(RRID:SCR_023697)及forestploter 1.1.2版本完成。MR分析基于三个核心假设:第一,相关性假设要求所选SNPs必须与暴露因素显著相关;第二,独立性假设要求所选SNPs不受暴露-结局关联中混杂因素的影响;第三,排他性限制假设要求所选SNPs只能通过暴露因素影响结局20(图1)。
SNP 选择方法:
工具变量的选择按以下顺序进行:(1)选择在 P < 5 × 10−8 水平上与暴露相关的 SNP;当工具变量数量不足时,采用 P < 5 × 10−6;(2)使用 clump_data() 函数进行连锁不平衡剔除,参数设置为 R2 < 0.001 且遗传距离为 10,000 kb;仅在必要时为保留可分析的工具变量集合,放宽标准至 5,000 kb 范围内 R2 < 0.01;(3)使用 P = 5 × 10−5 的阈值滤除与结局显著相关的 SNP;(4)使用 harmonise_data() 函数对暴露和结局的等位基因进行一致性校正,并排除回文或存在歧义的变异;(5)以 F = β2/SE2 计算工具变量强度,并剔除 F < 10 的 SNP;(6)在 PhenoScanner V2 中对保留的 SNP 进行筛查,以识别可能混杂系统性红斑狼疮(SLE)与妊娠丢失关系的表型21。抗磷脂抗体(aPL)可能是 SLE 与自然流产次数的共同风险因素。针对各个候选 SNP,在 PhenoScanner V2 中使用默认的 GWAS 目录检索所有已报道的全基因组关联研究(GWAS)关联结果。显著性阈值设为 P < 1 × 10⁻5,参考基因组版本采用默认的 GRCh37。由于研究人群为欧洲血统,启用基于欧洲参考面板(proxies = "EUR")的代理变异搜索,连锁不平衡阈值设为 R2 > 0.8,搜索窗口为 1,000 kb。其余所有搜索参数均保持默认设置。若 SNP 在预设的混杂因素——抗磷脂抗体(aPL)——上显示显著关联,则视为可能存在多效性,并从最终的工具变量集合中剔除,以尽量减少对孟德尔随机化排他性限制假设的违背。用于正向和反向 MR 分析的工具 SNP 分别列于补充表 1 和 2中。
统计分析:
在完成仪器选择和等位基因一致性校正后,使用 TwoSampleMR 软件包中的 mr() 函数计算因果估计值。分析流程按以下顺序进行。首先,采用四种孟德尔随机化(MR)方法估计总体因果效应:逆方差加权法(IVW)、MR-Egger 回归、加权中位数法和加权众数法。自发性流产数量与系统性红斑狼疮(SLE)的效应估计值以比值比(OR)及其相应的 95% 置信区间和 P 值表示。逆方差加权法(IVW)被指定为主要分析方法,因其在所有纳入的 SNP 均为有效工具变量且不存在水平多效性时具有较高的统计效能。然而,当存在水平多效性时,IVW 的估计结果可能产生偏倚22。MR-Egger 回归主要用于在可能存在水平多效性的情况下评估因果推断23。加权中位数法要求至少 50% 的分析权重来源于有效的工具变量(IVs),该方法在存在异质性但无水平多效性时表现最优24。加权众数法通过识别具有相似因果效应的工具变量簇,并基于最大簇估计效应值25。四种 MR 方法所得的效应估计结果见图 2。第二,使用 mr_heterogeneity() 函数实施的 Cochran's Q 检验评估各 SNP 特异性因果估计值之间的异质性。Q 统计量表示各 SNP 估计值相对于总体因果估计值的加权平方偏差之和。若 Q 检验的 P 值 < 0.05,则认为存在异质性,此时采用随机效应 IVW 模型;若无显著异质性,则采用固定效应 IVW 模型26。第三,使用 MRPRESSO 软件包(RRID:SCR_023697)中的 mr_pleiotropy_test() 函数执行 MR-Egger 截距检验,以评估方向性水平多效性。若截距值在 P < 0.05 水平上显著偏离零,则认为存在方向性水平多效性。第四,使用 mr_presso() 函数执行 MR-PRESSO 程序,以检测具有异常多效性效应的 SNP27。当检测到异常值时,将其剔除,并使用剩余的工具变量重新进行因果分析。采用 MR-PRESSO 全局检验评估总体水平多效性,并通过偏差检验判断剔除异常值是否对因果估计值产生实质性影响。第五,使用 mr_leaveoneout() 函数进行留一法敏感性分析。在该分析中,依次排除每个 SNP,并利用剩余的 SNP 重新计算合并的因果估计值。通过 mr_leaveoneout_plot() 可视化结果,以判断总体关联是否主要由单个工具变量驱动。第六,使用 mr_singlesnp() 函数生成各 SNP 特异性的估计值。这些估计值用于构建漏斗图(通过 mr_funnel_plot()),以直观评估可能由方向性水平多效性引起的不对称性。使用 forestploter(版本 1.1.2)生成汇总森林图,以展示不同 MR 方法所得的效应估计值及其置信区间。正向 MR 散点图、SNP 特异性森林图、留一法分析图和漏斗图分别见补充图 1–4。

图2。双向孟德尔随机化分析结果。
(A)正向孟德尔随机化(MR)分析的森林图,以系统性红斑狼疮(SLE)为暴露因素,自然流产次数为结局。(B)反向MR分析的森林图,以自然流产次数为暴露因素,SLE为结局。效应估计值以比值比(ORs)及其95%置信区间(CIs)表示,采用的方法包括逆方差加权法、MR-Egger法、加权中位数法和加权众数法。SNP:单核苷酸多态性。 请点击此处查看本图的高清版本。
当逆方差加权(IVW)估计值在 P < 0.05 水平上具有统计学显著性,且 MR-Egger、加权中位数和加权众数估计值的方向与 IVW 估计值一致,同时异质性、多效性、MR-PRESSO 或逐一剔除法敏感性分析未对结果产生实质性影响时,认为支持存在因果关联。所有统计检验均为双侧检验。
生物信息学分析
微阵列数据:
转录组数据集来自基因表达综合数据库(Gene Expression Omnibus, GEO;RRID:SCR_005012)28。下载了GSE61635、GSE165004、GSE50772和GSE198700的已处理Series Matrix文件、样本元数据以及平台注释文件。各数据集的平台、组织来源、样本量和分析类别汇总于表2。由于这些数据集来源于不同的组织和微阵列平台,因此每个数据集均被独立进行预处理和分析。不同数据集的表达矩阵未直接合并,也未进行跨平台批次校正。只有在各个发现数据集中独立完成差异表达分析后,才在基因符号水平上进行跨数据集整合。
| GEO 数据集 | 疾病 | 平台 | 组织(智人) | 病例数 | 对照数 | 实验类型 | 贡献者 | 数据集类别 |
| GSE61635 | 系统性红斑狼疮(SLE) | GPL570 | 全血 | 99 | 30 | 表达谱微阵列 | Greidinger EL | 发现数据集 |
| GSE165004 | 复发性流产(RPL) | GPL16699 | 子宫内膜 | 24 | 24 | 表达谱微阵列 | Keleş ID29 | 发现数据集 |
| GSE50772 | 系统性红斑狼疮(SLE) | GPL570 | 外周血单个核细胞(PBMCs) | 61 | 20 | 表达谱微阵列 | Kennedy WP30 | 验证数据集 |
| GSE198700 | 复发性流产(RPL) | GPL13534 | 绒毛膜绒毛 | 5 | 5 | 表达谱微阵列 | Li Y31 | 验证数据集 |
表2: 用于生物信息学分析的转录组学数据集。
该表格总结了纳入发现和验证分析的基因表达综合数据库(GEO)转录组学数据集,包括疾病类型、微阵列平台、组织来源、样本量、实验类型、原始研究贡献者以及数据集类别。
GSE61635 数据集使用 Affymetrix Human Genome U133 Plus 2.0 Array 平台(GPL570)生成,包含来自系统性红斑狼疮(SLE)患者的 99 个全血芯片数据,其中包括部分患者的多次随访样本,以及来自独立健康对照的 30 个芯片数据。已提交的表达矩阵由原始研究者进行了稳健多阵列平均背景校正、分位数归一化、探针集汇总以及 log2 转换。因此,未再进行二次背景校正或分位数归一化处理。患者标识符从 GEO 元数据中提取,并保留在重复测量模型分析中使用。
GSE165004 数据集使用 Agilent SurePrint G3 Human Gene Expression v2 8×60K 微阵列平台(GPL16699)生成。完整数据集包含 24 例生育正常对照、24 例复发性流产(RPL)患者和 24 例原因不明不孕症患者。仅纳入在月经周期第 19–21 天采集的 24 例 RPL 样本和 24 例生育正常对照样本;排除 24 例原因不明不孕症样本,因其超出预定义的比较范围29。采用数据提交者已归一化的表达矩阵,通过箱线图和密度图确认样本分布具有可比性后,未再进行额外的芯片间归一化处理。
GSE50772 被用作独立的系统性红斑狼疮(SLE)验证数据集,包含使用 GPL57030 生成的 61 例 SLE 患者和 20 例健康对照者的外周血单个核细胞样本。GSE198700 使用 GPL13534 生成,包含来自 5 例复发性流产(RPL)患者和 5 例选择性人工流产对照者的绒毛膜绒毛样本31。已提交的表达矩阵被完整导入,并进行一次 log2(x + 1) 转换,因为提交的表达值是以非对数尺度提供的。该转换在进行样本水平质量控制、探针注释、基因水平汇总、候选基因验证、差异表达分析、组间比较检验和 ROC 分析之前,应用于完整的表达矩阵。候选基因未单独进行转换,后续验证分析中也未进行额外的对数转换。对于所有数据集,在分析前均根据相应的 GEO 元数据对样本标识、疾病状态、组织来源和组别标签进行了交叉核对。质量控制包括对文库大小或表达分布、样本箱线图、主成分分析、层次聚类以及样本距离热图的评估。经过质量控制评估后,未排除任何额外样本。
差异表达分析:
使用 limma 版本 3.60.6(RRID:SCR_010943)分别对 GSE61635 和 GSE165004 进行差异表达分析。所有表达矩阵均以基因为行、样本为列进行组织。差异表达的阈值设定为 GSE61635 的 |log₂ 倍数变化| > 1 和 GSE165004 的 |log₂ 倍数变化| > 0.5,并且 Benjamini–Hochberg(BH)校正后的 P < 0.05。火山图使用 ggplot2 版本 3.5.1(RRID:SCR_014601)生成。根据校正后 P 值排序,选取前 50 个最显著的差异表达基因(DEGs)绘制热图,使用 pheatmap 版本 1.0.12(RRID:SCR_016418)完成。通过使用基础 R 的 intersect() 函数对显著的 SLE 和 RPL 差异表达基因列表中的官方基因符号取交集,识别共有差异表达基因,并利用 ggvenn 版本 0.1.16(RRID:SCR_025300)进行可视化。SLE 与 RPL 差异表达基因列表的热图、火山图及交集结果见图 3。

图3.系统性红斑狼疮与复发性流产中差异表达的基因。
(A)GSE61635数据集中系统性红斑狼疮(SLE)患者与健康对照之间50个最显著差异表达基因(DEGs)的热图。(B)GSE165004数据集中复发性流产(RPL)患者与生育正常对照之间50个最显著差异表达基因(DEGs)的热图。(C)GSE61635中差异基因表达的火山图。(D)GSE165004中差异基因表达的火山图。(E)SLE与RPL发现数据集中显著差异表达基因列表重叠的维恩图。DEGs,差异表达基因。 请点击此处查看该图的放大版本。
交集差异表达基因的功能富集分析:
为在分子水平上分析差异表达基因(DEG)的功能,使用 DAVID 在线工具(2021 版;RRID:SCR_001881)32 进行基因本体(GO)功能和京都基因与基因组百科全书(KEGG)通路富集分析。以官方人类基因符号作为标识符类型,并选择 Homo sapiens 为物种。自定义背景基因集由在 GSE61635 和 GSE165004 两个数据集中均通过探针注释和质控且可检测到的所有基因的交集构成。最小基因计数阈值设为 2,最大 EASE 分数(代表 DAVID 修正的单侧 Fisher 精确 P 值)设为 0.05。采用 DAVID 输出结果中“Benjamini”列所提供的 Benjamini–Hochberg 方法对多重比较进行校正。当 EASE 分数 < 0.05 且 Benjamini 校正后的 P 值 < 0.05 时,认为该功能条目具有统计学显著性。完整的 DAVID 输出结果,包括条目名称、基因数量、EASE 分数、Benjamini 校正后的 P 值、输入基因映射及背景基因映射,均以制表符分隔文件格式导出。使用 CNSknowall 网站对筛选后的 DAVID 结果进行可视化展示。GO 和 KEGG 富集分析结果见图 4A。

图4共享差异表达基因的功能富集分析及蛋白质-蛋白质相互作用网络
(A)59个共有差异表达基因(DEGs)的基因本体(GO)与京都基因与基因组百科全书(KEGG)富集分析。桑基图展示了基因与富集GO条目之间的关联关系,配套的气泡图根据富集因子、基因数量和统计显著性总结了富集的GO与KEGG条目。(B)基于59个共有DEGs利用STRING构建并在Cytoscape中可视化的蛋白质-蛋白质相互作用(PPI)网络。节点大小和颜色反映网络连接度,边表示预测的蛋白质-蛋白质关联。BP,生物过程;CC,细胞组分;MF,分子功能。 请点击此处以查看此图的放大版本。
蛋白质相互作用网络及核心基因的鉴定:
将共有差异表达基因(DEGs)上传至 STRING 数据库版本 11.0(RRID:SCR_005223)33,选择物种为智人(Homo sapiens,分类编号:9606)。采用完整的 STRING 网络,包含功能关联和物理相互作用的蛋白关系。启用所有可用的证据通道,包括实验验证数据、人工审编数据库、共表达、文本挖掘、基因邻域、基因融合以及基因共现信息。
最低要求的相互作用得分设定为 0.400,对应中等置信度。未添加任何额外的第一层或第二层互作蛋白;因此,该网络仅包含由提交的共有差异表达基因(DEGs)编码的蛋白质。网络边采用置信度模式显示,并导出为包含相互作用蛋白及其综合 STRING 得分的制表符分隔值文件。STRING 置信度得分表示存在关联关系的可信程度,而非相互作用的强度或结合强度。
将STRING网络文件导入Cytoscape 3.10.0版本(RRID:SCR_003032)34。在进行网络拓扑分析之前,移除了与其他提交蛋白无相互作用的节点35。剩余网络被视为无向网络。STRING综合得分作为边属性保留用于可视化,而cytoHubba的排序则基于默认的非加权拓扑定义生成。所得的蛋白质-蛋白质相互作用(PPI)网络如图4B所示。
使用 cytoHubba 版本 0.1(RRID:SCR_017677)中的六种算法对枢纽基因进行排序:最大团中心性(Maximal Clique Centrality, MCC)、最大邻域组分(Maximum Neighborhood Component, MNC)、边渗透组分(Edge Percolated Component, EPC)、度值(Degree)、接近性(Closeness)和径向性(Radiality)36。对于每种算法,基因按降序排列,并保留排名前 10 的基因。网络枢纽候选基因通过六个前 10 名基因列表的严格交集确定。因此,只有当某个基因在所有六种算法产生的前 10 名基因中均出现时,才被保留为网络枢纽基因。排序结果及交集分析过程已导出并存档。每种 cytoHubba 算法识别出的前 10 位基因列于表 3中。
| 排名 | 最大团中心性 (MCC) | 最大邻域组分 (MNC) | 边渗透组分 (EPC) | 度 | 接近中心性 | 径向性 |
| 1 | RSAD2 | RSAD2 | RSAD2 | RSAD2 | RSAD2 | RSAD2 |
| 2 | RTP4 | RTP4 | RTP4 | RTP4 | RTP4 | RTP4 |
| 3 | IFIT3 | IFIT3 | IFIT3 | IFIT3 | IFIT3 | IFIT3 |
| 4 | IFI27 | IFI27 | IFI27 | IFI27 | IFI27 | IFI27 |
| 5 | IFI44 | IFI44 | IFI44 | IFI44 | IFI44 | IFI44 |
| 6 | GBP1 | GBP1 | GBP1 | GBP1 | GBP1 | GBP1 |
| 7 | MX1 | MX1 | MX1 | MX1 | MX1 | MX1 |
| 8 | OAS1 | OAS1 | OAS1 | OAS1 | OAS1 | OAS1 |
| 9 | IFIT1 | IFIT1 | IFIT1 | IFIT1 | IFIT1 | IFIT1 |
| 10 | CXCL11 | CXCL11 | CXCL11 | CXCL11 | CXCL11 | CXCL11 |
表3: 通过六种cytoHubba排序算法识别出的前10个枢纽基因。
使用Cytoscape软件中cytoHubba插件实现的六种网络拓扑算法对共有差异表达基因进行排序。每种算法得出的排名前10位的基因分别列出,以便在最大团中心性(MCC)、最大邻域组分(MNC)、边渗透组分(EPC)、度(Degree)、接近性(Closeness)和径向性(Radiality)之间进行比较。
使用LASSO回归识别核心基因:
在系统性红斑狼疮(SLE)和复发性流产(RPL)的发现数据集中,分别使用glmnet 4.1-8版本(RRID:SCR_015505)独立进行LASSO逻辑回归分析。预测变量矩阵由网络枢纽候选基因的标准化表达值构成,其中行代表样本,列代表基因。疾病状态编码为1,对照状态编码为0。采用family = "binomial"和alpha = 1,拟合带有纯LASSO惩罚项的二项分布广义线性模型。预测变量通过standardize = TRUE在内部进行标准化处理,并包含截距项。使用自定义的基础R代码,分别为SLE和RPL数据集生成按类别分层的10折交叉验证分配方案。在每个疾病状态分层内,样本索引被随机置换,并尽可能均匀地分配到10个折中,分配方式为sample(rep(seq_len(10), length.out = n))。在生成每个数据集的折分配方案前,设定随机种子为123,以确保结果可重复。由于每个疾病状态组均包含超过10个样本,因此每次交叉验证的折中均包含病例和对照样本。生成的整数向量(foldid_sle和foldid_rpl)被提供给cv.glmnet()函数的foldid参数,并在同一数据集内对所有评估的λ值使用相同的折分配方案。
模型拟合时采用 family = "binomial"、alpha = 1、nfolds = 10、type.measure = "deviance"、standardize = TRUE、intercept = TRUE、nlambda = 100、thresh = 1 × 10⁻7 和 maxit = 100000。惩罚参数通过 lambda.min 选择,该值定义为产生最小平均交叉验证二项偏差的 lambda。更为保守的 lambda.1se 定义为在最小交叉验证误差一个标准误范围内的最大 lambda,作为敏感性分析结果予以记录。LASSO 方法分别独立应用于 GSE61635 和 GSE165004 数据集。在两种疾病特异性模型中系数均非零的基因被定义为共有的 LASSO 筛选候选基因。通过引入 L1 正则化项,该方法有效将信息量较低基因的系数压缩至零,从而实现特征筛选37。系统性红斑狼疮(SLE)和复发性流产(RPL)发现数据集的系数图谱及十折交叉验证曲线见图5。

图5网络枢纽基因的最小绝对收缩与选择算子回归分析。
(A)系统性红斑狼疮(SLE)发现数据集(GSE61635)通过最小绝对收缩与选择算子(LASSO)逻辑回归生成的系数图谱。(B)用于确定最优惩罚参数的十折交叉验证曲线λ) 用于系统性红斑狼疮(SLE)模型。(C)基于复发性流产(RPL)发现数据集(GSE165004)的LASSO逻辑回归生成的系数图谱。(D)用于确定最优惩罚参数的十折交叉验证曲线λ用于RPL模型。上方x轴上的数值表示在每个λ值下非零回归系数的数量 λ垂直虚线表示 λ分钟 λ_1se 请点击此处以查看该图的放大版本。
核心基因诊断价值的验证:
在独立的系统性红斑狼疮(SLE)数据集GSE50772和独立的狼疮性肾炎(RPL)数据集GSE198700中评估了LASSO筛选出的候选基因的表达模式。仅在GSE61635和GSE165004数据集中完成候选基因筛选后,才使用外部数据集。在验证数据集中未进行额外的特征筛选或模型拟合。使用双侧Wilcoxon秩和检验比较病例组与对照组之间的候选基因表达水平。当在同一数据集中检测多个候选基因时,采用Benjamini–Hochberg方法对所得的P值进行多重检验校正。若某候选基因在校正后其病例组与对照组间的表达差异具有统计学显著性,且其变化方向与相应的发现数据集一致,则认为该基因在外部数据集中得到成功复制。候选基因在发现数据集与验证数据集中的表达模式详见图6。

图6发现和验证数据集中 IFI27 和 CXCL11 的表达
(A,B) 系统性红斑狼疮(SLE)发现数据集(GSE61635)中 IFI27 和 CXCL11 的表达情况。(C,D) 独立 SLE 验证数据集(GSE50772)中 IFI27 和 CXCL11 的表达情况。(E,F) 复发性流产(RPL)发现数据集(GSE165004)中 IFI27 和 CXCL11 的表达情况。(G) 独立 RPL 验证数据集(GSE198700)中 IFI27 的表达情况。采用双侧 Wilcoxon 秩和检验比较各组间的基因表达水平。 P 当在同一数据集中检测多个候选基因时,采用Benjamini–Hochberg方法对数值进行校正。 P < 0.05; **** P < 0.0001;ns,无显著性差异。 请点击此处查看此图的放大版本。
采用 pROC 版本 1.18.5(RRID:SCR_024286)38 进行受试者工作特征(ROC)分析。在每个发现数据集中,针对每个候选基因分别生成 ROC 曲线。使用 DeLong 法计算 ROC 曲线下面积(AUC)及其双侧 95% 置信区间。探索性诊断阈值通过最大 Youden 指数确定。阈值、敏感性和特异性的置信区间通过设置随机种子为 123 的 2,000 次分层自助抽样法计算。AUC 被用作一种与阈值无关的区分能力度量39。由于这些数据集为回顾性数据,且使用不同的组织、平台和标准化方法生成,因此 Youden 指数导出的阈值在每个数据集中独立计算,并视为探索性的数据集特异性阈值。这些阈值不被视为标准化的临床阈值,也不在不同平台之间直接转移。外部 ROC 结果代表转录组学验证,而非前瞻性临床验证。pROC 通过 coords() 函数支持 AUC 的 DeLong 置信区间计算和 Youden 指数优化,同时可通过分层自助重采样法估计 ROC 坐标点的置信区间。SLE 和 RPL 发现数据集中候选基因区分能力的 ROC 曲线及汇总结果见图 7。

图7.候选基因的受试者工作特征分析。
(A)干扰素诱导蛋白27(IFI27)在系统性红斑狼疮(SLE)发现数据集中的受试者工作特征(ROC)曲线。(B)趋化因子(C-X-C基序)配体11(CXCL11)在SLE发现数据集中的ROC曲线。(C)IFI27和CXCL11在SLE发现数据集中的诊断效能汇总。(D)IFI27在复发性妊娠丢失(RPL)发现数据集中的ROC曲线。(E)CXCL11在RPL发现数据集中的ROC曲线。(F)IFI27和CXCL11在RPL发现数据集中的诊断效能汇总。曲线下面积(AUC)值以95%置信区间(CIs)表示。 请点击此处查看该图的高清版本。
ssGSEA免疫浸润分析:
鉴于免疫细胞失调在系统性红斑狼疮(SLE)和复发性流产(RPL)发病机制中的作用40,41,本研究在GSE61635和GSE165004发现数据集中通过计算方法推断了免疫细胞富集情况。分析在每个数据集中独立进行,未对数据集进行合并。免疫细胞基因特征集合包含了Charoentong et al.42所描述的28种免疫细胞群体的标志基因集。原始补充性基因特征表被转换为使用官方人类基因符号命名的基因集列表。每个基因集内的重复基因符号被去除。在对应表达矩阵中缺失的基因被剔除,且在标识符映射后包含少于五个匹配基因的基因集也被排除在该数据集分析之外。单样本基因集富集分析(ssGSEA)使用GSVA 1.52.3版本(RRID:SCR_021058)和GSEABase 1.66.0版本完成。在GSVA 1.52.3版本中,需提供方法特异的参数对象。所使用的参数如下:minSize = 5,maxSize = 500,alpha = 0.25,normalize = TRUE,checkNA = "yes”。
alpha 参数设置为 0.25,并启用了最终的 ssGSEA 评分归一化。在与表达矩阵匹配后,基因集被限制为包含 5 至 500 个基因。采用单线程执行方式,以确保在不同系统间计算结果的一致性。未使用 kcdf 参数,因为它不是 GSVA 1.52.3 版本中 ssgseaParam() 过程的参数。ssGSEA 生成的是相对的样本水平基因集富集评分,而非实验测定的免疫细胞计数或绝对细胞比例43。GSVA 1.52.3 工作流程需要一个方法特定的参数对象,ssGSEA 的参数包括 alpha、评分归一化和基因集大小限制。
针对每种免疫细胞特征,使用双侧 Wilcoxon 秩和检验比较疾病组与对照组之间的 ssGSEA 评分。28 种细胞类型比较的 P 值在每个数据集中分别采用 Benjamini–Hochberg 方法进行校正。校正后 P 值 < 0.05 的免疫细胞特征被视为富集差异显著。
计算了每个数据集中候选基因表达与每种免疫细胞特征的ssGSEA评分之间的Spearman秩相关系数。相关性 P 采用 Benjamini–Hochberg 方法对数据集中所有候选基因–免疫细胞组合进行多重检验校正。经校正后,相关性在调整后的显著性水平下被认为具有统计学意义 P 值 < 0.05。使用 ggcorrplot 版本 0.1.4.1 可视化相关性矩阵,并使用 ggplot2 版本 3.5.1 生成组间比较图。
所有统计检验均使用原始归一化的ssGSEA评分。热图和堆叠可视化仅用于描述性展示。这些评分并不表示免疫细胞的直接比例,所观察到的关联性应解释为计算得出的相关性,而非经实验验证的细胞-基因相互作用。免疫细胞特征富集谱、组间比较以及与候选基因表达的相关性结果见图8。

图8系统性红斑狼疮与复发性流产中共有候选基因的免疫细胞特征富集及其相关性分析
(A)系统性红斑狼疮(SLE)发现队列中28种免疫细胞特征的单样本基因集富集分析(ssGSEA)得分的层次聚类热图。(B)SLE患者与健康对照者之间免疫细胞特征ssGSEA得分的比较。(C)Spearman相关性热图,显示SLE发现队列中IFI27和CXCL11表达水平与28种免疫细胞特征ssGSEA得分之间的关联。(D)复发性流产(RPL)发现队列中28种免疫细胞特征的ssGSEA得分的层次聚类热图。(E)RPL患者与生育正常对照者之间免疫细胞特征ssGSEA得分的比较。(F)Spearman相关性热图,显示RPL发现队列中IFI27和CXCL11表达水平与28种免疫细胞特征ssGSEA得分之间的关联。相关性采用Spearman秩相关分析计算, P 使用Benjamini–Hochberg方法对数值进行校正。 P < 0.05; ** P < 0.01; *** P < 0.001;ns,无显著性差异。 请点击此处以查看此图的放大版本。
MR 分析
在工具变量筛选和数据整合后,保留了 16 个 SNP 用于正向 MR 分析,其中系统性红斑狼疮(SLE)被视为暴露因素,自然流产次数被视为结局。工具变量的详细信息见补充表 1。所有保留的 SNP 的 F 统计量均大于 10,表明弱工具变量偏倚的可能性较小。每个保留的 SNP 还通过 PhenoScanner V2 进行了筛查,未发现与抗磷脂抗体(aPL)相关的 SNP。MR-PRESSO 分析未检测到异常值。Cochran’s Q 检验显示 SNP 特异性估计值之间无显著异质性(Q = 16.12,P = 0.31);因此采用固定效应 IVW 模型。MR-Egger 截距检验未提示存在方向性水平多效性(P = 0.69)。IVW 分析显示,遗传预测的 SLE 与自然流产次数之间存在具有统计学意义但效应量较轻微的正向关联(比值比 [OR] = 1.01,95% 可信区间 [CI] = 1.00–1.02,P < 0.01;图 2A)。使用 MR-Egger 回归(OR = 1.01,95% CI = 1.00–1.03,P = 0.16)、加权中位数法(OR = 1.01,95% CI = 1.00–1.02,P = 0.17)和加权众数法(OR = 1.01,95% CI = 0.99–1.03,P = 0.42)获得的效应估计值与 IVW 估计值方向一致,尽管各自未达到统计学显著性。逐一剔除分析显示,剔除任一单个 SNP 均未对合并估计值产生实质性影响,且大致对称的漏斗图未提供视觉证据表明结果受显著的方向性多效性驱动。相应的散点图、SNP 特异性森林图、逐一剔除分析图和漏斗图见补充图 1–4。
在反向 MR 分析中,经过工具变量筛选后保留了 16 个 SNP,所有 SNP 的 F 统计量均大于 10(补充表 2)。MR-PRESSO 分析未发现任何离群值。Cochran’s Q 检验显示无显著异质性(Q = 13.41,P = 0.50),MR-Egger 截距检验也未提示存在定向水平多效性(P = 0.41)。IVW 估计结果不支持遗传预测的自然流产次数与系统性红斑狼疮(SLE)风险之间存在关联(OR = 0.93,95% CI = 0.21–4.23,P = 0.93;图 2B)。综上,MR 结果支持在正向分析中存在较弱的关联,即从遗传预测的 SLE 指向自然流产次数,而反向分析未支持从遗传预测的自然流产次数指向 SLE 风险的关联。
生物信息学分析
差异表达分析:
对GSE61635数据集的差异表达分析共鉴定出976个系统性红斑狼疮(SLE)组与健康对照组之间的差异表达基因(DEGs),其中包括678个上调基因和298个下调基因(图3C)。对GSE165004数据集的分析共鉴定出1,249个复发性流产(RPL)组与对照组之间的差异表达基因,其中包括578个上调基因和671个下调基因(图3D)。两个发现数据集中最显著的50个差异表达基因的热图分别展示于图3A和图3B。此外,在两个数据集中共鉴定出59个共同差异表达基因(图3E)。这些共同差异表达基因构成了后续功能富集分析和网络分析所使用的基因集合。
交集差异表达基因的功能富集分析:
使用 DAVID 对 59 个共有差异表达基因(DEGs)进行 GO 和 KEGG 通路富集分析。在生物过程类别中,这些共有差异表达基因显著富集于抗病毒防御反应、对病毒的反应、病毒基因组复制的负调控、抗病毒先天免疫反应、细胞凋亡过程的负调控以及细胞黏附等过程。在细胞组分方面,富集的条目包括细胞外区域、内质网膜、肌动蛋白细胞骨架和膜。在分子功能方面,富集的条目包含钙离子结合。KEGG 分析显示,这些基因在丙型肝炎和甲型流感相关通路中显著富集(图 4A)。这些结果表明,共有差异表达基因主要参与抗病毒和免疫相关的生物过程,为系统性红斑狼疮(SLE)和复发性种植失败(RPL)发现数据集中共有的基因提供了功能背景。
蛋白质-蛋白质相互作用网络及枢纽基因的鉴定
将59个共有差异表达基因(DEG)上传至STRING数据库,以0.400的最低相互作用置信度评分构建蛋白质-蛋白质相互作用(PPI)网络。所得网络包含59个节点和80条边。将该网络导入Cytoscape软件(版本3.10.0)进行可视化分析,并在进行拓扑分析前移除孤立节点(图4B)。使用cytoHubba插件对枢纽基因进行排序,共应用了六种算法:最大团中心性(Maximal Clique Centrality, MCC)、最大邻域成分(Maximum Neighborhood Component, MNC)、边渗透组分(Edge Percolated Component, EPC)、度值(Degree)、接近性(Closeness)和放射性(Radiality)。六种算法均一致识别出排名前十的相同基因:RSAD2、RTP4、IFIT3、IFI27、IFI44、GBP1、MX1、OAS1、IFIT1 和 CXCL11(表3)。因此,这些基因被保留作为后续LASSO回归分析的候选网络枢纽基因。
LASSO 回归分析鉴定出 IFI27 和 CXCL11 为共有的候选基因
在系统性红斑狼疮(SLE)和复发性流产(RPL)的发现数据集中,对10个候选枢纽基因进行了LASSO回归分析。在SLE数据集中,在选定的lambda值下,有四个基因保留了非零系数:IFIT3、IFI27、IFI44 和 CXCL11,其系数分别为 2.575、0.057、2.359 和 0.307(图 5A、B)。在RPL数据集中,四个基因保留了非零系数:IFI27、GBP1、OAS1 和 CXCL11,其系数分别为 −0.897、0.167、−1.007 和 −0.519(图 5C、D)。比较两种疾病特异性模型所筛选的基因后,确定 IFI27 和 CXCL11 为共有的LASSO筛选候选基因。随后,这两个基因在发现数据集和外部验证数据集中进行了进一步评估。
IFI27 和 CXCL11 表达的外部验证
使用从 GEO 数据库获取的独立验证数据集 GSE50772 和 GSE198700 评估了两个候选基因的表达模式。在 GSE61635 数据集中,与健康对照组相比,SLE 组中 IFI27 和 CXCL11 均显著上调(图 6A、B)。在独立的 SLE 验证数据集(GSE50772)中,IFI27 仍显著上调(图 6C),而 CXCL11 在两组间的表达差异无统计学意义(图 6D)。在 RPL 发现数据集(GSE165004)中,与对照组相比,RPL 组中 IFI27 和 CXCL11 均显著下调(图 6E、F)。在独立的 RPL 验证数据集(GSE198700)中,RPL 组中 IFI27 仍显著下调(图 6G),而 CXCL11 未被检测到。总体而言,IFI27 在 SLE 和 RPL 的发现数据集及验证数据集中均表现出一致的差异表达。相比之下,CXCL11 在外部验证数据集中未能得到一致重复。因此,在后续分析中优先选择 IFI27 作为共有的候选生物标志物。
诊断区分能力的探索性评估
采用受试者工作特征(ROC)分析,评估 IFI27 和 CXCL11 表达水平在已分析的回顾性转录组数据集中区分疾病样本与对照样本的能力。对于系统性红斑狼疮(SLE),IFI27 的 ROC 曲线下面积(AUC)为 0.822(95% CI = 0.752–0.892;图 7A),而 CXCL11 的 AUC 为 0.852(95% CI = 0.786–0.917;图 7B)。SLE 数据集中两个候选基因 ROC 曲线的比较见 图 7C。对于复发性流产(RPL),IFI27 的 AUC 为 0.872(95% CI = 0.773–0.970;图 7D),而 CXCL11 的 AUC 为 0.668(95% CI = 0.513–0.882;图 7E)。RPL 数据集中两个候选基因 ROC 曲线的比较见 图 7F。IFI27 在两个疾病数据集中 AUC 值均超过 0.80,且在不同表达数据集中表现出比 CXCL11 更一致的外部验证结果。这些结果支持 IFI27 作为进一步评估的候选生物标志物。然而,由于 ROC 分析使用的是回顾性公共转录组数据集,因此结果应被视为转录组水平区分能力的探索性证据,而非前瞻性临床诊断验证。
免疫浸润的计算评估
采用ssGSEA分析GSE61635和GSE165004发现数据集中28种免疫细胞特征的富集情况。SLE和RPL数据集的免疫细胞富集热图分别如图8A、D所示,相应的ssGSEA评分组间比较结果如图8B、E所示。在SLE数据集中,SLE患者与健康对照之间多种免疫细胞特征存在显著差异,包括代表CD8+ T细胞、CD4+ T细胞、B细胞、树突状细胞、1型辅助性T(Th1)细胞、2型辅助性T(Th2)细胞、17型辅助性T(Th17)细胞、自然杀伤细胞、巨噬细胞、嗜酸性粒细胞、肥大细胞、单核细胞和中性粒细胞的特征(图8B)。在RPL数据集中,RPL组中活化CD8+ T细胞、活化CD4+ T细胞、效应记忆CD4+ T细胞、Th17细胞和单核细胞的ssGSEA评分高于对照组;相反,调节性T细胞(Treg)和巨噬细胞的ssGSEA评分在RPL组中低于对照组(图8E)。相关性分析显示,在SLE数据集中,IFI27和CXCL11的表达水平与活化CD4+ T细胞、自然杀伤细胞、Th2细胞及中央记忆CD8+ T细胞的ssGSEA评分呈正相关,而与Th1细胞的ssGSEA评分呈负相关(图8C)。在RPL数据集中,IFI27的表达水平与Treg和Th2细胞的ssGSEA评分呈正相关,而CXCL11的表达水平与嗜酸性粒细胞的ssGSEA评分呈正相关(图8F)。
数据可用性:
本研究未生成任何新的主要人类受试者数据。所有分析均完全基于公开可用的全基因组关联研究(GWAS)汇总统计和转录组数据集。系统性红斑狼疮(SLE)的GWAS汇总统计数据来自FinnGen第11版(编号:finngen_R11_L12_LUPUS)。自然流产次数的汇总统计数据来自IEU OpenGWAS资源(编号:ukb-b-419),该数据基于英国生物银行(UK Biobank)数据。转录组数据集来自美国国家生物技术信息中心基因表达综合数据库(NCBI Gene Expression Omnibus,GEO),编号分别为GSE61635、GSE165004、GSE50772和GSE198700。这些公开可用的数据集可通过以下数据库获取:
-- FinnGen Release 11: https://r11.finngen.fi/
-- IEU OpenGWAS: https://gwas.mrcieu.ac.uk/
-- 基因表达综合数据库(GEO): https://www.ncbi.nlm.nih.gov/geo/
支持本研究结果的处理后数据已包含在文章及其补充材料中。未获取或保存任何个体层面或可识别个人身份的参与者数据。分析工作流程使用公开可获得的软件和程序包完成,具体如方案部分所述。
补充文件1. 已完成的STROBE-MR报告清单
已完成 使用孟德尔随机化加强流行病学观察性研究报告规范(STROBE-MR) 清单,标明稿件中每个推荐报告项目所对应的位置。 请点击此处下载此文件。
补充图1. 正向孟德尔随机化分析的散点图。
该散点图显示了工具性单核苷酸多态性(SNPs)对系统性红斑狼疮(SLE)的遗传效应与自然流产次数之间的关联。每个数据点代表一个SNP,水平和垂直误差线表示SNP效应估计值的标准误。回归线分别对应逆方差加权法、MR-Egger法、加权中位数法和加权众数法这几种孟德尔随机化方法。 请点击此处下载该文件。
补充图2. 前向孟德尔随机化分析中单核苷酸多态性特异性的因果估计。
森林图显示了每个工具性单核苷酸多态性(SNP)在系统性红斑狼疮(SLE)与自然流产次数之间关联中的因果效应估计值。黑色点表示具有95%置信区间的SNP特异性效应估计值。红色点表示使用逆方差加权法和MR-Egger方法获得的总体因果效应估计值。垂直虚线表示无效应线。 请点击此处下载该文件。
补充图3. 正向孟德尔随机化分析的留一法敏感性分析。
森林图展示了留出一个工具单核苷酸多态性(SNP)后,系统性红斑狼疮(SLE)与自然流产次数之间关联的敏感性分析结果。每个黑点表示依次排除一个工具SNP后得到的总体逆方差加权因果效应估计值,水平线表示相应的95%置信区间。红点表示使用所有工具SNP得到的总体逆方差加权估计值。垂直虚线表示无效应线。 请点击此处下载该文件。
补充图4. 正向孟德尔随机化分析的漏斗图。
该漏斗图显示了系统性红斑狼疮(SLE)与自然流产次数之间关联的单核苷酸多态性(SNP)特异性因果效应估计值的分布情况。每个点代表一个工具性单核苷酸多态性(SNP)。垂直线表示使用逆方差加权法和MR-Egger法获得的总体因果效应估计值。y轴表示标准误的倒数(1/SE)。 请点击此处下载该文件。
补充表 1. 前向孟德尔随机化分析中选定的工具性单核苷酸多态性位点。
该表格列出了用于系统性红斑狼疮与自然流产次数前向孟德尔随机化分析的工具性单核苷酸多态性位点(SNPs),包括最近的注释基因、染色体、基因组位置、效应等位基因、效应等位基因频率、效应量(Beta)、标准误(SE)、P 值以及 F 统计量。染色体位置基于全基因组关联研究中所使用的参考基因组组装。F 统计量的计算公式为 Beta2/SE2。 请点击此处下载该文件。
补充表 2. 用于反向孟德尔随机化分析的工具性单核苷酸多态性位点。
该表格列出了以自然流产次数为暴露因素、系统性红斑狼疮为结局进行反向孟德尔随机化分析所使用的工具性单核苷酸多态性位点(SNPs),包括最近注释的基因、染色体、基因组位置、效应等位基因、效应等位基因频率、效应大小(Beta)、标准误(SE)、P 值以及 F 统计量。染色体位置基于全基因组关联研究中所使用的参考基因组组装。F 统计量的计算公式为 Beta2/SE2。 请点击此处下载该文件。
通过双向孟德尔随机化(MR)分析结合多维度生物信息学分析,本研究发现系统性红斑狼疮(SLE)与复发性流产(RPL)之间存在正向的因果关联,并系统筛选了二者共有的转录组生物标志物。据我们所知,这是首个整合双向MR、转录组分析和免疫浸润分析来探讨这一关联的研究。尽管观察到的MR效应量较小(IVW OR = 1.01),但该关联在多种互补的MR方法和敏感性分析中均得到一致支持,且无显著异质性、水平多效性或强影响异常值的证据,提示所观察到的关系在统计学上稳健但效应量较小。因此,当前研究结果应被解读为支持SLE对RPL易感性具有轻微的遗传贡献,而非显著的临床效应。通过整合因果推断、转录组验证和免疫浸润分析,本研究为优先筛选复杂免疫介导的生殖系统疾病的候选生物标志物提供了可重复的研究框架。一项在埃及于2007年至2021年间开展的研究纳入了123名SLE女性患者,共201次妊娠,报道其中20.4%的妊娠以胎儿丢失告终44。先前的研究同样表明,SLE是RPL的重要风险因素,因为免疫失调可能增加妊娠丢失的风险12。
生物信息学分析显示,59个共同差异表达基因(DEGs)主要富集于抗病毒免疫应答、细胞黏附以及细胞凋亡调控相关的通路。病毒感染可能参与系统性红斑狼疮(SLE)的发病机制。SLE患者常表现出先天性和适应性免疫应答功能障碍45,46,使其更易发生病毒感染。这种易感性的增加可能通过胎盘炎症和滋养层细胞损伤等机制导致妊娠丢失47。作为胎盘的重要组成部分,滋养层细胞的自噬及生物学行为改变也与复发性流产(RPL)的发生相关48,49。综合来看,这些观察结果提示SLE相关的免疫失调可能通过影响滋养层细胞功能而影响妊娠结局。进一步分析鉴定出IFI27和CXCL11为SLE与RPL共有的候选枢纽基因。然而,由于IFI27在多个独立数据集中表现出更高的生物学一致性,因此被优先用于后续分析。尽管这两个基因均被LASSO模型筛选出,但仅IFI27在发现集和外部验证集中均表现出一致的差异表达,而CXCL11在验证数据集中未能稳定重复。此外,IFI27对RPL具有更强的诊断区分能力,并且在血液和生殖组织数据集中均持续显著失调。综上所述,这些结果支持IFI27相较于CXCL11是一个更为稳健的候选生物标志物,但仍需进一步实验验证。利用SLE数据集GSE50772和RPL数据集GSE198700进行验证表明,IFI27的表达在各验证数据集中均持续失调。值得注意的是,IFI27在SLE患者的血液样本中呈高表达,与既往研究结果一致50,但在RPL患者的子宫内膜和绒毛组织样本中则呈低表达。这种相反的表达模式可能反映了SLE合并RPL时,全身性免疫失调与母胎界面局部免疫微环境之间的差异。
IFI27 是一种干扰素刺激基因,参与抗病毒免疫、干扰素信号传导以及病毒感染后的宿主免疫应答51,52。在正常妊娠过程中,IFI27 在滋养层细胞中的表达显著上调53,提示其在维持滋养层功能中具有重要的生理作用。相比之下,我们的分析显示,复发性流产(RPL)患者子宫内膜和绒毛组织中 IFI27 的表达降低。尽管该发现与部分先前报道不同54,但应谨慎解读,因为本研究整合了来自不同组织的转录组数据集,而非配对的母-胎样本。一种可能的解释是,系统性红斑狼疮(SLE)中慢性的I型干扰素系统激活可导致循环免疫细胞中持续的干扰素信号传导,同时在母-胎界面诱导受体脱敏、免疫耗竭或代偿性负反馈机制。另一种可能是,外周血与生殖组织之间在细胞组成或组织特异性的表观遗传调控方面的差异,可能抑制局部 IFI27 的表达,即使存在全身性干扰素激活。这些假说仍属推测,需通过配对的母体血液、子宫内膜组织和滋养层样本进行机制验证,理想情况下应在单细胞水平上开展研究,以区分组织特异性和细胞类型特异性的调控机制55。
免疫浸润分析表明,系统性红斑狼疮(SLE)和复发性流产(RPL)中免疫细胞特征存在显著差异,主要表现为CD4+ T细胞相关群体的改变。IFI27的表达在两种疾病中均与Th2细胞富集呈正相关;然而,这些发现代表的是通过ssGSEA计算得出的相关性,而非经实验验证的生物学相互作用。既往研究表明,SLE患者外周血中Th1细胞和Treg细胞比例降低,而Th2细胞比例升高56,57,这与我们的研究结果一致。在正常妊娠过程中,母胎界面的Th1/Th2免疫平衡会向Th2优势状态偏移58。因此,生殖组织中IFI27表达的下调可能反映了与母胎耐受受损相关的局部免疫稳态改变,但IFI27是否直接调控该过程仍有待实验验证。
应当承认存在若干局限性。首先,尽管孟德尔随机化分析支持因果关联,但所估计的遗传效应相对较小,提示系统性红斑狼疮(SLE)仅代表复发性流产(RPL)多因素发病机制中的一个组成部分。其次,转录组整合分析纳入了来自不同组织(外周血、子宫内膜和绒毛组织)、不同微阵列平台以及独立队列的数据集,尽管IFI27的验证结果一致,但仍可能引入生物学和技术创新异质性。第三,用于绒毛组织的外部验证队列样本量有限,可能降低了统计效能和结果的普适性。第四,由于公开可用的数据集所含临床信息有限,一些重要因素如疾病活动度、抗磷脂抗体状态、药物暴露、妊娠阶段及其他临床协变量无法得到充分评估。最后,尽管PhenoScanner筛查已最大限度减少了孟德尔随机化分析中潜在的多效性混杂,但仍不能完全排除残余混杂因素的影响。
从转化医学的角度来看,IFI27 目前应被视为一种候选生物标志物,而非经过临床验证的诊断标志物。在进入临床应用之前,尚需开展前瞻性多中心研究,以验证其在不同人群中的诊断效能,建立标准化的检测平台和诊断阈值,并明确妊娠阶段、疾病活动度以及免疫抑制治疗如何影响 IFI27 的表达水平。此外,功能实验结合空间转录组学和单细胞转录组学分析,对于阐明 IFI27 是主动参与母胎免疫调节,还是仅仅反映了干扰素驱动的免疫激活,也将至关重要。
利益冲突:
作者声明不存在任何竞争性财务或非财务利益。
本研究由北京市中医管理局中西医结合重大疑难疾病重点攻关项目(2023BJSZDYNJBXTGG-003)、国家级公益科研院所基本科研业务费专项资金(ZZ16-XRZ-038)以及高水平中医医院提升项目(HLCMHPP2023087)资助。资助方在研究设计、数据收集、数据分析、数据解释、论文撰写以及决定投稿发表等过程中均未参与。作者感谢FinnGen研究、英国生物样本库(UK Biobank)以及美国国家生物技术信息中心基因表达综合数据库(GEO)的研究人员和参与者公开其数据集。作者同时感谢FinnGen联盟,该联盟通过芬兰各研究机构、生物样本库及国际合作伙伴之间的协作,将芬兰生物样本库样本与全国健康登记数据进行了整合。
| 姓名 | 公司 | 目录编号 | 评论 |
|---|---|---|---|
| 28个免疫细胞基因特征集合 | 已发表的补充基因特征资源 | Charoentong 等人(参考文献 42)描述的补充性免疫细胞标志基因列表 | 不适用 RRID: 暂无内容 目的/注意事项: 用于ssGSEA的免疫细胞特征谱 |
| 中枢神经系统知识库 | CNSknowall 网络平台 | DAVID 输出文件 | 不适用 RRID: 暂无内容 目的 / 说明: 功能富集分析结果的可视化 |
| 计算机工作站 | 机构计算环境 | 不适用 | 不适用 RRID: 不适用 目的 / 说明: 计算分析。 |
| 耗材 | 不适用 | 不适用 | 不适用 RRID: 不适用 目的 / 说明: 未使用湿实验耗材。 |
| Cytoscape | Cytoscape 联盟 | 不适用 | 3.10.0 RRID: SCR_003032 目的 / 说明: 蛋白质–蛋白质相互作用网络可视化与拓扑学分析 |
| cytoHubba | Cytoscape 应用商店 | 不适用 | 0.1 RRID: SCR_017677 目的 / 说明: 使用MCC、MNC、EPC、Degree、Closeness和Radiality进行枢纽基因排序。 |
| DAVID 功能注释工具 | 美国国立卫生研究院 / 国立癌症研究所 | 上传共享的差异表达基因列表及背景基因列表 | 2021 RRID: SCR_001881 目的 / 备注: 基因本体与KEGG通路富集分析 |
| 欧洲LD参考面板 | 1000基因组计划 / IEU OpenGWAS | 3 期欧洲人群组(与 GRCh37 兼容的变异) | 3 期 RRID: 未报告 目的 / 说明: 通过 OpenGWAS/TwoSampleMR 工作流程进行连锁不平衡剪枝 |
| FinnGen | FinnGen 联盟 | finngen_R11_L12_LUPUS | 发布 11 RRID: SCR_022254 目的 / 备注: 系统性红斑狼疮的全基因组关联研究汇总统计学数据 |
| 森林图绘制工具 | CRAN | 不适用 | 1.1.2 RRID: 暂无内容 目的 / 说明: 孟德尔随机化估计的森林图可视化 |
| 基因表达综合数据库(GEO) | 国家生物技术信息中心 | GSE61635;GSE165004;GSE50772;GSE198700 | 不适用 RRID: SCR_005012 目的 / 说明: 发现与验证转录组数据集的来源。 |
| 基因本体 | 基因本体联盟 | 通过 DAVID 获取的 GO 术语 | DAVID 2021 注释 RRID: SCR_002811 目的 / 备注: 生物过程、细胞组分和分子功能注释。 |
| ggcorrplot | CRAN | 不适用 | 0.1.4.1 RRID: 暂无内容 目的 / 说明: 候选基因的可视化–免疫细胞相关性矩阵 |
| ggplot2 | CRAN | 不适用 | 3.5.1 RRID: SCR_014601 目的 / 说明: 火山图、箱线图及其他统计图形。 |
| ggvenn | CRAN | 不适用 | 0.1.16 RRID: SCR_025300 目的 / 说明: 共享差异表达基因的可视化 |
| glmnet | CRAN | 不适用 | 4.1-8 RRID: SCR_015505 目的 / 说明: LASSO 逻辑回归与交叉验证 |
| GSE165004 | NCBI GEO | GSE165004 / GPL16699 | 处理后的系列矩阵 RRID: SCR_005012 目的 / 说明: RPL 子宫内膜发现数据集 |
| GSE198700 | NCBI GEO | GSE198700 / GPL13534 | 处理后的系列矩阵 RRID: SCR_005012 目的 / 备注: 独立的RPL绒毛膜绒毛验证数据集 |
| GSE50772 | NCBI GEO | GSE50772 / GPL570 | 处理后的系列矩阵 RRID: SCR_005012 目的 / 说明: 独立的系统性红斑狼疮外周血单个核细胞验证数据集 |
| GSE61635 | NCBI GEO | GSE61635 / GPL570 | 处理后的系列矩阵 RRID: SCR_005012 目的 / 说明: SLE全血发现数据集 |
| GSEABase | Bioconductor | 不适用 | 1.66.0 RRID: 暂无内容 目的 / 说明: ssGSEA中免疫细胞基因集的管理 |
| GSVA | Bioconductor | 不适用 | 1.52.3 RRID: SCR_021058 目的 / 备注: 单样本基因集富集分析(ssGSEA) |
| IEU OpenGWAS | MRC 整合流行病学单位 | ukb-b-419;finngen_R11_L12_LUPUS | 不适用 RRID: 未报告 目的 / 说明: 全基因组关联研究摘要统计结果与标准化遗传关联数据的获取 |
| 京都基因与基因组百科全书(KEGG) | Kanehisa 实验室 | 通过 DAVID 访问的 KEGG 通路 | DAVID 2021 注释 RRID: SCR_012773 目的 / 说明: 通路富集注释 |
| limma | Bioconductor | 不适用 | 3.60.6 RRID: SCR_010943 目的 / 说明: 差异表达分析 |
| MRPRESSO | Verbanck 等人 | 不适用 | 1 RRID: SCR_023697 目的/注意事项: 水平多效性及异常工具变量的检测 |
| pheatmap | CRAN | 不适用 | 1.0.12 RRID: SCR_016418 目的 / 说明: 表达热图 |
| PhenoScanner V2 | PhenoScanner 联盟 | SNP水平表型查询 | 版本 2 RRID: 暂无内容 目的 / 说明: 筛选保留的SNP以排除潜在的混杂表型关联。 |
| pROC | CRAN | 不适用 | 1.18.5 RRID: SCR_024286 目的 / 说明: ROC曲线、AUC、DeLong置信区间、Youden指数截断值和自助法置信区间 |
| R | R 统计计算基金会 | 不适用 | 4.4.2 RRID: SCR_001905 目的 / 说明: 统计计算环境 |
| 试剂 | 不适用 | 不适用 | 不适用 RRID: 不适用 目的 / 说明: 未使用湿实验试剂。 |
| STRING | STRING 联盟 | 智人(分类号 9606);最低相互作用评分 0.400 | 11 RRID: SCR_005223 目的 / 说明: 蛋白质–蛋白质相互作用网络构建 |
| TwoSampleMR | MRC综合流行病学研究组 | 不适用 | 0.6.6 RRID: SCR_019010 目的 / 备注: 双向双样本孟德尔随机化、数据提取、数据整合、因果效应估计及敏感性分析。 |
| 英国生物样本库 | 英国生物样本库 | ukb-b-419 | 2018年汇总数据集 RRID: SCR_012815 目的 / 备注: 自发性流产次数的全基因组关联研究汇总统计结果 |