问题描述
在生物信息分析过程中，常常需要修改 vcf 染色体名称，即删除染色体名中的chr字符或添加chr字符到vcf文件中。
这里使用shell脚本和bcftools两种方式实现。大致比较了一下，shell脚本在速度上稍微快一些。
1. 通过 shell 脚本实现
下面的shell脚本要求输入vcf为解压缩（因为我的数据是没有压缩的…），如果是压缩文档，自行修改一下代码即可。
通过 awk 添加 chr 字符
awk '{
if($0 !~ /^#/)
print "chr"$0;
else if(match($0,/(##contig=<ID=)(.*)/,m))
print m[1]"chr"m[2];
else print $0
}' no_chr.vcf > with_chr.vcf
通过 sed 删除 chr 字符
sed 's/##contig=<ID=chr/##contig=<ID=/g' with_chr.vcf | sed 's/^chr//g' > no_chr.vcf
代码封装
因为需要经常用到上面的两段代码，所以将它们封装到一个名为change_vcf_name.sh的脚本中。
#! /user/bin/bash
echo "usage: change_vcf_name.sh a|r input_vcf_file output_vcf_file"
echo "a (add) : 添加chr字符到vcf文件中"
echo "r (remove): 删除vcf文件中的chr字符"
if [ $1 = "a" ]
then
awk '{
if($0 !~ /^#/)
print "chr"$0;
else if(match($0,/(##contig=<ID=)(.*)/,m))
print m[1]"chr"m[2];
else print $0
}' $2 > $3
elif [ $1 = "r" ]
then
sed 's/##contig=<ID=chr/##contig=<ID=/g' $2 | sed 's/^chr//g' > $3
else
echo "输入的第一个参数不合适，请输入a或r。"
fi
代码调用：
# 删除vcf中的chr字符
bash change_vcf_name.sh r input.vcf output.vcf
# 加入chr字符到vcf中
bash change_vcf_name.sh a input.vcf output.vcf
2. 通过 bcftools 程序实现
需先准备一个bcftools的配置文件 （这里命名为chr_name_change.txt，内容如下，由tab分隔成两列。
bcftools会根据这个配置文件修改vcf的名称（将左边的名称替换为右边，例如这里会将chr1替换为1、chr2替换为2，…）
chr1 1
chr2 2
chr3 3
chr4 4
chr5 5
chr6 6
chr7 7
chr8 8
chr9 9
chr10 10
chr11 11
chr12 12
chr13 13
chr14 14
chr15 15
chr16 16
chr17 17
chr18 18
chr19 19
chr20 20
chr21 21
chr22 22
chrX X
chrY Y
chrMT M
程序运行代码
- 调用bcftools的annotate子程序
- 调用annotate子程序的rename-chrs模块
- 调用bgzip进行压缩
bcftools annotate --rename-chrs chr_name_change.txt cohart.filter.vcf.gz | bgzip -c > cohart.filter.nochr.vcf.gz
版权声明：本文为liangbilin原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/liangbilin/article/details/108991510
网友评论