Alternative Depth (AD) field in VCF created with samtools mpileup and bcftools call

Question

[dutoit.ludovic]

Hi,

we created a VCF file for multiple samples "using samtools mpileup and bcftools call". Those VCF files only have GT and PL field while we need DP and AD for downstream analyses.

We need Alternative depth AD and Depth DP as a field in our VCF. It does not seem possible to get it out of bcftools or samtools.

using bcfview on any single sample bcf, one can obtain DP or DP4 from the info field. The only way I can think of extracting AD is by summing the fields of DP4 2 by 2 and them putting them back in the vcf using my own custom scripting.

Is that methodologically correct or am I missing the point here and that would not really be AD? Does anyone has a better way of doing it or any suggestions?

Many many thanks for your help

Ludo

[dutoit.ludovic]

That is extremely helpful, thank you both for your time!

Ludo

Answer 1

[finswimmer]

Hello dutoit.ludovic,

using samtools mpileup with the aim to call variants is deprecated. You should use bcftools mpileup instead. Here you have the -a parameter to add additional annotations like AD (which means Allelic Depth by the way :) )

Have a look at the "output options" for mpileup in the docs for the fields available.

fin swimmer

Answer 2

[Kevin Blighe]

You can add tags, including DP, during bcftools mpileup via the --annotate parameter - take a look HERE in the BCFtools manual page. I also show an example of using this, here:

"${BCFtools_root}"bcftools mpileup \
  --redo-BAQ \
  --min-BQ 30 \
  --per-sample-mF \
  --annotate FORMAT/AD,FORMAT/ADF,FORMAT/ADR,FORMAT/DP,FORMAT/SP,INFO/AD,INFO/ADF,INFO/ADR \
  -f "${Ref_FASTA}" \
  "${Results_root}"/"${RunNumber}"/"${PatientID}"/tmp/"${PatientID}"_Aligned_Sorted_PCRDuped_FiltMAPQ.bam \
  "${Results_root}"/"${RunNumber}"/"${PatientID}"/tmp/"${PatientID}"_Aligned_Sorted_PCRDuped_FiltMAPQ_75pcReads.bam \
  "${Results_root}"/"${RunNumber}"/"${PatientID}"/tmp/"${PatientID}"_Aligned_Sorted_PCRDuped_FiltMAPQ_50pcReads.bam \
  "${Results_root}"/"${RunNumber}"/"${PatientID}"/tmp/"${PatientID}"_Aligned_Sorted_PCRDuped_FiltMAPQ_25pcReads.bam | \
    "${BCFtools_root}"bcftools call \
      --multiallelic-caller \
      --variants-only \
      -Ob > "${Results_root}"/"${RunNumber}"/"${PatientID}"/tmp/"${PatientID}"_Aligned_Sorted_PCRDuped_FiltMAPQ.bcf` ;

[ source: https://github.com/kevinblighe/ClinicalGradeDNAseq/blob/master/AnalysisMasterVersion1.sh#L302-L316 ]

----------------------

Extra tags can also be added 'after the fact' via BCFtools plugins, as I show here, in this previous question: A: How to use bcftools to calculate AF INFO field from AC and AN in VCF?

Kevin