端粒长度与甲状腺毒症的关联：基于两项样本孟德尔随机化研究的洞察

摘要

本方案描述了一种双样本孟德尔随机化分析流程，旨在利用公开的遗传学汇总数据评估端粒长度是否对甲状腺毒症风险具有因果影响。内容涵盖工具变量选择、数据一致性处理、主要效应估计、敏感性分析，以及可重复的R语言代码，并提供可直接用于绘图的结果输出，以支持透明化报告。

摘要

甲状腺毒症是一种以甲状腺激素过多为特征的内分泌疾病，但其与全身性衰老生物学之间的病因学联系尚未完全明确。端粒长度（TL）反映了细胞衰老和基因组稳定性，并已被认为与多种复杂疾病相关。我们开展了一项两样本孟德尔随机化（MR）研究，以评估遗传预测的端粒长度对甲状腺毒症风险的因果效应。用于TL的遗传工具变量来自一项大规模的欧洲人群全基因组关联研究（GWAS）（n > 470,000）。甲状腺毒症的汇总统计数据取自最新的FinnGen发布版本（≈4,000例病例和>210,000例对照）。主要的逆方差加权分析表明，遗传预测较长的TL与较低的甲状腺毒症风险相关，该效应的方向和大小在多种互补估计方法（MR-Egger、加权中位数/最大似然法、MR-PRESSO和MR-RAPS）中保持一致。敏感性分析未发现定向多效性的证据，异质性通过Cochran's Q检验进行评估。Steiger方向性检验支持从TL到甲状腺毒症的因果路径。

据我们所知，本研究是首批利用当前全基因组关联研究资源评估总体甲状腺毒症风险与端粒长度（TL）之间因果关系的孟德尔随机化分析之一，拓展了既往聚焦于甲状腺功能亢进相关表型的研究证据。这些发现提示，由TL反映的细胞衰老过程可能影响甲状腺毒症的易感性，并为进一步开展关于端粒生物学在甲状腺功能障碍中作用的纵向及机制性研究提供了动力。

引言

甲状腺毒症是一种常见的内分泌疾病，其特征是循环中甲状腺激素水平过高。典型表现包括体重减轻、心动过速、易激惹、焦虑和震颤，严重病例易发生心律失常、心力衰竭和骨量丢失¹。尽管临床特征和治疗策略已相对明确，但影响个体间易感性和疾病进展的生物学机制仍不完全清楚。病因包括自身免疫性功能亢进（如格雷夫斯病）、毒性结节、甲状腺炎、外源性激素暴露、药物反应以及肿瘤相关性原因²。在这些不同病因中，免疫激活、氧化应激以及线粒体生物合成异常反复出现，成为连接激素过量与组织损伤及全身并发症的可能机制主题。

端粒长度（TL）是染色体末端的串联重复DNA-蛋白质结构，反映了细胞衰老、复制历史和基因组稳定性。随着细胞分裂，TL会逐渐缩短，而氧化应激和慢性炎症可加速这一过程，这些过程在甲状腺功能障碍中同样显著^3,4。较短的TL与心脑血管代谢性疾病、某些癌症以及神经退行性疾病相关⁵，新兴证据还提示其与肥胖和2型糖尿病等内分泌表型存在关联，这些疾病中激素和免疫-代谢通路与复制应激相互交汇⁶。鉴于甲状腺激素调控基础代谢率、线粒体功能和增殖信号通路，TL动态变化在甲状腺毒症状态下可能尤为相关，因为全身性高代谢和免疫失调可能与端粒生物学相互影响^7,8。

然而，端粒长度（TL）在甲状腺毒症中是否起因果作用仍不明确。传统的观察性研究容易受到混杂因素（如生活方式、合并症）和反向因果关系（疾病影响TL）的干扰。为解决这些局限性，我们采用两样本孟德尔随机化（MR）方法，利用与TL稳健相关的生殖系遗传变异作为无混杂的工具变量，检验其对疾病风险的影响⁹。全基因组关联研究现已提供了足够强效的TL工具变量以及大规模、表型详尽的结果数据集，从而能够进行严格的因果评估¹⁰。本研究通过整合遗传流行病学与端粒生物学，评估基因预测的TL是否影响甲状腺毒症的发病风险，进而阐明其潜在的机制通路，并探索TL作为风险分层及未来干预开发的候选生物标志物的可能性。

据我们所知，本分析是首批利用现代大规模全基因组关联研究（GWAS）资源来探讨端粒长度（TL）与甲状腺毒症总体风险关系的孟德尔随机化研究之一，拓展了以往主要关注甲状腺功能亢进相关表型或自身免疫性甲状腺疾病的研究¹¹。

