This is a test version of Biostars. For the public version, visit https://www.biostars.org.
VCF file REF column only contains N

I am using samtools mpileup to call variant.

Command used:

samtools mpileup -u -f "/gpfs/home/g/o/Homo_sapiens.GRCh37.75.dna.primary_assembly.fa" "/gopisiva/Mason/Single_cell/SC_Stressed/*.bam" | bcftools view -bvcg - > SC_Stressed.raw.bcf"

bcftools view SC_Stressed.raw.bcf -v vcfutils.pl varFilter -D100 > SC_Stressed.flt100.vcf

OUTPUT:

INFO= <ID=PR,Number=1,Type=Integer,Description="# permutations yielding a smaller PCHI2.">
CHROM   POS ID  REF ALT QUAL    FILTER  INFO    FORMAT SRR1481944_Aligned.sorted.bam

    chr1    16961   .   N   T   39.8    .   DP=2;VDB=5.960000e-02;AF1=1;AC1=2;DP4=0,0,2,0;MQ=60;FQ=-33;AF1=1;AC1=2;FQ=-33   GT:PL:GQ    1/1:71,6,0:10
    chr1    16962   .   N   C   39.8    .   DP=2;VDB=5.960000e-02;AF1=1;AC1=2;DP4=0,0,2,0;MQ=60;FQ=-33;AF1=1;AC1=2;FQ=-33   GT:PL:GQ    1/1:71,6,0:10
    chr1    16963   .   N   G   39.8    .   DP=2;VDB=5.960000e-02;AF1=1;AC1=2;DP4=0,0,2,0;MQ=60;FQ=-33;AF1=1;AC1=2;FQ=-33   GT:PL:GQ    1/1:71,6,0:10
    chr1    16964   .   N   T   39.8    .   DP=2;VDB=5.960000e-02;AF1=1;AC1=2;DP4=0,0,2,0;MQ=60;FQ=-33;AF1=1;AC1=2;FQ=-33   GT:PL:GQ    1/1:71,6,0:10
vcf ref mpileup samtools
• 3,704 views • link
updated by John Torcivia • John Torcivia • written by genie66 • genie66

1 answer

If you run your bcftools line like:

bcftools view SC_Stressed.raw.bcf -v | vcfutils.pl varFilter -D100 > SC_Stressed.flt100.vcf

Do you get a better result? It looks like the pipe to vcfutils.pl might have been missed. Otherwise everything looks fine to me so its strange it is reporting things with a read depth of 2 and not a minimum of 100.

• 1 views • link
written by John Torcivia • John Torcivia

When I add pipe, I am getting nothing after the header description. Do you find any mistake in the command?

• 1 views • link
written by genie66 • genie66

I think everything looks good. What is the size of "SC_Stressed.raw.bcf" after you do the mpileup stage but before the bcftools?

• 0 views • link
written by John Torcivia • John Torcivia

Thanks @jackfrost this following command worked for me. Hope it help someone

samtools mpileup -uf hg19.fa SRR14819$i.sorted.bam | bcftools call -mv > SRR14819$i.Control.vcf
• 0 views • link
written by genie66 • genie66

Log in to answer this question.

More posts like this