Hi there,
I am trying to get ‘samtools mpileup’ to output info fields like AD and DP from
bamfile for a list of genomic positions that I provide with the –l option. I
also would like mpileup to output a position even if with zero depth.
When I use the command below, it works fine:
samtools mpileup -OsB --min-BQ 0 -a -f ${genome_fasta} -r chr4 -l ${loci_red}
$bam
[mpileup] 1 samples in 1 input files
<mpileup> Set max per-file depth to 8000
chr4 106158247 T 0 * * * *
However, I would like the output in VCF format so that I can specify the output
tags with --output-tags. When I try to run:
samtools mpileup -OuvsB --min-BQ 0 -a --output-tags AD,DP -f ${genome_fasta} -r
chr1 -l ${loci_red} $bam
##fileformat=VCFv4.2
##FILTER=<ID=PASS,Description="All filters passed">
##samtoolsVersion=1.6+htslib-1.6
##samtoolsCommand=samtools mpileup -OuvsB --min-BQ 0 -a --output-tags AD,DP -f
/home/users/allstaff/quaglieri.a/PHD_project/GEO_Leucegene_data/genomes/hg19/Homo_sapiens/UCSC/hg19/Sequence/WholeGenomeFasta/genome.fa
-r chr1 -l
/wehisan/general/user_managed/grpu_majewski_3/quaglieri.a/GEO_Leucegene_data/scripts/02-variant_loss_downsampling/mutations_loci_red.bed
G6.6248889.Cyc1.Resp_Recal.reordered_DeDupl.rg.split.bam
##reference=file:///home/users/allstaff/quaglieri.a/PHD_project/GEO_Leucegene_data/genomes/hg19/Homo_sapiens/UCSC/hg19/Sequence/WholeGenomeFasta/genome.fa
##contig=<ID=chr1,length=248956422>
##contig=<ID=chr2,length=242193529>
##contig=<ID=chr3,length=198295559>
##contig=<ID=chr4,length=190214555>
##contig=<ID=chr5,length=181538259>
##contig=<ID=chr6,length=170805979>
##contig=<ID=chr7,length=159345973>
##contig=<ID=chr8,length=145138636>
##contig=<ID=chr9,length=138394717>
(…..)
##FORMAT=<ID=AD,Number=R,Type=Integer,Description="Allelic depths">
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT G6
<mpileup> Set max per-file depth to 8000
I get no locations as output. It looks like the option –a does not work with –v
–r –g output.
Here is my samtools version
Samtools --version
samtools 1.6
Using htslib 1.6
Copyright (C) 2017 Genome Research Ltd.
Do you have any suggestion of why this happens?
Thanks for your support.
All the best,
Anna
--
Anna Quaglieri
PhD student
Division of Bioinformatics - Speed lab
Walter & Eliza Hall Institute of Medical Research
1G Royal Parade, Parkville VIC 3052, Australia
Contact: +614 6892 5003
---I'm an Athena SWAN Advocate. I call out sexism and promote equality.---
------------------------------------------------------------------------------
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