方案

本研究采用双样本孟德尔随机化（MR）设计，评估基因预测的端粒长度（TL）是否对甲状腺毒症风险具有因果影响。分析仅使用去标识化的公开可用的全基因组关联研究（GWAS）汇总统计数据，未访问任何个体水平数据。青岛市立医院机构审查委员会认定，由于本分析完全依赖于公开的汇总数据，因此无需进一步审查。所有提供数据的全基因组关联研究在其原始研究方案中均已获得知情同意和伦理批准。本分析通过使用与TL相关的生殖系变异作为无混杂的工具变量，并实施敏感性分析程序以检验多效性、异质性及效应方向，从而最大限度地减少混杂因素和反向因果关系的影响。

数据来源
端粒长度（TL）的暴露数据集来自 IEU OpenGWAS，标识符为 ieu-b-4879，包含约 472,174 名欧洲血统参与者。甲状腺毒症结局数据集来自 FinnGen 联盟，采用 2021 年终点 finn-b-thyROTOXICOSIS，包含 4,142 例病例和 213,693 例对照。这些资源提供了构建单核苷酸多态性（SNP）水平工具、等位基因整合以及使用既定的孟德尔随机化（MR）估计方法评估因果效应所需的效应估计值和标准误。

软件与计算环境
所有分析均在 R（版本 4.3.1）中进行。主要分析软件包为 TwoSampleMR（版本 0.5.7），并结合使用 ieugwasr 以编程方式访问 GWAS 资源，MRPRESSO 用于异常值检测与偏差检验，mr.raps 在弱工具变量和特异性多效性条件下进行稳健估计，RadialMR 用于径向可视化，以及 tidyverse 和 data.table 等通用软件包。会话信息与软件包版本已保存至文件，以确保严格的可重复性。

孟德尔随机化的核心假设
孟德尔随机化（MR）框架基于以下前提：遗传工具变量与端粒长度（TL）密切相关，且独立于TL与甲状腺毒症之间关系的混杂因素，并仅通过TL影响甲状腺毒症，而非经由其他替代通路。分析计划通过量化每个变异体的效应强度、使用基于截距的方法和异常值筛查检测不平衡的水平多效性、评估变异体间的异质性，以及通过正式的方向性检验确认从TL到甲状腺毒症的因果方向，来实现对这些假设的验证¹²。

工具变量选择与质量控制
从TL全基因组关联研究中筛选遗传工具变量，筛选标准为全基因组显著性水平p值阈值 < 5×10^-8. 为确保独立性，通过连锁不平衡剔除（clumping）处理，设定 r² = 0.001 的范围内，使用1000基因组计划的欧洲参考数据，在10,000 kb 窗口内进行筛选；若存在相关性变异，则保留与端粒长度（TL）关联 p 值更小的变异。每个变异的工具强度通过统计量 F_i = \beta_{E,i}^{2}/\mathrm{SE}(\beta_{E,i})^{2} 进行汇总，并依据传统的 F 值标准进行评估。 > 10项标准^13,14. Where appropriate, the Sanderson-Windmeijer approximation was referenced to describe aggregate instrument strength across multiple variants in a sample of size N, with K instruments explaining exposure variance R^{2}: F_{\mathrm{SW}} = \{[R^{2}/(1-R^{2})]\,(N-K-1)\}/K ¹⁵等位基因频率接近0.5且为模糊回文变异的位点被排除，除非等位基因频率信息允许明确比对。

数据整合
针对每个工具变量，从 FinnGen 中提取其与甲状腺毒症的相应关联，并进行对齐，使效应大小在暴露和结局中代表相同的效应等位基因。数据整合过程去除了等位基因不匹配的情况，纠正了链方向问题，并排除了具有未解决歧义的回文型单核苷酸多态性（SNP），从而生成适用于有效构建 Wald 比值的数据库。

主要因果估计
主要分析采用逆方差加权（IVW）估计法，将各单核苷酸多态性（SNP）特异性的Wald比值合并为一个综合因果效应¹⁶。结果以端粒长度（TL）每增加一个标准差的比值比表示，95%置信区间基于固定效应和随机效应模型计算得出。模型选择由异质性诊断结果决定，同时呈现两种模型结果以支持稳健解读。

