genotype calling using samtools

Question

[Karim]

Dear All,

I am using samtools/bcftools (version 1.3) to call SNP genotypes from RNA-Seq data. I added some input and output options as well as some filter options to command line as below;

samtools mpileup -uf gg.fa -r 27:4875164-4886606 a.bam b.bam | bcftools view -v snps -g hom --no-header - > out.bcf

It seems that the -g option does not work appropriately because when replacing "hom" with "het" or removing it from the command line, nothing changes in the output file. My question is: Do I need some extra options to get the -g option working?

Any help appreciated

[Chris F.]

You're not actually calling variants here, just doing the pileup. You probably want to use the -g or -v options with mpileup and pipe it to bcftools call to create vcf file. Once that's generated, you can filter it through view with your criteria as above.