Question: samtools/bcftools strange (error ?) lack of call
0
gravatar for mbourgey
5.1 years ago by
mbourgey0
Canada
mbourgey0 wrote:

Hi all,

I got a strange lack of call when I used samtolls mpileup + bcftools to call my SNP (version 0.19 or 1.0 give me the same  outputs)

here as what I observe when I use samtools+bcftools view:

samtools mpileup -q 1 -u -D -S -g -f hg1k_v37.fasta -r 3:44775917-44775917 alignment/sample1.sorted.dup.recal.bam alignment/sample2.sorted.dup.recal.bam alignment/sample3.sorted.dup.recal.bam | bcftools view

##ALT=<ID=X,Description="Represents allele(s) other than observed.">
##INFO=<ID=INDEL,Number=0,Type=Flag,Description="Indicates that the variant is an INDEL.">
##INFO=<ID=IDV,Number=1,Type=Integer,Description="Maximum number of reads supporting an indel">
##INFO=<ID=IMF,Number=1,Type=Float,Description="Maximum fraction of reads supporting an indel">
##INFO=<ID=DP,Number=1,Type=Integer,Description="Raw read depth">
##INFO=<ID=VDB,Number=1,Type=Float,Description="Variant Distance Bias for filtering splice-site artefacts in RNA-seq data (bigger is better)",Version=3>
##INFO=<ID=RPB,Number=1,Type=Float,Description="Mann-Whitney U test of Read Position Bias (bigger is better)">
##INFO=<ID=MQB,Number=1,Type=Float,Description="Mann-Whitney U test of Mapping Quality Bias (bigger is better)">
##INFO=<ID=BQB,Number=1,Type=Float,Description="Mann-Whitney U test of Base Quality Bias (bigger is better)">
##INFO=<ID=MQSB,Number=1,Type=Float,Description="Mann-Whitney U test of Mapping Quality vs Strand Bias (bigger is better)">
##INFO=<ID=SGB,Number=1,Type=Float,Description="Segregation based metric.">
##INFO=<ID=MQ0F,Number=1,Type=Float,Description="Fraction of MQ0 reads (smaller is better)">
##INFO=<ID=I16,Number=16,Type=Float,Description="Auxiliary tag used for calling, see description of bcf_callret1_t in bam2bcf.h">
##INFO=<ID=QS,Number=R,Type=Float,Description="Auxiliary tag used for calling">
##FORMAT=<ID=PL,Number=G,Type=Integer,Description="List of Phred-scaled genotype likelihoods">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Number of high-quality bases">
##FORMAT=<ID=SP,Number=1,Type=Integer,Description="Phred-scaled strand bias P-value">
##bcftools_viewVersion=1.0+htslib-1.0
##bcftools_viewCommand=view
#CHROM    POS    ID    REF    ALT    QUAL    FILTER    INFO    FORMAT    sample1    sample2    sample3
3    44775917    .    A    C,<X>    0    .    DP=247;I16=207,1,36,0,6174,186100,1117,35331,12480,748800,2160,129600,3824,84932,643,14149;QS=2.65514,0.344865,0;VDB=0.986901;SGB=0.623532;RPB=0.19008;MQB=1;MQSB=1;BQB=2.58832e-10;MQ0F=0    PL:DP:SP    0,90,152,90,152,152:30:0    20,0,115,255,185,236:112:0    0,42,110,255,149,215:102:0

when I then use: bcftools call -m -v I didn't get any variant called whereas the PL field of the sample 2 (20,0,115,255,185,236) show an higher likelihood for the heterozygote call.

Any idea what's happens here ?

thanks in advance

Mathieu

snp bcftools ngs samtools • 2.6k views
ADD COMMENTlink modified 5.1 years ago • written 5.1 years ago by mbourgey0

Hello mbourgey!

It appears that your post has been cross-posted to another site: samtools-help list

This is typically not recommended as it runs the risk of annoying people in both communities.

ADD REPLYlink written 5.1 years ago by Devon Ryan92k
0
gravatar for mbourgey
5.1 years ago by
mbourgey0
Canada
mbourgey0 wrote:

yes it 's me

sorry about that !

I posted here then think perhaps it should be more appropriate to ask that on the mailing list.

Mathieu

ADD COMMENTlink written 5.1 years ago by mbourgey0
0
gravatar for Cytosine
5.1 years ago by
Cytosine440
Ljubljana, Slovenia
Cytosine440 wrote:

I don't believe it is an issue of genotype probabilites. What value did you feed to the option "-m"?

ADD COMMENTlink written 5.1 years ago by Cytosine440
0
gravatar for mbourgey
5.1 years ago by
mbourgey0
Canada
mbourgey0 wrote:

Hi Cytosine,

Thanks for your reply.

 

I din't feed any value for the -m option in bcftools call as it does not require anyone:

-m, --multiallelic-caller       alternative model for multiallelic and rare-variant calling (conflicts with -c)

 

and let the default value (not set) in samtools mpileup for the -m option:

 -m, --min-ireads INT    minimum number gapped reads for indel candidates [1]

 

Mathieu

 

ADD COMMENTlink written 5.1 years ago by mbourgey0
0
gravatar for Cytosine
5.1 years ago by
Cytosine440
Ljubljana, Slovenia
Cytosine440 wrote:

Strange... I have bcftools v0.1.19-44428cd and the option menu says this:

 -m FLOAT  alternative model for multiallelic and rare-variant calling, include if P(chi^2)>=FLOAT

Meaning filter according to the chi^2  statistic... I thought your output was empty since you did not provide a value to -m

 

On the other hand, your option does tell something interesting (conflicts with -c)

And you call bcftools view also with -v, which makes it imply -c as well:

-v        output potential variant sites only (force -c)

 

Could be that these two are interfering in your case, since -v forces -c, which conflicts with your -m.

ADD COMMENTlink written 5.1 years ago by Cytosine440
0
gravatar for mbourgey
5.1 years ago by
mbourgey0
Canada
mbourgey0 wrote:

just for information: it's the commands and options I used with samtools 1.0

I will check for the -c vs. -m conflict

Thanks

Mathieu

ADD COMMENTlink written 5.1 years ago by mbourgey0
0
gravatar for mbourgey
5.1 years ago by
mbourgey0
Canada
mbourgey0 wrote:

in bcftools 1.0 the -v option does not imply/force -c.

but just to be sure I tried using bcftools call -c -v instead and the result is the same no variant called

ADD COMMENTlink written 5.1 years ago by mbourgey0
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: 2018 users visited in the last hour