敏感性分析与稳健性检验
通过 MR-Egger 回归（包含截距项）评估稳健性，以检验是否存在不平衡的方向性多效性；采用最大似然估计法，在满足同质性的前提下提高估计效率并减轻测量误差；使用 MR-PRESSO 进行全局离群值检验，并估算偏差及经离群值校正后的效应；应用 MR-RAPS 提供对弱工具变量和特异性多效性具有稳健性的估计结果¹⁷。在 IVW 框架下，采用 Cochran's Q 统计量量化不同遗传变异之间的异质性¹⁸。利用 Steiger 检验分析效应方向，该检验通过比较暴露和结局中被解释的方差比例，判断数据更支持端粒长度（TL）导致甲状腺毒症，还是相反的因果关系。MR-Egger 截距被用作方向性多效性的正式检验方法。通过逐一剔除分析（leave-one-out analyses）检查整体关联是否由某个单一变异主导。

反向孟德尔随机化
为了探究潜在的反向因果关系，将分析流程重复进行，以甲状腺毒症作为暴露因素，端粒长度（TL）作为结局。采用相同的工具变量筛选标准、数据协调步骤、逆方差加权（IVW）主要估计方法及敏感性分析，以确保关于因果方向的结论均在相同的因果框架内得出。

多重检验、检验效能与结果报告
主要假设检验针对的是 TL 对 thyrotoxicosis 影响的 IVW 估计值。敏感性分析估计值及诊断性检验结果被视为支持性证据；p 值以科学计数法报告以确保清晰，结论强调各方法间结果的一致性，而非孤立的显著性阈值。工具变量强度的汇总统计及解释的变异比例，基于标准非中心性分布公式提供近似的检验效能，同时需认识到检验效能取决于样本量、工具变量强度以及真实效应大小。

计算可重复性与确切命令
通过提供完整的 R 命令序列，确保研究的可重复性，这些命令可重现工具变量筛选、结局提取、数据一致性匹配、主要分析与敏感性分析、反向 MR 分析、诊断结果输出，以及可直接用于分析的文件导出。该脚本生成稳定的、便于阅读的 CSV 文件，分别对应工具变量列表、一致性匹配后的数据集，以及 MR 估计结果和诊断信息的汇总。

# R 4.3.1；TwoSampleMR 0.5.7
# 可选安装：
# install.packages(c("TwoSampleMR","ieugwasr","MRPRESSO","mr.raps",
"RadialMR","tidyverse","data.table"))

library(TwoSampleMR)
library(ieugwasr)
library(MRPRESSO)
library(mr.raps)
library(RadialMR)
library(tidyverse)
library(data.table)

# 暴露因素：端粒长度（IEU OpenGWAS）
exposure_id <- "ieu-b-4879"

# 结果：所报告分析中使用的 FinnGen 2021 甲状腺毒症终点
outcome_id <- "finn-b-thyROTOXICOSIS"

# 使用全基因组阈值和严格LD剪枝进行工具变量筛选
exp <- extract_instruments(outcomes = exposure_id, p1 = 5e-8, clump = TRUE, r2 = 0.001, kb = 10000)
exp$F_stat <- (exp$beta.exposure^2) / (exp$se.exposure^2)
fwrite(exp, "S1_instruments_TL.csv")

# 结果提取与数据整合
out <- extract_outcome_data(snps = exp$SNP, outcomes = outcome_id)
dat <- harmonise_data(exp, out, action = 2)
fwrite(dat, "S2_harmonised_TL_vs_thyrotoxicosis.csv")

# 主要 MR 及敏感性分析估计方法
res <- mr(dat, method_list = c("mr_ivw","mr_ivw_fe","mr_egger_regression",
"mr_weighted_median","mr_raps","mr_maxlik"))

het <- mr_heterogeneity(dat) # Cochran's Q
pleio <- mr_pleiotropy_test(dat) # Egger 截距
steiger <- directionality_test(dat) # Steiger 方向性检验

# MR-PRESSO 全局估计值及异常值校正估计值
mrpresso <- mr_presso(BetaOutcome = "beta.outcome", BetaExposure = "beta.exposure",
SdOutcome = "se.outcome", SdExposure = "se.exposure",
OUTLIERtest = TRUE, DISTORTIONtest = TRUE,
data = dat, NbDistribution = 1000, SignifThreshold = 0.05) #

# 逆向 MR：甲状腺毒症（暴露因素） -> 端粒长度（结果）
rev_exp <- extract_instruments（outcomes = outcome_id, p1 = 5e-8, clump = TRUE, r2 = 0.001, kb = 10000）
反向输出 <- extract_outcome_data(snps = rev_exp$SNP, outcomes = exposure_id)
rev_dat <- harmonise_data(rev_exp, rev_out, action = 2)

rev_res <- mr(rev_dat, method_list = c("mr_ivw","mr_ivw_fe","mr_egger_regression",

"mr_weighted_median","mr_raps","mr_maxlik"))
rev_het <- mr_heterogeneity(rev_dat)
rev_pleio <- mr_pleiotropy_test(rev_dat)
rev_steiger <- directionality_test(rev_dat)

