我见过太多人做完GWAS拿到一张漂亮的曼哈顿图之后就卡在下一步了。最经典的操作是把显著位点拉到UCSC或Ensembl里找离它最近的基因然后宣布“找到了疾病相关基因”。但这个做法在遗传学圈子里经常翻车因为连锁不平衡LD会让信号“挂”在错误的基因旁边而eQTL就是用来解决这个问题的关键工具。把GWAS和eQTL结合起来做共定位或孟德尔随机化分析已经成为目前筛选疾病致病基因的主流策略也是我从统计信号走向生物学机制最实用的一条路径。这篇内容适合手里已经有一批GWAS显著位点、但不知道如何往下推到基因层面的研究者也适合刚接触遗传学分析、想搞懂这两个术语到底怎么用的学生。1. GWAS结果交到手上之后真正的麻烦才开始1.1 显著位点不等于致病位点更不等于致病基因GWAS输出的核心结果是一堆达到全基因组显著性水平的SNP比如p值小于5×10⁻⁸。很多初学者会本能地把这些SNP当作“致病变异”但实际情况要复杂得多。人类基因组在传递过程中不是整段完整遗传的而是以“块”为单位——这就是连锁不平衡。一个位点上的变异往往和周围一大段区域里的其他变异高度关联。这意味着GWAS捡到的“最显著SNP”未必是真正起功能作用的变异它可能只是因为和那个真正的因果变异靠得近、相关性高而被检测到。更有意思的是GWAS显著位点里有相当大比例位于基因间区或内含子区域它们本身不改变氨基酸序列却可能通过影响调控元件来改变基因表达水平。这就引出一个关键问题从统计关联到生物学机制中间隔着一条“功能注释”的鸿沟。如果你只停留在SNP层面你能回答的是“某个基因组区域与疾病相关”但无法回答“这个区域里到底是哪个基因在捣乱”。这时候就需要引入eQTL数据。eQTL做的就是给每个SNP和基因表达量之间搭一座桥某个SNP的基因型不同某个基因的表达水平是否随之显著改变。如果GWAS显著位点恰好也是一个eQTL位点那就意味着这个位点很可能通过调控特定基因的表达来影响疾病风险。1.2 “最近基因法”为什么会带偏方向我刚开始做这一行的时候也走过“离哪个基因最近就选哪个”的弯路。这种策略在某些情况下碰巧是对的但本质上没有理论依据。原因很简单三维基因组结构决定了基因表达的调控不一定遵守“距离最近”原则。一个位于基因A下游50kb的变异完全可能通过染色质环化作用与远处的基因B的启动子区域接触从而调控基因B的表达。在二维基因组坐标上看起来“近”的基因在三维空间里可能压根儿不相关。业内已经有不少案例展示了这种误导性。比如某些GWAS位点落在基因间区按最近距离归属于基因C但后续功能实验证明真正受调控的是位于数百kb之外的基因D。如果你只靠距离来判断后续所有实验敲除、过表达、动物模型都会押错筹码。所以GWAS之后加一步eQTL共定位分析本质上是在用功能证据校正统计信号的归属问题。1.3 关联、共定位与因果推断三者不是一回事在进入实操之前有几组概念必须分清楚不然你在阅读文献和跑工具的时候会很痛苦。关联GWAS发现的“这个SNP和这个疾病在统计上相关”不涉及方向性和机制。共定位在某个基因组区域内GWAS信号和eQTL信号是否由同一个因果变异驱动。如果共定位后验概率高说明这个位点的基因型和疾病之间的关系很可能通过基因表达来介导。因果推断更进一步用孟德尔随机化的思想以基因型为工具变量评估基因表达水平对疾病风险的因果效应。从关联到共定位再到因果推断是证据强度逐级提升的过程。我见过有的人拿到coloc结果后直接说“证实了XX基因导致XX疾病”其实coloc只是支持“存在共同因果变异”的证据真正的因果结论还需要结合MR结果、功能实验甚至动物模型才能站得住脚。2. eQTL数据的底层逻辑为什么它能帮GWAS“指认”基因2.1 从DNA到RNA之间的调控层eQTL的全称是expression Quantitative Trait Loci翻译过来是“表达数量性状基因座”。它研究的是一个非常简单的问题不同个体的基因组序列有差异这些差异是否会导致基因表达水平的差异。举个例子假设有1000个人某个SNP位点上有的人带A等位基因有的人带G等位基因。如果带A的人群体中某个基因的平均表达量显著高于带G的人群那么这个SNP就是这个基因的一个eQTL。这里的“表达量”指的是RNA水平通常用RNA-seq数据来量化。也就是说eQTL分析需要两类数据基因型数据SNP分型和基因表达数据RNA-seq而且这两类数据必须来自同一批个体。这一层机制非常重要因为很多疾病相关变异并不会改变蛋白质的氨基酸序列而是通过影响基因的转录调控来发挥作用。比如改变了转录因子结合位点、影响了启动子活性、或者干扰了增强子与启动子的相互作用。eQTL分析正是把这类调控效应转化为可量化的统计关联。2.2 cis-eQTL和trans-eQTL为什么优先看ciseQTL按照作用距离分为两类cis-eQTLSNP调控的是它自身附近通常定义为1Mb范围内的基因表达。这种调控通常比较直接比如SNP位于某个基因的启动子区域直接影响该基因的转录效率。trans-eQTLSNP调控的是距离很远的基因甚至跨染色体。这种调控往往通过间接机制实现比如SNP影响了某个转录因子的表达这个转录因子又去调控下游一大批基因。在实际的GWAS-eQTL联合分析中绝大多数人只看cis-eQTL。原因有几个方面。第一cis-eQTL的效应通常更强、更稳定跨数据集的可重复性更高。第二cis-eQTL信号的LD结构相对清晰容易做共定位分析。第三从生物学解释的角度看一个位点顺式调控它旁边的基因这个故事讲起来更合理后续功能验证也更容易设计。相比之下trans-eQTL存在严重得多重检验问题一个SNP要测试对所有基因的影响而且远距离调控关系容易受混杂因素干扰。GTEx项目的统计结果显示cis-eQTL在几乎所有组织中都能检测到大量信号而trans-eQTL就稀少得多。因此除非你有非常明确的理由要研究trans调控否则建议把精力集中在cis-eQTL上。2.3 主流eQTL数据资源及选择标准目前eQTL分析可用的公开数据主要来自几个大型项目每个项目有自己的人群特征和组织覆盖范围选择时不能闭着眼睛抓一个就用。数据资源覆盖组织样本规模适合场景GTExV849个组织~838个供体多组织差异分析器官特异性疾病首选eQTLGen Consortium全血~31,684人血液相关性状或大样本验证DGN全血~922人有RNA-seq原始数据可自定义分析TCGA多种癌症组织数百至上千例癌症eQTL分析结合肿瘤表达谱Brain eQTL如MetaBrain、ROSMAP脑组织多个区域数千例神经精神类疾病选择标准我个人的经验是病在哪个组织就优先用哪个组织的eQTL数据。比如研究冠心病用GTEx的冠状动脉组织或心脏组织研究阿尔茨海默病一定要找脑组织来源的eQTL数据而不是偷懒用全血的。因为eQTL的组织特异性非常强一个SNP在肝脏里可能是某个基因的eQTL在脑子里可能完全没有调控效应。如果你用了错误的组织数据共定位结果大概率是阴性。如果疾病相关组织没有足够的公开eQTL数据退而求其次可以用全血但必须明确说明这是替代方案并且最好做一个跨组织的敏感性分析在不同组织的数据里看结果是否稳定。3. 把两套信号放在一起共定位与SMR的本质区别和执行细节3.1 coloc和SMR/HEIDI分别回答什么问题GWAS和eQTL都拿到之后联合分析的主流工具有两个流派Bayesian共定位代表工具是coloc和基于工具变量的回归方法代表工具是SMR/HEIDI。coloc的核心思想是在一个区域内同时观察GWAS关联信号和eQTL关联信号评估这两种信号是否共享同一个因果变异。它输出的是五个后验概率PP.H0该区域内没有任何因果变异PP.H1只有GWAS信号有因果变异PP.H2只有eQTL信号有因果变异PP.H3两个信号有各自的因果变异但它们不共享PP.H4两个信号共享同一个因果变异我们通常关注PP.H4如果它大于0.75或0.8说明共定位证据充分。**SMRSummary-data-based Mendelian Randomization**走的是另一条思路。它把eQTL的效应量SNP对表达的影响当作暴露把GWAS的效应量SNP对疾病的影响当作结局用类似孟德尔随机化的框架去估算“基因表达每变化一个单位疾病风险变化多少”。但SMR存在一个问题两个信号可能只是LD造成的关系而不是真正的因果介导。因此SMR内置了一个HEIDI过滤步骤用来检测关联信号是否是由于LD造成的而不是同一个因果变异驱动的。这两者不是竞争关系而是互补关系。我常用的策略是先用SMR做全基因组范围内的快速扫描得到一个候选基因列表然后用coloc对每个候选位点做精细的共定位验证。3.2 coloc操作的完整流程与参数设置关于coloc我需要给你一个可以直接上手的实操路径。假设你已经准备好了GWAS的sumstats文件包含SNP、效应等位基因、效应量、标准误、p值、等位基因频率和eQTL的sumstats文件流程如下。第一步是数据格式统一。coloc要求输入每个SNP的效应量beta、方差通常用标准误平方、MAF、SNP ID。两个数据集的SNP ID必须一致建议都使用chr_pos_allele格式避免rsID在不同版本注释中的混乱问题。第二步是区域划分。coloc假设每个区域只有一个因果变异所以不能把整条染色体丢进去。一般做法是以GWAS最显著SNP为中心取前后500kb或1Mb的区域。如果你用的是PLINK格式数据可以用--ld-snp配合--r2来界定独立信号区间。第三步是运行coloc。核心代码我这里给一个R语言的骨架library(coloc) # 准备数据 dataset1 - list(beta gwas$beta, varbeta gwas$se^2, snp gwas$snp, type cc, s 0.3, MAF gwas$maf) # s参数代表病例占比如果是定量性状就设置typequant dataset2 - list(beta eqtl$beta, varbeta eqtl$se^2, snp eqtl$snp, type quant, MAF eqtl$maf) # 运行共定位 my_coloc - coloc.abf(dataset1, dataset2, p12 1e-5) print(my_coloc$summary)这里有一个容易踩的坑p12参数的设置。p12是先验概率表示该区域同时存在GWAS和eQTL因果变异的概率。默认值是1e-5但如果你筛查的区域很多适当调低一点可以控制假阳性。我一般跑批量分析时用1e-6单个候选位点验证时用默认值。还有一个细节需要注意如果GWAS是case-control设计必须正确设置s参数人群中的患病率它会直接影响PP.H4的计算。很多人忽略这个参数得到的后验概率会偏高或偏低。3.3 SMR操作的参数细节mBIM和HEIDI阈值的正确理解SMR的操作相对简洁。下载SMR软件之后通常用以下命令格式smr --bfile gwas_plink --gwas-summary gwas.ma --beqtl-summary eqtl_merged --out smr_result --thread-num 10--bfileGWAS数据对应的参考面板PLINK文件用于计算LD--gwas-summaryGWAS汇总统计文件--beqtl-summaryeQTL数据的二进制格式SMR输出的每个探针通常对应一个基因会有一个p值。判断显著位点时常用的阈值是p小于0.05除以有效探针数Bonferroni校正。但真正需要关注的是HEIDI检验。HEIDI的p值如果大于0.05说明没有显著证据表明该关联是由LD造成的可以认为SMR信号是真实的因果效应。如果HEIDI p值很小小于0.01那就要警惕了——这个信号很可能只是LD混淆造成的假象。我见过很多人在SMR结果里挑p值很显著的基因却不看HEIDI的过滤结果结果选了一堆假阳性。SMR官方文档的建议是优先保留HEIDI p 0.05且SMR p达到显著阈值的基因这个组合才是比较可信的候选。3.4 从SMR到coloc的筛选流水线实际操作中我不建议把coloc直接跑在几百万个位点上计算量大且多重检验问题很棘手。推荐的做法是先粗后精全基因组范围内运行SMR得到初步候选基因列表。对每个候选基因的所在区域以最显著SNP为中心±500kb提取GWAS和eQTL数据。对每个区域运行coloc记录PP.H4。保留同时满足HEIDI p 0.05且PP.H4 0.8的基因作为最终的“优先级基因”。这套组合拳的好处是SMR负责快速筛选并给出效应方向和大小coloc负责在局部区域仔细核查两种信号是否真的共享同一变异。两者都通过的话候选基因的可信度就高多了。4. 一次完整的实战从零到候选基因的逐步拆解4.1 场景设定和初始数据准备为了让你对上面的流程有更清晰的感知我来模拟一个分析场景。假设你研究的是炎症性肠病IBD手里有一批来自公开数据库的GWAS汇总统计包含约800万个SNP。你想找到IBD相关位点中究竟哪些基因是真正的致病基因。第一步拿到GWAS数据后先看看样本信息种族构成欧洲人群就用欧洲人群的eQTL数据、病例对照比例、使用了哪些协变量。然后下载GTEx V8的结肠组织的eQTL数据因为IBD主要是肠道疾病以及一份欧洲人群的LD参考面板比如1000 Genomes的EUR队列。数据都到齐之后预处理阶段有两个细节格外重要。一是检查GWAS数据中的效应等位基因是否和eQTL数据一致不一致的话把所有等位基因都翻转到正链格式。二是确认SNP的位置使用的是哪个基因组版本GRCh37还是GRCh38GTEx V8用的是GRCh38而很多经典GWAS数据是GRCh37需要做坐标转换。我在第一次做这个转换的时候因为没有仔细检查基因组版本导致几千个SNP的坐标对不上白白排查了大半天。后来学乖了所有分析前先抽样10个位点人工对齐一遍。4.2 用SMR做全基因组扫描预处理完成后用SMR对全基因组范围内约800万个SNP扫描一遍肠组织的eQTL数据。扫描结果通常会产生几千个“探针-基因”水平的关联测试。经过Bonferroni校正后以有效的探针数量为基准大概会有几十个基因达到显著水平。在结果表里每一行展示一个基因的SMR效应量、SE、p值和HEIDI检验结果。在这一步还要关注效应量的方向。比如某个基因的SMR beta值为正说明该基因表达量升高会增加IBD风险为负则表达量降低会增加风险。这个方向信息在你后面解读生物学意义时很有用——如果直观上觉得某个基因应该是高表达致病但SMR结果显示低表达致病那就要回过去检查数据处理或者考虑是否存在反馈调控等复杂机制。4.3 对SMR显著基因做coloc精细验证以SMR扫描出的显著基因列表为基础提取这些基因所在区域的GWAS原始信号和eQTL原始信号然后用coloc逐个区域跑共定位。举个简化示例。假设SMR扫描发现GATA3基因附近区域信号显著那么以该区域最显著的GWAS SNP为中心取前后500kb的所有SNP提取这些SNP在结肠eQTL数据中对应的associations然后运行coloc。运行后结果可能显示PP.H4 0.92PP.H3 0.04说明这个区域有极强的共定位证据——GWAS的信号和eQTL的信号很可能共享同一个因果变异。这个时候基因GATA3就从“统计关联区域附近的一个基因”升级为“具有表达调控证据支持的IBD候选基因”。4.4 可视化验证和最终优先基因排序不能只看数字就下结论。我习惯把共定位区域做两个额外的可视化检查。一个是LocusZoom风格的区域关联图把GWAS的-log10(p)和eQTL的-log10(p)画在同一张图上同时标记LD颜色。如果两个信号的峰值位置几乎重合并且LD模式一致共定位的证据就非常稳固。另一个是效果量散点图横轴是各SNP对基因表达的效应量纵轴是对疾病的效应量。如果这些点沿着一条直线分布说明两种关联很可能由同一个因果变异驱动。如果散点图是杂乱的一团云即使统计上的PP.H4碰巧较高我也要留个心眼。最终把通过SMR和coloc双重验证的基因按照证据强度HEIDI p值、PP.H4、效应方向一致性排序输出一个优先级列表。这个列表就是你后续做功能实验比如在细胞系中敲低或过表达目标基因后观察表型变化的最佳起点。5. 实操中绕不开的坑以及我的应对经验5.1 LD参考面板选错全盘皆输无论是SMR的HEIDI检验还是coloc的共定位计算都需要LD信息。LD参考面板的选择直接决定了你在某个区域内考虑哪些SNP作为“相关变异”以及它们在LD结构中的关系。常用的参考面板包括1000 Genomes约2504个样本五个超级人群、UK10K、以及更大规模的TopMed。选择的原则是GWAS样本是什么人群就用什么人群的参考面板。欧洲人群GWAS就用EUR参考面板东亚人群就用EAS。一个常见的误区是嫌1000 Genomes样本量太小干脆用全人类混在一起的参考面板结果LD结构被稀释和扭曲有些真实的LD关系被抹平有些虚假的LD关系被制造出来最终HEIDI的过滤效果大打折扣。另外如果GWAS数据本身是来自千人基因组时代的经典汇总统计比如某些早期发表的IBD或类风湿关节炎研究那1000 Genomes参考面板完全够用。如果是近年的大样本GWAS可以考虑使用UK Biobank衍生的LD矩阵但要注意版本兼容问题。5.2 组织特异性用错组织的eQTL等于白跑前面提到过组织特异性的重要性这里我更详细地拆解一下为什么“随便用一个组织”的危害这么大。同一套基因型数据在不同组织里eQTL信号的强度差异可以非常悬殊。一个SNP在肝脏组织中调控某基因表达的效果量可能是0.5在脑组织中这个效应可能完全消失。这意味着如果你研究的是神经性疾病却用了GTEx的肝脏eQTL数据即使某个位点在基因组坐标上正好位于某个基因附近共定位分析也无法检测到信号。应对方案是如果没有目标组织的eQTL数据可以考虑用与该组织细胞类型相似的替代组织比如用全血替代免疫相关组织的部分场景或者使用基于多个组织整合分析的算法如Probabilistic Integration of multiple eQTL datasetsPMC等这些方法可以按组织权重整合多个来源的信号。还有一个进阶思路是使用单细胞eQTL数据。近年来像GTEx的单细胞数据扩展项目和DICEDatabase of Immune Cell Expression等资源可以提供特定免疫细胞类型的eQTL。对于免疫介导的疾病这远比全组织混匀的数据更有分辨率。5.3 多重检验与结果过拟合GWAS的检验规模是百万级别的eQTL分析的检验规模也是百万级别的两者联合分析后检验次数会叠加。如果不在每个环节做多重检验校正你的候选基因列表大概率充满假阳性。我的一般做法是SMR阶段用Bonferroni校正显著p值阈值设为0.05除以有效探针数。有效探针数可以在SMR输出的.esi文件中查到——注意不是原始探针总数而是经过QC后实际参与分析的探针数量。到了coloc阶段因为只对那些已经通过了SMR筛选的基因运行检验次数已经大幅减少但也不能掉以轻心。如果筛选出的基因太多可以考虑对PP.H4再做一次FDR控制。在实际的审稿过程中很多评审都会问“你们做了多重检验校正吗”。如果能在方法部分明确写出各阶段的校正方案整个分析的严谨性会显著提升。5.4 阴性结果也可能是正确答案最后一个想提醒的是不要为了找“显著候选基因”而硬凑结果。如果你在一个GWAS显著区域内跑了coloc发现PP.H4很低而PP.H3很高这意味着GWAS信号和eQTL信号各有自己的因果变异它们只是“碰巧住得近”。这种情况下该基因其实不太可能通过表达调控来介导疾病风险。它或许通过改变蛋白结构、影响剪接、改变非编码RNA等机制起作用但至少不是eQTL介导的。我见过有人在这种情况下偷偷把eQTL数据换成另一批直到跑出一个PP.H40.8的结果才罢休。这种做法是典型的p-hacking最后只会浪费后续功能实验的时间和经费。阴性结果也是信息它告诉你这个位点的机制不在转录调控层面可能在剪接调控、翻译调控甚至更下游的层面。5.5 一个额外的判断技巧跨组织一致性验证在完成了标准流程之后我建议额外做一步跨组织验证。如果你在GTEx的结肠组织中找到了一个IBD相关候选基因把这个基因放到GTEx的全血、肝脏、肌肉等组织中再看一遍共定位结果。如果只有结肠组织显著而其他组织都没有信号这不仅符合IBD的组织特异性预期还增强了结果的可信度。相反如果所有组织都显示共定位信号那就要考虑是不是某种普遍的调控机制或者是否存在全血细胞比例等混杂因素在驱动。这一步虽然简单但在我参与的项目里它经常是说服评审最有力的一张图。6. 从筛选到机制验证拿到候选基因之后还能做什么走到这一步你手里应该有一批“通过了SMRcoloc双重筛选”的候选基因。但这些基因只是统计层面的“嫌疑人”要真正定罪还需要机制的交叉验证。最直接的下一步是在公开的单细胞转录组数据中查看候选基因的细胞类型表达分布。比如你筛出一个免疫相关基因可以看看它是不是主要在T细胞或巨噬细胞亚群中表达这能为了解它在疾病中的角色提供第一层线索。更硬核的做法是回到分子实验。在疾病相关的细胞模型中用CRISPR干扰或小干扰RNA降低候选基因的表达观察细胞表型是否发生变化增殖、凋亡、迁移、炎症因子分泌等。如果降低该基因表达可以改变与疾病相关的细胞表型那就为基因的致病作用提供了直接的功能证据。有些人可能觉得从生信到湿实验跨度太大但现在的模式已经越来越普遍。你可以先发表方法学层面筛选的结果把功能验证作为后续工作的基础。还有一个成本更低、速度更快的验证方向是检查该基因是否已经被药物靶点数据库收录。比如在DrugBank、Open Targets等数据库中查询你的候选基因如果它已经是某些在研药物的靶点那意味着你发现的基因不仅有生物学意义还可能具备转化价值。这在写基金申请和论文引言的时候都是很好的加分项。关于工具版本和运行环境的一点提示最后补充一个实操层面的具体建议安装环境和使用版本问题。SMR、coloc这类工具虽然不像深度学习框架那么吃资源但不同版本的结果可能会有差异。比如R包coloc从3.x升级到5.x之后内部算法细节有一些调整输出结果中部分字段的命名也发生了改变。如果你的分析延续了好几个月中间升级过R包那我建议把环境信息R版本、coloc版本、GTEx数据版本全部记录在分析日志里。这些信息在论文方法部分通常也要求详细描述。另外建议所有分析脚本用sessionInfo()输出运行环境日志随结果文件一起存档。这些看似琐碎的习惯在审稿人索要数据复现说明时就是救命稻草。我自己的体会是GWAS和eQTL的整合分析并不复杂真正拉开差距的地方在于对数据的理解深度和对细节的把控能力。把每个参数的含义搞清楚把每条结果的生物学含义想明白再配合严谨的过滤条件你就能从这个经典流程里得到真正有价值的信息。 
长河君 · 2026/10/4 7:03:58 · 阅读 1.2k GWAS后如何锁定真正致病基因？eQTL共定位与SMR实战指南
我见过太多人做完GWAS拿到一张漂亮的曼哈顿图之后就卡在下一步了。最经典的操作是把显著位点拉到UCSC或Ensembl里找离它最近的基因然后宣布“找到了疾病相关基因”。但这个做法在遗传学圈子里经常翻车因为连锁不平衡LD会让信号“挂”在错误的基因旁边而eQTL就是用来解决这个问题的关键工具。把GWAS和eQTL结合起来做共定位或孟德尔随机化分析已经成为目前筛选疾病致病基因的主流策略也是我从统计信号走向生物学机制最实用的一条路径。这篇内容适合手里已经有一批GWAS显著位点、但不知道如何往下推到基因层面的研究者也适合刚接触遗传学分析、想搞懂这两个术语到底怎么用的学生。1. GWAS结果交到手上之后真正的麻烦才开始1.1 显著位点不等于致病位点更不等于致病基因GWAS输出的核心结果是一堆达到全基因组显著性水平的SNP比如p值小于5×10⁻⁸。很多初学者会本能地把这些SNP当作“致病变异”但实际情况要复杂得多。人类基因组在传递过程中不是整段完整遗传的而是以“块”为单位——这就是连锁不平衡。一个位点上的变异往往和周围一大段区域里的其他变异高度关联。这意味着GWAS捡到的“最显著SNP”未必是真正起功能作用的变异它可能只是因为和那个真正的因果变异靠得近、相关性高而被检测到。更有意思的是GWAS显著位点里有相当大比例位于基因间区或内含子区域它们本身不改变氨基酸序列却可能通过影响调控元件来改变基因表达水平。这就引出一个关键问题从统计关联到生物学机制中间隔着一条“功能注释”的鸿沟。如果你只停留在SNP层面你能回答的是“某个基因组区域与疾病相关”但无法回答“这个区域里到底是哪个基因在捣乱”。这时候就需要引入eQTL数据。eQTL做的就是给每个SNP和基因表达量之间搭一座桥某个SNP的基因型不同某个基因的表达水平是否随之显著改变。如果GWAS显著位点恰好也是一个eQTL位点那就意味着这个位点很可能通过调控特定基因的表达来影响疾病风险。1.2 “最近基因法”为什么会带偏方向我刚开始做这一行的时候也走过“离哪个基因最近就选哪个”的弯路。这种策略在某些情况下碰巧是对的但本质上没有理论依据。原因很简单三维基因组结构决定了基因表达的调控不一定遵守“距离最近”原则。一个位于基因A下游50kb的变异完全可能通过染色质环化作用与远处的基因B的启动子区域接触从而调控基因B的表达。在二维基因组坐标上看起来“近”的基因在三维空间里可能压根儿不相关。业内已经有不少案例展示了这种误导性。比如某些GWAS位点落在基因间区按最近距离归属于基因C但后续功能实验证明真正受调控的是位于数百kb之外的基因D。如果你只靠距离来判断后续所有实验敲除、过表达、动物模型都会押错筹码。所以GWAS之后加一步eQTL共定位分析本质上是在用功能证据校正统计信号的归属问题。1.3 关联、共定位与因果推断三者不是一回事在进入实操之前有几组概念必须分清楚不然你在阅读文献和跑工具的时候会很痛苦。关联GWAS发现的“这个SNP和这个疾病在统计上相关”不涉及方向性和机制。共定位在某个基因组区域内GWAS信号和eQTL信号是否由同一个因果变异驱动。如果共定位后验概率高说明这个位点的基因型和疾病之间的关系很可能通过基因表达来介导。因果推断更进一步用孟德尔随机化的思想以基因型为工具变量评估基因表达水平对疾病风险的因果效应。从关联到共定位再到因果推断是证据强度逐级提升的过程。我见过有的人拿到coloc结果后直接说“证实了XX基因导致XX疾病”其实coloc只是支持“存在共同因果变异”的证据真正的因果结论还需要结合MR结果、功能实验甚至动物模型才能站得住脚。2. eQTL数据的底层逻辑为什么它能帮GWAS“指认”基因2.1 从DNA到RNA之间的调控层eQTL的全称是expression Quantitative Trait Loci翻译过来是“表达数量性状基因座”。它研究的是一个非常简单的问题不同个体的基因组序列有差异这些差异是否会导致基因表达水平的差异。举个例子假设有1000个人某个SNP位点上有的人带A等位基因有的人带G等位基因。如果带A的人群体中某个基因的平均表达量显著高于带G的人群那么这个SNP就是这个基因的一个eQTL。这里的“表达量”指的是RNA水平通常用RNA-seq数据来量化。也就是说eQTL分析需要两类数据基因型数据SNP分型和基因表达数据RNA-seq而且这两类数据必须来自同一批个体。这一层机制非常重要因为很多疾病相关变异并不会改变蛋白质的氨基酸序列而是通过影响基因的转录调控来发挥作用。比如改变了转录因子结合位点、影响了启动子活性、或者干扰了增强子与启动子的相互作用。eQTL分析正是把这类调控效应转化为可量化的统计关联。2.2 cis-eQTL和trans-eQTL为什么优先看ciseQTL按照作用距离分为两类cis-eQTLSNP调控的是它自身附近通常定义为1Mb范围内的基因表达。这种调控通常比较直接比如SNP位于某个基因的启动子区域直接影响该基因的转录效率。trans-eQTLSNP调控的是距离很远的基因甚至跨染色体。这种调控往往通过间接机制实现比如SNP影响了某个转录因子的表达这个转录因子又去调控下游一大批基因。在实际的GWAS-eQTL联合分析中绝大多数人只看cis-eQTL。原因有几个方面。第一cis-eQTL的效应通常更强、更稳定跨数据集的可重复性更高。第二cis-eQTL信号的LD结构相对清晰容易做共定位分析。第三从生物学解释的角度看一个位点顺式调控它旁边的基因这个故事讲起来更合理后续功能验证也更容易设计。相比之下trans-eQTL存在严重得多重检验问题一个SNP要测试对所有基因的影响而且远距离调控关系容易受混杂因素干扰。GTEx项目的统计结果显示cis-eQTL在几乎所有组织中都能检测到大量信号而trans-eQTL就稀少得多。因此除非你有非常明确的理由要研究trans调控否则建议把精力集中在cis-eQTL上。2.3 主流eQTL数据资源及选择标准目前eQTL分析可用的公开数据主要来自几个大型项目每个项目有自己的人群特征和组织覆盖范围选择时不能闭着眼睛抓一个就用。数据资源覆盖组织样本规模适合场景GTExV849个组织~838个供体多组织差异分析器官特异性疾病首选eQTLGen Consortium全血~31,684人血液相关性状或大样本验证DGN全血~922人有RNA-seq原始数据可自定义分析TCGA多种癌症组织数百至上千例癌症eQTL分析结合肿瘤表达谱Brain eQTL如MetaBrain、ROSMAP脑组织多个区域数千例神经精神类疾病选择标准我个人的经验是病在哪个组织就优先用哪个组织的eQTL数据。比如研究冠心病用GTEx的冠状动脉组织或心脏组织研究阿尔茨海默病一定要找脑组织来源的eQTL数据而不是偷懒用全血的。因为eQTL的组织特异性非常强一个SNP在肝脏里可能是某个基因的eQTL在脑子里可能完全没有调控效应。如果你用了错误的组织数据共定位结果大概率是阴性。如果疾病相关组织没有足够的公开eQTL数据退而求其次可以用全血但必须明确说明这是替代方案并且最好做一个跨组织的敏感性分析在不同组织的数据里看结果是否稳定。3. 把两套信号放在一起共定位与SMR的本质区别和执行细节3.1 coloc和SMR/HEIDI分别回答什么问题GWAS和eQTL都拿到之后联合分析的主流工具有两个流派Bayesian共定位代表工具是coloc和基于工具变量的回归方法代表工具是SMR/HEIDI。coloc的核心思想是在一个区域内同时观察GWAS关联信号和eQTL关联信号评估这两种信号是否共享同一个因果变异。它输出的是五个后验概率PP.H0该区域内没有任何因果变异PP.H1只有GWAS信号有因果变异PP.H2只有eQTL信号有因果变异PP.H3两个信号有各自的因果变异但它们不共享PP.H4两个信号共享同一个因果变异我们通常关注PP.H4如果它大于0.75或0.8说明共定位证据充分。**SMRSummary-data-based Mendelian Randomization**走的是另一条思路。它把eQTL的效应量SNP对表达的影响当作暴露把GWAS的效应量SNP对疾病的影响当作结局用类似孟德尔随机化的框架去估算“基因表达每变化一个单位疾病风险变化多少”。但SMR存在一个问题两个信号可能只是LD造成的关系而不是真正的因果介导。因此SMR内置了一个HEIDI过滤步骤用来检测关联信号是否是由于LD造成的而不是同一个因果变异驱动的。这两者不是竞争关系而是互补关系。我常用的策略是先用SMR做全基因组范围内的快速扫描得到一个候选基因列表然后用coloc对每个候选位点做精细的共定位验证。3.2 coloc操作的完整流程与参数设置关于coloc我需要给你一个可以直接上手的实操路径。假设你已经准备好了GWAS的sumstats文件包含SNP、效应等位基因、效应量、标准误、p值、等位基因频率和eQTL的sumstats文件流程如下。第一步是数据格式统一。coloc要求输入每个SNP的效应量beta、方差通常用标准误平方、MAF、SNP ID。两个数据集的SNP ID必须一致建议都使用chr_pos_allele格式避免rsID在不同版本注释中的混乱问题。第二步是区域划分。coloc假设每个区域只有一个因果变异所以不能把整条染色体丢进去。一般做法是以GWAS最显著SNP为中心取前后500kb或1Mb的区域。如果你用的是PLINK格式数据可以用--ld-snp配合--r2来界定独立信号区间。第三步是运行coloc。核心代码我这里给一个R语言的骨架library(coloc) # 准备数据 dataset1 - list(beta gwas$beta, varbeta gwas$se^2, snp gwas$snp, type cc, s 0.3, MAF gwas$maf) # s参数代表病例占比如果是定量性状就设置typequant dataset2 - list(beta eqtl$beta, varbeta eqtl$se^2, snp eqtl$snp, type quant, MAF eqtl$maf) # 运行共定位 my_coloc - coloc.abf(dataset1, dataset2, p12 1e-5) print(my_coloc$summary)这里有一个容易踩的坑p12参数的设置。p12是先验概率表示该区域同时存在GWAS和eQTL因果变异的概率。默认值是1e-5但如果你筛查的区域很多适当调低一点可以控制假阳性。我一般跑批量分析时用1e-6单个候选位点验证时用默认值。还有一个细节需要注意如果GWAS是case-control设计必须正确设置s参数人群中的患病率它会直接影响PP.H4的计算。很多人忽略这个参数得到的后验概率会偏高或偏低。3.3 SMR操作的参数细节mBIM和HEIDI阈值的正确理解SMR的操作相对简洁。下载SMR软件之后通常用以下命令格式smr --bfile gwas_plink --gwas-summary gwas.ma --beqtl-summary eqtl_merged --out smr_result --thread-num 10--bfileGWAS数据对应的参考面板PLINK文件用于计算LD--gwas-summaryGWAS汇总统计文件--beqtl-summaryeQTL数据的二进制格式SMR输出的每个探针通常对应一个基因会有一个p值。判断显著位点时常用的阈值是p小于0.05除以有效探针数Bonferroni校正。但真正需要关注的是HEIDI检验。HEIDI的p值如果大于0.05说明没有显著证据表明该关联是由LD造成的可以认为SMR信号是真实的因果效应。如果HEIDI p值很小小于0.01那就要警惕了——这个信号很可能只是LD混淆造成的假象。我见过很多人在SMR结果里挑p值很显著的基因却不看HEIDI的过滤结果结果选了一堆假阳性。SMR官方文档的建议是优先保留HEIDI p 0.05且SMR p达到显著阈值的基因这个组合才是比较可信的候选。3.4 从SMR到coloc的筛选流水线实际操作中我不建议把coloc直接跑在几百万个位点上计算量大且多重检验问题很棘手。推荐的做法是先粗后精全基因组范围内运行SMR得到初步候选基因列表。对每个候选基因的所在区域以最显著SNP为中心±500kb提取GWAS和eQTL数据。对每个区域运行coloc记录PP.H4。保留同时满足HEIDI p 0.05且PP.H4 0.8的基因作为最终的“优先级基因”。这套组合拳的好处是SMR负责快速筛选并给出效应方向和大小coloc负责在局部区域仔细核查两种信号是否真的共享同一变异。两者都通过的话候选基因的可信度就高多了。4. 一次完整的实战从零到候选基因的逐步拆解4.1 场景设定和初始数据准备为了让你对上面的流程有更清晰的感知我来模拟一个分析场景。假设你研究的是炎症性肠病IBD手里有一批来自公开数据库的GWAS汇总统计包含约800万个SNP。你想找到IBD相关位点中究竟哪些基因是真正的致病基因。第一步拿到GWAS数据后先看看样本信息种族构成欧洲人群就用欧洲人群的eQTL数据、病例对照比例、使用了哪些协变量。然后下载GTEx V8的结肠组织的eQTL数据因为IBD主要是肠道疾病以及一份欧洲人群的LD参考面板比如1000 Genomes的EUR队列。数据都到齐之后预处理阶段有两个细节格外重要。一是检查GWAS数据中的效应等位基因是否和eQTL数据一致不一致的话把所有等位基因都翻转到正链格式。二是确认SNP的位置使用的是哪个基因组版本GRCh37还是GRCh38GTEx V8用的是GRCh38而很多经典GWAS数据是GRCh37需要做坐标转换。我在第一次做这个转换的时候因为没有仔细检查基因组版本导致几千个SNP的坐标对不上白白排查了大半天。后来学乖了所有分析前先抽样10个位点人工对齐一遍。4.2 用SMR做全基因组扫描预处理完成后用SMR对全基因组范围内约800万个SNP扫描一遍肠组织的eQTL数据。扫描结果通常会产生几千个“探针-基因”水平的关联测试。经过Bonferroni校正后以有效的探针数量为基准大概会有几十个基因达到显著水平。在结果表里每一行展示一个基因的SMR效应量、SE、p值和HEIDI检验结果。在这一步还要关注效应量的方向。比如某个基因的SMR beta值为正说明该基因表达量升高会增加IBD风险为负则表达量降低会增加风险。这个方向信息在你后面解读生物学意义时很有用——如果直观上觉得某个基因应该是高表达致病但SMR结果显示低表达致病那就要回过去检查数据处理或者考虑是否存在反馈调控等复杂机制。4.3 对SMR显著基因做coloc精细验证以SMR扫描出的显著基因列表为基础提取这些基因所在区域的GWAS原始信号和eQTL原始信号然后用coloc逐个区域跑共定位。举个简化示例。假设SMR扫描发现GATA3基因附近区域信号显著那么以该区域最显著的GWAS SNP为中心取前后500kb的所有SNP提取这些SNP在结肠eQTL数据中对应的associations然后运行coloc。运行后结果可能显示PP.H4 0.92PP.H3 0.04说明这个区域有极强的共定位证据——GWAS的信号和eQTL的信号很可能共享同一个因果变异。这个时候基因GATA3就从“统计关联区域附近的一个基因”升级为“具有表达调控证据支持的IBD候选基因”。4.4 可视化验证和最终优先基因排序不能只看数字就下结论。我习惯把共定位区域做两个额外的可视化检查。一个是LocusZoom风格的区域关联图把GWAS的-log10(p)和eQTL的-log10(p)画在同一张图上同时标记LD颜色。如果两个信号的峰值位置几乎重合并且LD模式一致共定位的证据就非常稳固。另一个是效果量散点图横轴是各SNP对基因表达的效应量纵轴是对疾病的效应量。如果这些点沿着一条直线分布说明两种关联很可能由同一个因果变异驱动。如果散点图是杂乱的一团云即使统计上的PP.H4碰巧较高我也要留个心眼。最终把通过SMR和coloc双重验证的基因按照证据强度HEIDI p值、PP.H4、效应方向一致性排序输出一个优先级列表。这个列表就是你后续做功能实验比如在细胞系中敲低或过表达目标基因后观察表型变化的最佳起点。5. 实操中绕不开的坑以及我的应对经验5.1 LD参考面板选错全盘皆输无论是SMR的HEIDI检验还是coloc的共定位计算都需要LD信息。LD参考面板的选择直接决定了你在某个区域内考虑哪些SNP作为“相关变异”以及它们在LD结构中的关系。常用的参考面板包括1000 Genomes约2504个样本五个超级人群、UK10K、以及更大规模的TopMed。选择的原则是GWAS样本是什么人群就用什么人群的参考面板。欧洲人群GWAS就用EUR参考面板东亚人群就用EAS。一个常见的误区是嫌1000 Genomes样本量太小干脆用全人类混在一起的参考面板结果LD结构被稀释和扭曲有些真实的LD关系被抹平有些虚假的LD关系被制造出来最终HEIDI的过滤效果大打折扣。另外如果GWAS数据本身是来自千人基因组时代的经典汇总统计比如某些早期发表的IBD或类风湿关节炎研究那1000 Genomes参考面板完全够用。如果是近年的大样本GWAS可以考虑使用UK Biobank衍生的LD矩阵但要注意版本兼容问题。5.2 组织特异性用错组织的eQTL等于白跑前面提到过组织特异性的重要性这里我更详细地拆解一下为什么“随便用一个组织”的危害这么大。同一套基因型数据在不同组织里eQTL信号的强度差异可以非常悬殊。一个SNP在肝脏组织中调控某基因表达的效果量可能是0.5在脑组织中这个效应可能完全消失。这意味着如果你研究的是神经性疾病却用了GTEx的肝脏eQTL数据即使某个位点在基因组坐标上正好位于某个基因附近共定位分析也无法检测到信号。应对方案是如果没有目标组织的eQTL数据可以考虑用与该组织细胞类型相似的替代组织比如用全血替代免疫相关组织的部分场景或者使用基于多个组织整合分析的算法如Probabilistic Integration of multiple eQTL datasetsPMC等这些方法可以按组织权重整合多个来源的信号。还有一个进阶思路是使用单细胞eQTL数据。近年来像GTEx的单细胞数据扩展项目和DICEDatabase of Immune Cell Expression等资源可以提供特定免疫细胞类型的eQTL。对于免疫介导的疾病这远比全组织混匀的数据更有分辨率。5.3 多重检验与结果过拟合GWAS的检验规模是百万级别的eQTL分析的检验规模也是百万级别的两者联合分析后检验次数会叠加。如果不在每个环节做多重检验校正你的候选基因列表大概率充满假阳性。我的一般做法是SMR阶段用Bonferroni校正显著p值阈值设为0.05除以有效探针数。有效探针数可以在SMR输出的.esi文件中查到——注意不是原始探针总数而是经过QC后实际参与分析的探针数量。到了coloc阶段因为只对那些已经通过了SMR筛选的基因运行检验次数已经大幅减少但也不能掉以轻心。如果筛选出的基因太多可以考虑对PP.H4再做一次FDR控制。在实际的审稿过程中很多评审都会问“你们做了多重检验校正吗”。如果能在方法部分明确写出各阶段的校正方案整个分析的严谨性会显著提升。5.4 阴性结果也可能是正确答案最后一个想提醒的是不要为了找“显著候选基因”而硬凑结果。如果你在一个GWAS显著区域内跑了coloc发现PP.H4很低而PP.H3很高这意味着GWAS信号和eQTL信号各有自己的因果变异它们只是“碰巧住得近”。这种情况下该基因其实不太可能通过表达调控来介导疾病风险。它或许通过改变蛋白结构、影响剪接、改变非编码RNA等机制起作用但至少不是eQTL介导的。我见过有人在这种情况下偷偷把eQTL数据换成另一批直到跑出一个PP.H40.8的结果才罢休。这种做法是典型的p-hacking最后只会浪费后续功能实验的时间和经费。阴性结果也是信息它告诉你这个位点的机制不在转录调控层面可能在剪接调控、翻译调控甚至更下游的层面。5.5 一个额外的判断技巧跨组织一致性验证在完成了标准流程之后我建议额外做一步跨组织验证。如果你在GTEx的结肠组织中找到了一个IBD相关候选基因把这个基因放到GTEx的全血、肝脏、肌肉等组织中再看一遍共定位结果。如果只有结肠组织显著而其他组织都没有信号这不仅符合IBD的组织特异性预期还增强了结果的可信度。相反如果所有组织都显示共定位信号那就要考虑是不是某种普遍的调控机制或者是否存在全血细胞比例等混杂因素在驱动。这一步虽然简单但在我参与的项目里它经常是说服评审最有力的一张图。6. 从筛选到机制验证拿到候选基因之后还能做什么走到这一步你手里应该有一批“通过了SMRcoloc双重筛选”的候选基因。但这些基因只是统计层面的“嫌疑人”要真正定罪还需要机制的交叉验证。最直接的下一步是在公开的单细胞转录组数据中查看候选基因的细胞类型表达分布。比如你筛出一个免疫相关基因可以看看它是不是主要在T细胞或巨噬细胞亚群中表达这能为了解它在疾病中的角色提供第一层线索。更硬核的做法是回到分子实验。在疾病相关的细胞模型中用CRISPR干扰或小干扰RNA降低候选基因的表达观察细胞表型是否发生变化增殖、凋亡、迁移、炎症因子分泌等。如果降低该基因表达可以改变与疾病相关的细胞表型那就为基因的致病作用提供了直接的功能证据。有些人可能觉得从生信到湿实验跨度太大但现在的模式已经越来越普遍。你可以先发表方法学层面筛选的结果把功能验证作为后续工作的基础。还有一个成本更低、速度更快的验证方向是检查该基因是否已经被药物靶点数据库收录。比如在DrugBank、Open Targets等数据库中查询你的候选基因如果它已经是某些在研药物的靶点那意味着你发现的基因不仅有生物学意义还可能具备转化价值。这在写基金申请和论文引言的时候都是很好的加分项。关于工具版本和运行环境的一点提示最后补充一个实操层面的具体建议安装环境和使用版本问题。SMR、coloc这类工具虽然不像深度学习框架那么吃资源但不同版本的结果可能会有差异。比如R包coloc从3.x升级到5.x之后内部算法细节有一些调整输出结果中部分字段的命名也发生了改变。如果你的分析延续了好几个月中间升级过R包那我建议把环境信息R版本、coloc版本、GTEx数据版本全部记录在分析日志里。这些信息在论文方法部分通常也要求详细描述。另外建议所有分析脚本用sessionInfo()输出运行环境日志随结果文件一起存档。这些看似琐碎的习惯在审稿人索要数据复现说明时就是救命稻草。我自己的体会是GWAS和eQTL的整合分析并不复杂真正拉开差距的地方在于对数据的理解深度和对细节的把控能力。把每个参数的含义搞清楚把每条结果的生物学含义想明白再配合严谨的过滤条件你就能从这个经典流程里得到真正有价值的信息。