修改vcf文件中的染色体名

问题描述

在生物信息分析过程中,常常需要修改 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
  • 7
    点赞
  • 21
    收藏
    觉得还不错? 一键收藏
  • 5
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 5
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值