# 用于归档及生成图表的导出命令
write.csv(bind_rows(res), "S3_mr_results_primary.csv", row.names = FALSE)
write.csv(het, "S3_mr_heterogeneity.csv", row.names = FALSE)
write.csv(pleio, "S3_mr_pleiotropy_egger.csv",row.names = FALSE)
write.csv(steiger, "S3_mr_steiger.csv", row.names = FALSE)
write.csv(bind_rows(rev_res), "S3_reverse_mr_results.csv", row.names = FALSE)
write.csv(rev_het, "S3_reverse_heterogeneity.csv", row.names = FALSE)
write.csv(rev_pleio, "S3_reverse_pleiotropy_egger.csv", row.names = FALSE)
write.csv(rev_steiger, "S3_reverse_steiger.csv", row.names = FALSE)

结果

仪器特征与假设检验
从端粒长度全基因组关联研究（GWAS）中，保留了129个全基因组显著且连锁不平衡（LD）独立的工具变量（p < 5×10^-8；r² < 0.001；10,000 kb 窗口）。绝大多数单核苷酸多态性（SNP）的F统计量均超过常规阈值，表明弱工具偏倚的风险有限。通过数据协调去除了意义不明确的回文变异，并对等位基因进行对齐，使暴露和结局效应均指向相同的效应等位基因，从而得到适用于Wald比率构建和合并估计的数据集。

主要因果估计
所有估计值在方向上一致，表明遗传预测的较长端粒长度与较低的甲状腺毒症风险相关。主要的IVW分析结果为OR = 0.689，95% CI 0.564–0.842，p < 1.0×10^-3。最大似然估计法得出OR = 0.689，95% CI 0.576–0.823，p < 1.0×10^-3。MR-PRESSO（原始）分析结果为OR = 0.689，9...

讨论

这项孟德尔随机化（Mendelian randomization, MR）研究提供了证据，表明遗传预测的较长端粒可降低甲状腺毒症的风险¹⁹。在逆方差加权法（IVW）、MR-Egger、最大似然法、MR-PRESSO 和 MR-RAPS 多种方法中，方向一致的估计值表明存在稳健的保护性关联，而方向性检验支持端粒长度（TL）→ 甲状腺毒症的因果路径，而非反向因果²⁰。其生物学合理性与端粒在维持染色体完整性和调控细胞寿命中的作用相符。进行性端粒缩短会加速细胞衰老和功能障碍²¹，而甲状腺毒症的特征则是代谢转换增强和增殖信号上调²²。在此高代谢需求状态下，较长的端粒可能缓冲氧化应激并减少DNA损伤，从而维持上皮细胞和滤泡细胞的完整性，并延缓病理性重塑²³。类似地，较长端粒与其他年龄相关疾病（包括心血管疾病和某些癌症）风险降低之间的保护性关联，进一步支持了端粒动态在疾病易感性中的系统性作用

材料

本文使用的材料清单

姓名	公司	目录编号	评论
1000 Genomes Project	连锁不平衡（LD）clumping 的参考数据。	1000 Genomes Project	-
Cochran’s Q statistic	用于孟德尔随机化（MR）中异质性检验的统计方法。	N/A	-
FinnGen 数据集	甲状腺毒症（病例与对照）的全基因组关联研究（GWAS）汇总统计量。	FinnGen Consortium	2021 年版本，finn-b-thyROTOXICOSIS
IEU OpenGWAS 数据集	端粒长度的全基因组关联研究（GWAS）汇总统计量。	IEU OpenGWAS	数据集编号：ieu-b-4879
ieugwasr 包	用于程序化访问 GWAS 数据的 R 软件包。	CRAN	版本 0.1.5
mr.raps 包	在工具变量较弱的情况下，用于孟德尔随机化中稳健估计的 R 软件包。	CRAN	版本 0.1.2
MR-Egger 回归	用于检测方向性多效性的敏感性分析方法。	N/A	-
MR-PRESSO 异常值检验	用于识别并校正孟德尔随机化分析中异常值的统计检验方法。	N/A	-
MRPRESSO 包	用于孟德尔随机化分析中异常值检测与偏差检验的 R 软件包。	CRAN	版本 1.0
R 软件	用于所有分析的统计计算与图形环境。	统计计算基金会 R Foundation for Statistical Computing	版本 4.3.1
RadialMR 包	用于孟德尔随机化结果径向可视化的 R 软件包。	CRAN	版本 1.1.2
TwoSampleMR 包	用于进行孟德尔随机化分析的 R 软件包。	CRAN	版本 0.5.7