Question: Remove not polymorphic sites from VCF
0
gravatar for GabrielMontenegro
21 months ago by
United Kingdom
GabrielMontenegro560 wrote:

I have a VCF file and I want to remove sites that are not polymorphic sites i.e. sites where all individuals are REF/REF or ALT/ALT. Is there a quick way to do this using bcftools or vcftools?

Thanks.

genome vcf • 1.9k views
ADD COMMENTlink modified 12 months ago by amy.roberts20 • written 21 months ago by GabrielMontenegro560
2
gravatar for finswimmer
21 months ago by
finswimmer13k
Germany
finswimmer13k wrote:

Now I hopefully understand your goal. This will remove side where all samples are either hom ref or hom alt:

$ bcftools view -e 'COUNT(GT="AA")=N_SAMPLES || COUNT(GT="RR")=N_SAMPLES' input.vcf

Have a look at the manual to explore more ways to filter vcf files using bcftools.

I've deleted my old answer, because it hasn't answer your question.

ADD COMMENTlink written 21 months ago by finswimmer13k

Thanks! Yes, I think that answers my question. I got this error message though: the tag "INFO/COUNT" is not defined in the VCF header. Do you know how to add this info on my vcf?

ADD REPLYlink written 21 months ago by GabrielMontenegro560
1

I think you're using an outdated version of bcftools. The method COUNT was introduced in v1.7. The current version is v1.9.

Which version are you using? Please upgrade.

ADD REPLYlink written 21 months ago by finswimmer13k

Using v1.9 solved it - Thanks!

ADD REPLYlink written 21 months ago by GabrielMontenegro560

@finswimmer have a quick question, will bcftools get the COUNT of the AA or RR genotypes based on the INFO field? Or will it compute the COUNT for all sites? I am applied this filter after merging my dataset with another dataset and I think I am getting weird results.

ADD REPLYlink modified 20 months ago • written 20 months ago by GabrielMontenegro560

For each line in the vcf file bcftools will have a look at the genotype of each sample. In the INFO column there is no information about the counts.

Why do you think you are getting a weird result? As always, a small example dataset which shows your problem would be helpful.

ADD REPLYlink written 20 months ago by finswimmer13k
2
gravatar for Pierre Lindenbaum
21 months ago by
France/Nantes/Institut du Thorax - INSERM UMR1087
Pierre Lindenbaum131k wrote:

using vcffilterjdk: http://lindenb.github.io/jvarkit/VcfFilterJdk.html

 java -jar dist/vcffilterjdk.jar -e 'return !(variant.getGenotypes().stream().allMatch(G->G.isHomVar()) || variant.getGenotypes().stream().allMatch(G->G.isHomRef()));' input.vcf
ADD COMMENTlink written 21 months ago by Pierre Lindenbaum131k
2
gravatar for amy.roberts
12 months ago by
amy.roberts20
amy.roberts20 wrote:

You can also use the -c function in view

-c, --min-ac INT[:nref|:alt1|:minor|:major|:'nonmajor']
minimum allele count (INFO/AC) of sites to be printed. Specifying the type of allele is optional and can be set to non-reference (nref, the default), 1st alternate (alt1), the least frequent (minor), the most frequent (major) or sum of all but the most frequent (nonmajor) alleles.

By setting -c 1, you will only keep sites with at least one nonref allele (=polymorphic in your data)

bcftools view -c 1 input.vcf.gz -o output.vcf.gz -Oz
ADD COMMENTlink modified 12 months ago by RamRS30k • written 12 months ago by amy.roberts20
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 2315 users visited in the last hour