Entering edit mode
9.7 years ago
mbourgey
•
0
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
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.
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
I don't believe it is an issue of genotype probabilites. What value did you feed to the option
-m
?Hi Cytosine,
Thanks for your reply.
I didn't feed any value for the
-m
option in bcftools call as it does not require anyone:and let the default value (not set) in
samtools mpileup
for the-m
option:Mathieu
Strange... I have bcftools v0.1.19-44428cd and the option menu says this:
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:Could be that these two are interfering in your case, since
-v
forces-c
, which conflicts with your-m
.Just for your information: it's the commands and options I used with samtools 1.0
I will check for the
-c
vs.-m
conflictThanks
Mathieu
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