Dear James,
I have installed the latest github htslib, samtools and bcftools
Program: bcftools (Tools for variant calling and manipulating VCFs and BCFs)
Version: 1.5-13-gab7b4c1 (using htslib 1.5-5-gda5c0c7)
and still I get not the full number of sites.
However, it seems that these are all the positions that have ZERO
coverage throughout all BAM files. However, it was stated that with the
-a -a option one would also retain these.
So it is likely not a BUG but than it would be a request how to also
retain these sites as DP=0 in the mpileup file?
Unfortunately the INDELs end up as seperate lines, which can be
partially fixed with bcftools norm:
bcftools mpileup -r chr1:50000001-50025000 --fasta-ref REF.fasta -a
AD,DP -Q 0 -A -B -O z -b BAM.list | bcftools norm -O v -m +any | awk
'{if(substr($1,1,1)!="#") print $2}' | sort | uniq | wc -l
[mpileup] 6 samples in 6 input files
Lines total/split/realigned/skipped: 25140/0/0/0
24923
I have also checked if there might be something wrong in that reference
location, but this region consists of:
A C G T
7014 4852 4668 8466
So no gaps or 'N' chars in it.
Thank you in anticipation
Kristian
On 07/12/2017 06:08 PM, James Bonfield wrote:
On Wed, Jul 12, 2017 at 05:54:06PM +0200, Kristian Ullrich wrote:
HOWEVER, if I add the options '-u -v' so that the output should be
uncompressed and in VCF format, than the output has duplicated
POSITIONs and in total less than the expected 25000 sites, if I
remove the duplicated. So there might be an issue in reporting the
VCF output or in translating the compressed to uncompressed output.
It's unlikely to be a compression issue. Is this the latest samtools
and htslib release?
I don't see how specifying a single region could ever give duplicated
positions (and some missing) in VCF though, hence this is likely a
bug!
It is worth testing bcftools mpileup, just to be sure it's not been
fixed. The samtools mpileup -v code has migrated to bcftools now and
it's deprecated in samtools. There are some minor command line
differences and maybe other tweaks, so it's worth checking it there.
If it still fails with the latest bcftools code then can you please
file a ticket on github?
https://github.com/samtools/bcftools
Apologies for the problems.
James
--
Kristian Ullrich, Ph.D.
Department of Evolutionary Genetics
Max Planck Institute for Evolutionary Biology
August Thienemann Str. 2
24306 Ploen
Germany
+49 4522 763 274
ullr...@evolbio.mpg.de
------------------------------------------------------------------------------
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
_______________________________________________
Samtools-help mailing list
Samtools-help@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/samtools-help