HIC测序数据生信分析——第三节，HIC数据挂载实战：ALLHiC与3D-DNA双路径解析

1. 从Hi-C数据到染色体为什么需要“挂载”你好我是老张在基因组组装这个行当里摸爬滚打了十来年。今天咱们接着聊Hi-C数据分析的硬核实战部分——数据挂载。你可能已经完成了Hi-C数据的预处理拿到了一堆比对好的文件比如那个关键的.bam文件。但看着这些文件你可能会有点懵这离拿到染色体级别的基因组序列还差多远呢简单来说数据挂载就是利用Hi-C数据中蕴含的“空间临近”信息把前期组装出来的一堆碎片化的contig你可以理解为基因组拼图的小碎片按照它们在真实染色体上的物理位置和顺序给“挂”起来、连起来最终拼成一条条完整的染色体。这个过程也叫染色体构象捕获辅助组装或者支架构建。为什么这步这么关键因为二代或三代测序技术组装出来的基因组往往是由成千上万个contig组成的它们就像一本被撕成无数碎片的书你知道每一片上写了什么字碱基序列但不知道这些碎片在原书里的顺序和位置。Hi-C技术提供的“碎片A和碎片B在细胞核里经常挨在一起”的信息就是帮助我们恢复这本书原貌的“胶水”。目前主流的“胶水”工具有两个ALLHiC和3D-DNA。它们思路不同适用场景也各有侧重。这篇文章我就带你手把手走通这两条路结合我踩过的坑和实战经验帮你根据手头的数据情况选出最适合你的那把“瑞士军刀”。2. 路径一ALLHiC——应对复杂基因组的利器ALLHiC 这个名字直译过来就是“等位基因Hi-C”它的设计初衷就是为了解决多倍体或高杂合度基因组的组装难题。比如很多农作物、鱼类它们的基因组里包含多套非常相似的染色体等位基因普通的挂载工具很容易把这些相似序列的错误连接当成正确的空间临近信号导致组装混乱。ALLHiC 通过一个巧妙的“修剪”步骤专门来处理这种复杂性。2.1 核心思想与六大模块ALLHiC 把整个挂载过程拆解成六个逻辑清晰的模块像一条流水线prune修剪这是ALLHiC的特色步骤目的是利用近缘物种的基因注释信息识别并“剪除”contig之间因等位基因相似性而非真实三维空间临近产生的虚假Hi-C连接信号。这一步需要你提供一个近缘物种的CDS序列和GFF3注释文件。partition分区基于Hi-C互作信号将所有的contig初步划分到不同的组里理想情况下每个组对应一条染色体。你需要告诉它你预计的染色体条数-k参数。rescue救援在分区后有些contig可能因为信号弱没有被划分到任何组。这个模块会再次扫描数据尝试将这些“孤儿”contig归入合适的染色体组。optimize优化在每个染色体组内部调整contig的排列顺序和方向正链或反链使得最终的Hi-C互作矩阵看起来最“顺滑”、最符合预期。build构建将优化好顺序和方向的contig连接起来生成最终的染色体序列文件.fasta和描述连接关系的AGP文件。plot绘图生成可视化图让你直观地检查每条染色体的Hi-C互作热图评估组装质量。是不是感觉有点抽象别急我们直接上实战命令边操作边理解。2.2 实战演练基于Hicup预处理结果的ALLHiC全流程假设你已经用Hicup软件处理好了数据得到了clean.bam文件并且有一个初步组装的genome.facontig序列。我们的物种染色体数是10-k 10使用的限制性内切酶是DpnII识别位点GATC。第一步环境与数据准备我强烈推荐使用Docker来运行ALLHiC能避免一大堆依赖环境的问题。下面是我常用的命令# 拉取ALLHiC镜像并启动容器 docker run -it -d -v /本地/数据路径:/容器内/工作路径 --name allhic_work wangnan9394/allhic /bin/bash # 进入容器内部 docker exec -it allhic_work /bin/bash # 在容器内设置ALLHiC脚本路径假设已安装在/home/ALLHiC-master/ export PATH/home/ALLHiC-master/scripts/:/home/ALLHiC-master/bin/:$PATH进入工作目录准备好输入文件hicup.bamHicup输出的比对文件。genome.fa你的contig基因组文件。如果需要运行prune近缘物种的reference.cds.fasta和reference.gff3。第二步运行prune模块针对复杂基因组如果你研究的物种是多倍体或杂合度很高强烈建议走这一步。这里我演示基于GMAP的方法它比基于BLAST的方法更简便。# 1. 使用GMAP将近缘物种CDS序列比对到你的contig上生成GFF3 gmap_build -D . -d gmap_db genome.fa gmap -D . -d gmap_db -t 20 -f 2 -n 2 reference.cds.fasta gmap_aligned.gff3 # 注意-n 参数指定基因组倍性二倍体为2请根据实际情况调整。 # 2. 利用自定义脚本如修改自gmap2AlleleTable.pl的脚本生成等位基因关联表 perl my_allele_table.pl reference.gff3 gmap_aligned.gff3 Allelic.ctg.table # 3. 运行ALLHiC_prune ALLHiC_prune -i Allelic.ctg.table -b hicup.bam -r genome.fa运行成功后你会得到pruned.bam和pruned.clean.bam等文件其中虚假的等位基因互作信号已被削弱。第三步分区与优化核心步骤接下来无论是否进行了prune核心流程都是一样的。我们以使用prune后生成的pruned.bam为例。# 1. partition 分区 ALLHiC_partition -b pruned.bam -r genome.fa -e GATC -k 10 -m 25 # -e: 酶切类型必须与实验一致 # -k: 染色体数目务必设置准确。 # -m: 每个contig最少酶切位点数过滤掉位点太少的不可靠contig。运行后会生成pruned.counts_GATC.txt、pruned.clm以及一系列.counts_GATC.10gX.txt文件X从1到10每个文件对应一条初步划分的染色体组。# 2. rescue 救援如果运行了prune则必须运行rescue ALLHiC_rescue -b hicup.bam -r genome.fa -c pruned.clusters.txt -i pruned.counts_GATC.txt这一步会尝试“拯救”未被分区的contig生成最终的group1.txt到group10.txt。# 3. optimize 优化对每条染色体组单独进行 for i in {1..10}; do allhic optimize group${i}.txt pruned.clm done优化后每个groupX.txt文件中的contig顺序和方向就被调整到最优状态。第四步构建与可视化# 4. build 构建染色体 ALLHiC_build genome.fa最关键的两个结果文件诞生了groups.asm.fasta染色体级别的组装序列。groups.agp详细记录了每个contig在染色体中的位置、方向和间隙。# 5. plot 绘图检查效果 # 先获取染色体长度列表 seqkit fx2tab -n -l groups.asm.fasta chr_length.txt grep group chr_length.txt chrn.list # 生成Hi-C热图 ALLHiC_plot pruned.bam groups.agp chrn.list 500k pdf打开生成的500K_all_chrs.pdf你会看到每条染色体的Hi-C热图。理想情况下对角线上的信号应该最强自身互作且整个矩阵看起来干净、有条理。如果出现明显的错位或强烈的非对角线信号块可能意味着有错误连接。注意原始文章里提到运行某些步骤后“没有拿到文件说明Hi-C数据不好”。这确实是常见问题。数据质量是根基如果测序深度不足、酶切效率低或背景噪音高再好的工具也难为无米之炊。遇到这种情况你需要回头检查预处理步骤或者考虑增加测序数据量。3. 路径二3D-DNA——快速高效的“组装抛光器”如果说ALLHiC是精雕细琢的“手术刀”那么3D-DNA更像是高效快捷的“自动组装流水线”。它由Aiden实验室就是开发Hi-C技术的那个实验室开发与Juicer预处理流程是天作之合特别适合在已有较好参考基因组或初步支架的基础上进行快速、自动化的染色体挂载和错误修正。3.1 工作原理与流程特点3D-DNA的核心算法基于一个叫做“迭代校正”的策略。它不严格区分等位基因而是专注于利用全基因组范围的Hi-C互作信号通过多次迭代来优化contig的顺序、方向和位置。它的输入输出非常直接输入Juicer软件输出的merged_nodups.txt文件记录了所有有效的Hi-C互作对以及你的contig序列genome.fa。输出最终的组装文件.fasta、描述文件.assembly以及可供Juicebox可视化软件手动编辑的.hic文件。最大的亮点在于它和Juicebox的深度整合。3D-DNA生成的初步结果可以导入Juicebox进行人工浏览和调整。你可以直观地看到Hi-C热图手动拆分错误的连接、调整contig顺序然后将修改保存后重新跑一遍3D-DNA的后续步骤得到最终结果。这种人机交互的纠错能力对于处理复杂区域或数据质量稍差的情况是无可替代的。3.2 实战演练基于Juicer结果的3D-DNA流程假设你已经用Juicer跑完了预处理在aligned目录下找到了宝贵的merged_nodups.txt文件。第一步软件安装3D-DNA的安装相对简单从GitHub克隆即可。cd /your/soft/path git clone https://github.com/theaidenlab/3d-dna.git cd 3d-dna # 尝试运行帮助命令检查是否安装成功 /bin/bash ./run-asm-pipeline.sh -h第二步一键式运行3D-DNA的魅力就在于其命令的简洁。一条命令即可启动整个自动化流程。cd /your/work/path /bin/bash /path/to/3d-dna/run-asm-pipeline.sh \ -r 2 \ ../Genome/genome.fa \ ../aligned/merged_nodups.txt这里解释一下关键参数-r 2指定运行模式。2是标准组装模式。如果是第一次运行或数据质量好可以试试0快速模式。如果组装效果不佳可以增加迭代次数比如-r 3。后面两个参数就是你的contig文件和Juicer输出文件。运行开始后你会看到屏幕上滚动大量的日志信息它正在自动进行比对矩阵构建、迭代校正、生成支架等步骤。这个过程可能需要一些时间取决于数据量大小。第三步结果解读与人工干预运行结束后当前目录下会生成大量文件我们重点关注这几个*.hic这是Hi-C互作矩阵文件可以直接用Juicebox软件打开进行可视化检查。*.assembly这是一个文本文件记录了最终的组装结构。这个文件是可以手动编辑的*.fasta最终的基因组序列文件。如何人工干预这是3D-DNA流程的精髓。用Juicebox打开.hic文件你会看到一个全基因组的Hi-C热图。仔细观察。正确的组装热图应该呈现出清晰的“棋盘格”状对角线明亮。如果你看到某些区域有强烈的“非对角线”信号块比如在热图左上角到右下角的一条斜线之外出现亮块这通常意味着两个本不临近的基因组区域被错误地连接在了一起。在Juicebox中你可以直接框选这些区域进行“拆分”、“移动”、“反转”等操作。所有的编辑操作都会同步修改.assembly文件。人工编辑满意后保存.assembly文件。然后重新运行3D-DNA的最后一步命令基于你修改后的.assembly文件来生成最终的.fasta序列文件。具体命令可以参考3D-DNA文档中的run-asm-pipeline-post.sh脚本。提示对于初学者我建议先让3D-DNA全自动跑一遍得到一个基础结果。然后花时间学习使用Juicebox查看结果。即使不进行手动编辑仅仅通过可视化你也能对基因组组装的质量有一个非常直观的认识这是任何文本报告都无法替代的。4. 双剑合璧如何根据你的数据选择最佳路径讲完了两条独立路径你可能会问我到底该选哪个别急我帮你列个表从几个关键维度对比一下特性维度ALLHiC3D-DNA核心优势专门处理高杂合、多倍体能有效区分等位基因组装精度高。流程自动化程度高与Juicer/Juicebox生态无缝集成支持人工可视化纠错。输入要求需要clean.bam来自Hicup或ALLHiC预处理和genome.fa。如需prune还需近缘物种注释。需要Juicer输出的merged_nodups.txt和genome.fa。适用场景植物多倍体常见、某些鱼类、高杂合度动植物基因组。数据质量要求相对宽容。哺乳动物、模式生物、或已有较好参考基因组的物种。追求快速、可交互的流程。操作复杂度步骤多模块化需手动按顺序运行对用户理解流程要求高。几乎一键运行但后期人工纠错环节有学习成本。结果可靠性算法针对复杂基因组设计在适用场景下结果更可靠。自动化结果可能包含错误但结合Juicebox手动修正后可以达到极高精度。可视化自带绘图模块生成静态PDF热图用于质量评估。生成.hic文件需用Juicebox进行动态、交互式查看和编辑功能强大。我的实战选择建议如果你的物种是已知的多倍体比如小麦、棉花或者杂合度特别高比如许多果树别犹豫优先选择ALLHiC。它的prune模块是解决这类问题的“特效药”。虽然步骤繁琐点但能从根本上减少错误连接。如果你做的是人类、小鼠、果蝇等模式生物或者你有非常近缘的高质量参考基因组那么3D-DNA会是更高效的选择。它的自动化流程快而且Juicebox的人工检查环节能让你心里特别有底。你可以先快速跑出一个草案然后用Juicebox细细调整。如果你的数据质量一般或者你是新手我推荐从3D-DNAJuicebox入手。因为可视化纠错的能力能帮你更直观地理解Hi-C数据的含义以及组装错误的表现形式。这个过程本身就是最好的学习。当你不确定时或者数据非常珍贵可以两条路都试试。用同样的预处理数据分别跑一遍ALLHiC不运行prune和3D-DNA对比它们自动组装的结果。看看在染色体分群、contig排序上是否一致。如果核心部分一致那结果就非常可靠如果有较大分歧就需要你深入查看Hi-C热图判断哪一种更符合生物学规律。最后无论选择哪条路都要记住Hi-C挂载不是魔法。它高度依赖于预处理数据的质量比对率、有效互作对比例以及你提供的contig的连续性N50值。在开始挂载前务必确保你的预处理结果是可靠的否则后续步骤都是空中楼阁。好了工具和路都指给你了接下来就打开终端用你的数据实战吧。遇到具体报错别慌多看看日志大部分问题在GitHub的Issues里都能找到答案。