ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

亲手跑通 Snippy 变异检测:三步把测序 reads 变成一张 SNP 变异表

亲手跑通 Snippy 变异检测:三步把测序 reads 变成一张 SNP 变异表 亲手跑通 Snippy 变异检测三步把测序 reads 变成一张 SNP 变异表【免费下载链接】snippy:scissors: :zap: Rapid haploid variant calling and core genome alignment项目地址: https://gitcode.com/gh_mirrors/sn/snippySnippy 是一款面向细菌、病毒等单倍体基因组的快速变异检测工具它能把测序 reads 与参考基因组逐位比对一次性输出 SNP、插入缺失和基因影响注释省去你手动串联七八个生信工具的时间。这篇文章不绕弯子先让你看清装完之后能拿到什么再带你亲手跑通第一个样本全程每一步都有可以对照的预期结果。先看一个真实卡壳现场手工变异检测到底卡在哪隔壁实验室的阿哲上周接了个活儿一批大肠杆菌的测序数据要找耐药相关突变。他的工具链其实挺全——bwa 比对、samtools 过滤、freebayes 找变异、snpEff 注释每个都单独装过。但真正跑起来才知道有多折磨每换一个样本都要重敲一大串管道命令参数记不住翻旧脚本抄还抄错比对结果、变异结果、注释结果散落在四五个目录里想回看某个位点得来回切同事分享的脚本用的是旧版语法跑一半就崩报错信息他看不太懂最后拿到手的 .vcf 文件注释列密密麻麻他只知道有个 missense_variant却不知道到底影响了哪个基因。他凌晨一点还在群里问就没有一个工具能把比对→找变异→注释打包成一条命令把所有结果整齐收进一个文件夹吗有就是Snippy。它最初为细菌基因组设计实际只要是单倍体每个位置只有一份拷贝的样本都能处理——细菌、病毒、质粒、线粒体都可以。接下来咱们就把它请进家门。装完之后你会多出一个变异车间先看产物长什么样先别急着敲安装命令。咱们先把目标画面立起来等 Snippy 装好、跑完一个样本之后你的输出目录里会躺着这么一套文件文件用途snps.tab人类最好读的变异汇总表一列一个字段snps.vcf标准 VCF 格式带 snpEff 注释可丢给下游工具snps.bam所有 reads 的比对结果可放进 IGV 里逐个位点肉眼核对snps.consensus.fa把变异贴回参考基因组得到的修正版基因组snps.raw.vcf / snps.filt.vcffreebayes 的原始调用与过滤后结果snps.log每一步跑了什么命令、输出了什么全程留痕想象你拥有了一座变异检测车间投进去两样原料——测序 reads 和一份参考基因组出来的是一件成品——一张写着哪个位点、什么类型、落在哪个基因、会不会改变氨基酸的变异表。单个样本能跑十几个样本批量也能跑最后还能自动算出一份所有样本共有变异位点的核心比对文件直接喂给建树软件画系统发育树。把目标记在心里接下来就按开车间→点人名→接订单三步走。车间运转原理为什么 Snippy 自己装好还不够先花两分钟看懂这个车间怎么转后面排错会省很多时间。把参考基因组想象成一张标准图纸你的测序 reads 是这台机器实际零件的照片碎片。Snippy 的车间里有五个工位备料读取 FASTQ必要时按比例抽样降低深度对应--subsample参数比对bwa 把每条 read 贴到图纸上samtools 负责排序、去重、过滤找差异freebayes 逐位对比找出样本与图纸不同的位点区分 snp、插入、缺失注释snpEff 对照 GenBank 注释告诉你变异落在哪个基因、会不会改变氨基酸打包bcftools 等把结果整理成 vcf、tab、consensus 一整套产物。Snippy 的角色是车间主任你只对主任下单一条命令排程、调用工位、验收成品都是它的事。这也解释了一个新手最容易踩的坑——把 Snippy 单独装好远远不够bwa、samtools、freebayes、snpEff 这些工人缺了任何一个订单跑到一半就会卡住。所以接下来的安装环节核心关注点只有一个依赖能不能一并解决。把车间开起来两条通路把 Snippy 本体请进家门装 Snippy 本体其实只有两种主流方式选哪种取决于你手里的环境。通路 A源码克隆追最新版。在你想放代码的目录里执行git clone https://gitcode.com/gh_mirrors/sn/snippy.git export PATH$PWD/snippy/bin:$PATH预期结果敲snippy --help能弹出完整的参数说明而不是 command not found。通路 Bconda 一键装齐最省心。只要机器上已经配好 bioconda 频道一条命令把 Snippy 和全部依赖一次性装上conda install -c conda-forge -c bioconda -c defaults snippy预期结果conda 打印一串待安装包列表里面除了 snippy还能看到 bwa、samtools、freebayes、snpEff 等工人的名字装完自动回到命令行提示符。两相对比源码通路永远拿到最新代码、方便二次开发但依赖要自己逐个补齐参考仓库根目录的environment.yml它把全部依赖列得清清楚楚conda 通路版本可能略旧但依赖自动解决对新手友好得多。嫌麻烦的直接选 B。开工前先点名用两条命令确认八个工位都在车间开了门别急着接真订单先点个名。第一声确认主任在不在。snippy --version预期结果打印类似snippy 5.0.0-dev的版本串这个版本号定义在仓库的perl5/Snippy/Version.pm里想确认版本信息可以直接翻它。如果这里报错说明 PATH 没配对回上一步检查。第二声确认每个工位都有人。snippy --check预期结果工具会逐个探测 bwa、minimap2、samtools、bcftools、bedtools、freebayes、snpEff、samclip 等组件可用的显示 OK 或版本号缺失的会明确标出来。哪一项标了缺失就单独补哪一项conda install -c bioconda samtools bcftools bwa freebayes snpeff samclip seqtk补完再跑一遍snippy --check直到全部通过。这一步值得认真对待——点名全勤再开工能避免真数据跑到一半才发现缺工具。第一张订单用仓库自带测试材料跑出 snps.tab仓库的test目录里有一套现成的测试材料example.gbk带注释的参考基因组、example.fna参考序列、example.bed掩蔽区域文件。咱们的第一张订单就用它们。真实测序 reads 文件太大这里先用 wgsim 按参考序列模拟一对双端 reads——仓库里test/Makefile的官方测试流程也是这么干的wgsim -S 1 -h -r 0.005 -N 12000 -1 100 -2 100 -d 200 example.fna R1.fq R2.fq预期结果生成R1.fq、R2.fq两个文件各约 12000 条 100 bp 的 reads其中约 0.5% 的位点被随机改成了变异。然后正式下单snippy --cpus 4 --outdir my_first_run --ref example.fna --R1 R1.fq --R2 R2.fq预期结果日志里依次出现 bwa 比对、freebayes 变异识别等字样最后打印Walltime used: 3 min, 42 sec Results folder: my_first_run Done.看成品head -5 my_first_run/snps.tab每一列的含义是CHROM变异所在的序列名POS位置TYPE变异类型snp 单碱基替换 / mnp 多碱基替换 / ins 插入 / del 缺失 / complex 复合REF参考碱基ALT样本碱基EVIDENCE支持各碱基的 reads 计数。如果你把参考换成example.gbk这种带注释的文件表格还会多出基因名、产物描述和影响效应列——直接告诉你变异落在哪个基因、会不会改变氨基酸这正是 Snippy 最贴心的设计之一。到这里第一张订单已经交付了。批量交付与三种省力打法单样本跑通后真实场景通常是十几个样本对同一个参考。一个个手敲命令太低效Snippy 准备了批量入口snippy-multi先写一个制表符分隔的清单文件每行一个样本Isolate1 /path/to/R1.fq.gz /path/to/R2.fq.gz Isolate2 /path/to/SE.fq.gz Isolate3 /path/to/contigs.fa然后生成并检查批量脚本snippy-multi input.tab --ref Reference.gbk --cpus 16 runme.sh less runme.sh sh runme.sh预期结果每个样本输出一个结果文件夹最后自动调用snippy-core打印类似Found 2814 core SNPs from 96615 SNPs.的汇总——意思是总共有 96615 个变异位点其中 2814 个是全体样本共有的核心 SNP。产物core.aln可以直接丢给 FastTree 等工具画树。跑批量的同时还有三个高频场景值得记住深度太高跑得慢上千倍深度时绝大多数变异在 50~100 倍就够检出了用--subsample 0.1按比例抽读日志里会出现 Sub-sampling reads at rate 0.1速度立竿见影只想看特定区域比如只关心耐药基因把区间写进 BED 文件用--targets sites.bed限定范围能省大量计算时间只有 contigs 没有原始 reads改用--ctgs contigs.faSnippy 会把 contigs 撕成 250 bp 的伪 reads 再比对日志会提示 Shredding ... into pseudo-reads输出目录与 reads 样本完全兼容可以混在一起参与批量分析。另外如果你的参考基因组里有重复区域容易产生假阳性比如结核分枝杆菌的 PE/PPE 家族可以用--mask传入掩蔽文件仓库的etc/Mtb_NC_000962.3_mask.bed就自带一份现成的结核菌掩蔽区域。验收清单达到这五条才算真正装好对照检查一下别急着收工snippy --version能打印出版本号snippy --check全部组件 OK无一缺失测试样本跑通my_first_run里出现 snps.vcf、snps.tab、snps.bam、snps.consensus.fa能看懂 snps.tab 每一列的含义知道 snp / ins / del 的区别进阶snippy-multi生成了 runme.sh批量结束后得到了 core.aln。装好之后接下来做这三件事用测试数据把单样本→批量完整走一遍把每一步的预期输出记在脑子里再碰真实数据真实样本跑之前备份好原始 reads并给样本取清晰规范的 ID——样本 ID 会写进 BAM 和 VCF 的 Read Group命名混乱以后追溯起来很痛苦记下版本号。变异检测结果跟工具版本强相关写文章或汇报时注明snippy --version的输出你的结果才可复现、可信。等这些基本功都扎实了还有进阶玩法等着你用--ctgs做 contig 纠错用--unmapped抓出没比上的 reads那里往往藏着质粒等新元件用snippy-core产出核心基因组比对后接上建树流程。Snippy 的价值就是把比对—变异识别—注释—输出这条长流水线收进一条命令当你熟练地把第一份 reads 变成一张清晰的变异表时后面的群体分析和系统发育研究就有了可靠的数据地基。【免费下载链接】snippy:scissors: :zap: Rapid haploid variant calling and core genome alignment项目地址: https://gitcode.com/gh_mirrors/sn/snippy创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考
返回列表