企业官网建设流程全解析

我见过太多人做完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：为什么优先看cis

eQTL按照作用距离分为两类：

• cis-eQTL：SNP调控的是它自身附近（通常定义为1Mb范围内）的基因表达。这种调控通常比较直接，比如SNP位于某个基因的启动子区域，直接影响该基因的转录效率。
• trans-eQTL：SNP调控的是距离很远的基因，甚至跨染色体。这种调控往往通过间接机制实现，比如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分析可用的公开数据主要来自几个大型项目，每个项目有自己的人群特征和组织覆盖范围，选择时不能闭着眼睛抓一个就用。

选择标准我个人的经验是：病在哪个组织，就优先用哪个组织的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，说明共定位证据充分。

**SMR（Summary-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参数代表病例占比，如果是定量性状就设置type="quant" 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

• --bfile：GWAS数据对应的参考面板PLINK文件，用于计算LD
• --gwas-summary：GWAS汇总统计文件
• --beqtl-summary：eQTL数据的二进制格式

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直接跑在几百万个位点上，计算量大且多重检验问题很棘手。推荐的做法是先粗后精：

1. 全基因组范围内运行SMR，得到初步候选基因列表。
2. 对每个候选基因的所在区域（以最显著SNP为中心±500kb）提取GWAS和eQTL数据。
3. 对每个区域运行coloc，记录PP.H4。
4. 保留同时满足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还是GRCh38），GTEx 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.92，PP.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 datasets，PMC等），这些方法可以按组织权重整合多个来源的信号。

还有一个进阶思路是使用单细胞eQTL数据。近年来，像GTEx的单细胞数据扩展项目和DICE（Database 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.H4>0.8的结果才罢休。这种做法是典型的p-hacking，最后只会浪费后续功能实验的时间和经费。阴性结果也是信息：它告诉你这个位点的机制不在转录调控层面，可能在剪接调控、翻译调控甚至更下游的层面。

5.5 一个额外的判断技巧：跨组织一致性验证

在完成了标准流程之后，我建议额外做一步跨组织验证。

如果你在GTEx的结肠组织中找到了一个IBD相关候选基因，把这个基因放到GTEx的全血、肝脏、肌肉等组织中再看一遍共定位结果。如果只有结肠组织显著，而其他组织都没有信号，这不仅符合IBD的组织特异性预期，还增强了结果的可信度。相反，如果所有组织都显示共定位信号，那就要考虑是不是某种普遍的调控机制，或者是否存在全血细胞比例等混杂因素在驱动。

这一步虽然简单，但在我参与的项目里，它经常是说服评审最有力的一张图。

6. 从筛选到机制验证：拿到候选基因之后还能做什么

走到这一步，你手里应该有一批“通过了SMR+coloc双重筛选”的候选基因。但这些基因只是统计层面的“嫌疑人”，要真正定罪还需要机制的交叉验证。

最直接的下一步是在公开的单细胞转录组数据中查看候选基因的细胞类型表达分布。比如你筛出一个免疫相关基因，可以看看它是不是主要在T细胞或巨噬细胞亚群中表达，这能为了解它在疾病中的角色提供第一层线索。

更硬核的做法是回到分子实验。在疾病相关的细胞模型中，用CRISPR干扰或小干扰RNA降低候选基因的表达，观察细胞表型是否发生变化（增殖、凋亡、迁移、炎症因子分泌等）。如果降低该基因表达可以改变与疾病相关的细胞表型，那就为基因的致病作用提供了直接的功能证据。

有些人可能觉得从生信到湿实验跨度太大，但现在的模式已经越来越普遍。你可以先发表方法学层面（筛选）的结果，把功能验证作为后续工作的基础。

还有一个成本更低、速度更快的验证方向是检查该基因是否已经被药物靶点数据库收录。比如在DrugBank、Open Targets等数据库中查询你的候选基因，如果它已经是某些在研药物的靶点，那意味着你发现的基因不仅有生物学意义，还可能具备转化价值。这在写基金申请和论文引言的时候都是很好的加分项。

关于工具版本和运行环境的一点提示

最后补充一个实操层面的具体建议：安装环境和使用版本问题。

SMR、coloc这类工具虽然不像深度学习框架那么吃资源，但不同版本的结果可能会有差异。比如R包coloc从3.x升级到5.x之后，内部算法细节有一些调整，输出结果中部分字段的命名也发生了改变。如果你的分析延续了好几个月，中间升级过R包，那我建议把环境信息（R版本、coloc版本、GTEx数据版本）全部记录在分析日志里。这些信息在论文方法部分通常也要求详细描述。

另外，建议所有分析脚本用sessionInfo()输出运行环境，日志随结果文件一起存档。这些看似琐碎的习惯，在审稿人索要数据复现说明时就是救命稻草。

我自己的体会是：GWAS和eQTL的整合分析并不复杂，真正拉开差距的地方在于对数据的理解深度和对细节的把控能力。把每个参数的含义搞清楚，把每条结果的生物学含义想明白，再配合严谨的过滤条件，你就能从这个经典流程里得到真正有价值的信息。