vcf文件合并

1、vcf文件合并

(1)同样本vcf合并,这里的相同样本合并指的是不同的染色体,或者染色体区间。

使用方法 perl merge.pl list merge.vcf.gz

这里的list就是不同的染色体的vcf文件,如:

chr1:0-2000.vcf.gz
chr1:2001-4000.vcf.gz
chr3.vcf.gz
……

merge.pl

 

use strict;
use warnings;
open IN,"< $ARGV[0]" or die $!;
open OUT,'>:gzip',"$ARGV[1]" or die $!;
my $n=0;
while(<IN>){
	chomp;
	my $line=$_;
	if($line=~m/\.gz/){
	open VCF,"gunzip -dc $line | " or die $!;
	}else{
	open VCF,"< $line" or die $!;
	}
	$n++;
	if($n==1){
		while(<VCF>){
			chomp;
				print OUT "$_\n";
			}
		}
	if($n>1){
		while(<VCF>){
			chomp;
			my @aa=split /\t/;
			next if($aa[0]=~/^#/);
			print OUT "$_\n";
			}
		}
	}
close VCF;
close IN;
close OUT;

这里也可以使用zcat或者cat的模式,但是我觉得这样太麻烦了。需要将head行进行处理。

 

(2)不同样本vcf合并

主要是的软件就是bcftools和gatk。

bcftools参考网址:http://vcftools.sourceforge.net/htslib.html#merge

一行简单的代码:

bcftools merge A.vcf.gz B.vcf.gz C.vcf.gz -Oz -o ABC.vcf.gz

 

这个软件使用起来较为麻烦,需要将原始的vcf转换为bgzip的压缩,之后再建立index,之后才能用于合并。

这里简单举一个例子,方便理解,如下

有两个文件,分别执行下面的命令:

gunzip file.vcf.gz 
bgzip file.vcf
tabix -p vcf file.vcf.gz 

这样两个文件的tbi索引就建立好了,直接执行:

bcftools merge file1.vcf.gz file2.vcf.gz -Oz -o merge.vcf.gz

容易出现的报错:vcf 没有进行排序就行进行index。必须是进行排序后的vcf才能进行。

除了使用bcftools,还可以使用gatk,但是gatk必须是4.0以下的版本才行。代码如下:

 

gatk -T CombineVariants -V file1.vcf.gz -V file2.vcf.gz -o merge.vcf.gz -R ref.fa

总结:这样相同样本合并和不同样本合并就都可以进行了,推荐使用GATK进行合并。

 

 

 

 

评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值