使用bcftools call 感兴趣的基因

最新推荐文章于 2024-05-10 09:38:51 发布

爱笑的小牙

最新推荐文章于 2024-05-10 09:38:51 发布

阅读量2k

点赞数

分类专栏：生物信息分析

本文链接：https://blog.csdn.net/Cassiel60/article/details/88998254

版权

看很多文章都是用samtools mpileup 来提取感兴趣的位置，其实现在使用这个命令时，会提示我们这个命令已经过时了，使用bcftools mpileup 和bcftools call 。

bcftools是附属于samtools的程序，大多数用法是相同的，只是一些参数的变化，可以用 bcftools mpileup来查看具体用法，或者看文档。https://samtools.github.io/bcftools/howtos/variant-calling.html

第一次使用时，-r 后面是1：226249552-226259702 ，总会报错，提示：

[E::faidx_adjust_position] The sequence not found; 
#数据库版本弄错

根据samtools view WH.sorted.markdup.BQSR.bam | head -n 10 时发现里面的染色体是chr，所以就加上了chr,如下

bcftools mpileup -r chr1:226249552-226259702 -Ou -f $reference  WH.sorted.markdup.BQSR.bam  | bcftools call -mv -Oz --ploidy 1 -o H3F3A.vcf.gz

此时运行ok

批量运

关注