Variant Calling With Bcftools Call

I have been trying to call variants using bcftools call after installing it as instructed here. The problem is that I am not getting any calls in the output vcf (it contains only the header).

My bam files were aligned with bwa-mem and processed with GATK pipeline. They have worked well with muTect and VarScan. I have tried many parameter combinations e.g. -DgB , for mpileup and -mv for bcftools, in vain. The versions I am using are samtools 0.2.0-rc7-54-ge35361b and bcftools 0.2.0-rc7-13-gca1a8b5

Does anyone know what the problem might be?

An example commandline I used is below:

samtools mpileup -u -f ref.fa -b sampleList | bcftools call -cv - > output.vcf

bcftools samtools mpileup variant-calling

• 26,390 views • link

updated 6.8 years ago by Ram

• written 12.4 years ago by tayebwajb

Without looking at things further, keep in mind that you're using a development version, so should presumably be reporting things back to the developers. You might get more stable results with the current release.

• 149 views • link

written 12.4 years ago by Devon Ryan

I think the problem is with bcftools since samtools mpileup looks fine to me. When I used the current release of bcftools, I got the error "invalid BCF2 magic string" which this link says was fixed in recent commits so that's why I tried the development version.

• 0 views • link

updated 6.8 years ago by Ram

• written 12.4 years ago by tayebwajb

What happens if you drop -v? Could the simple explanation be that there are no non-reference calls?

• 0 views • link

written 12.4 years ago by pd3

I tried without -v too but that didn't help either. It would be helpful to know which versions of samtools, htslib and bcftools are compatible since I've had a hard time figuring this out. Thanks!

• 0 views • link

updated 6.8 years ago by Ram

• written 12.4 years ago by tayebwajb

This does not seem to be a samtools/htslib/bcftools compatibility problem, as you say you followed instructions from here http://samtools.github.io/bcftools/. If you run mpileup with -uv, it will produce raw VCF output. Are any sites listed?

• 0 views • link

written 12.4 years ago by pd3

I used -uv on mpileup as you suggested and only -c on bcftools call and still got no sites listed.

• 0 views • link

written 12.4 years ago by tayebwajb

I meant if you could run without bcftools at all like this: samtools mpileup -uv -f ref.fa -b sampleList | less. Are there any sites listed? If not, then there is something odd with your BAMs. Could all reads be excluded because of low quality for example? Is the result any different with -B,-x,-Q0?

• 0 views • link

written 12.4 years ago by pd3

Yes, I see sites listed with just mpileup. The problem seems to be with bcftools.

• 0 views • link

written 12.4 years ago by tayebwajb

Can you please provide a sample of mpileup -uv output for testing?

• 0 views • link

written 12.4 years ago by pd3

I share tayebwajb issues and cannot have bcftools call to work on samtools mpileup data. I have samtools and bcftools/htslib fresh installed and also get the complain about BCF2.

samtools mpileup -u -f ${ref} ${bam} | bcftools call -c > raw_calls.bcf
[mpileup] 1 samples in 1 input files
<mpileup> Set max per-file depth to 8000
[vcf.c:643 bcf_hdr_read] invalid BCF2 magic string: only BCFv2.2 is supported.

while my samtools pileup seems OK

samtools mpileup -u -f ${ref} ${bam} | head
[mpileup] 1 samples in 1 input files
<mpileup> Set max per-file depth to 8000
?BC???BCF?chrMchr1chr2chr3chr4chr5chr6chr7chr8chr9chr10chr11chr12chr13chr14chr15chr16chr17chr18chr19chr20chr21chr22chrXchrYGAIIx-chr21-BWA.aln?##samtoolsVersion=0.1.19-44428cd
##reference=file://ref/HiSeq_UCSC_hg19.fa
##contig=<ID=chrM,length=16571>
##contig=<ID=chr1,length=249250621>
##contig=<ID=chr2,length=243199373>
##contig=<ID=chr3,length=198022430>
##contig=<ID=chr4,length=191154276>
##contig=<ID=chr5,length=180915260>
##contig=<ID=chr6,length=171115067>
##contig=<ID=chr7,length=159138663>

A full example command would be great with data we can download somewhere to test our executables. Any help most welcome

Stephane

• 0 views • link

updated 6.8 years ago by Ram

• written 12.3 years ago by stephane.plaisance •

Almost there (thanks pd3 for pointing the git branch error of mine)

Now I still have an issue building razip, likely a mac OSX maverick specific option missing in the makefile

make razip
gcc -g -Wall -O2 -D_FILE_OFFSET_BITS=64 -D_LARGEFILE64_SOURCE -D_CURSES_LIB=1 -I. -I../htslib   razip.c   -o razip
Undefined symbols for architecture x86_64:
  "_razf_close", referenced from:
      _main in ccPCU9N4.o
  "_razf_dopen", referenced from:
      _main in ccPCU9N4.o
  "_razf_open", referenced from:
      _main in ccPCU9N4.o
  "_razf_read", referenced from:
      _main in ccPCU9N4.o
  "_razf_seek", referenced from:
      _main in ccPCU9N4.o
  "_razf_write", referenced from:
      _main in ccPCU9N4.o
ld: symbol(s) not found for architecture x86_64
collect2: error: ld returned 1 exit status
make: *** [razip] Error 1

Thanks for help on this!

• 0 views • link

updated 6.8 years ago by Ram

• written 12.3 years ago by stephane.plaisance •

Hi Stephane,

I have the same problem with bcftools not recognizing mpileup output and reporting the same "invalid BCF2 magic string" error. I have updated all versions of samtools, bcftools and htslib but still the same error

Have you managed to get over this issue?

Thank you!

Sahar

• 0 views • link

updated 6.8 years ago by Ram

• written 12.3 years ago by sahargel •

5 answers

There is not enough information to answer the original question, but Stephane's problem is caused by a version mismatch. Note that samtools and bcftools went through major changes, the new versions of bcftools will not work with pre-HTSlib versions of samtools. For the latest installation instructions please see bcftools' wiki, the current commands (which are likely to change in near future) are:

git clone --branch=bcftools+calling git://github.com/samtools/htslib.git
git clone git://github.com/samtools/bcftools.git
git clone git://github.com/samtools/samtools.git
cd bcftools; make; cd ..
cd samtools; make;

• 0 views • link

updated 6.8 years ago by Ram

• written 12.3 years ago by pd3

You need to use samtools 1.3 if you're using bcftools 1.3.

• 0 views • link

written 10.5 years ago by Devon Ryan

How about using GATK's Unified Genotyper as an alternative to samtools mpileup? I have shared a tutorial here (C. Difficile 630 dataset that I have): http://userweb.eng.gla.ac.uk/umer.ijaz/bioinformatics/oneliners.html?#VARIANTCALLING

• 0 views • link

updated 6.8 years ago by Ram

• written 12.3 years ago by umer.zeeshan.ijaz

Thanks Umer but the point IS to use a second caller and see what is shared with GATK.

• 0 views • link

written 12.3 years ago by stephane.plaisance •

My apologies, I might have misinterpreted while answering in a hurry!

Best Wishes,

Umer

• 0 views • link

updated 6.8 years ago by Ram

• written 12.3 years ago by umer.zeeshan.ijaz

no problem and thanks for sharing your methods ++

• 0 views • link

written 12.3 years ago by stephane.plaisance •

Hi,

Has this been resolved. I have done fresh installs of samtools and bcftools and get:

[vcf.c:753 bcf_hdr_read] invalid BCF2 magic string: only BCFv2.2 is supported

Its very frustrating. Also I cannot use GATK, I can't get Java 1.7 on my Debian system (yet) and it always gives me numerous java errors.

Thanks.

• 0 views • link

updated 6.8 years ago by Ram

• written 11.8 years ago by susan.klein •

There have been two updates to samtools and bcftools since this post was made, so it's quite likely that the issue has been resolved. This is particularly so since tayebwajb was using a pre-release version.

• 0 views • link

written 11.8 years ago by Devon Ryan

Hi all!

After running:

samtools mpileup -uf ../samples/human_g1k_v37.fasta ../samples/sorted_sample.bam | bcftools call -c

I am also getting the same error:

[vcf.c:805 bcf_hdr_read] invalid BCF2 magic string: only BCFv2.2 is supported.
Failed to open -: could not parse header

Latest versions of tools (bcftools 1.3, samtools 0.1.19-96b5f2294a). Does anyone has any news on this one nowadays?

• 0 views • link

updated 6.8 years ago by Ram

• written 10.5 years ago by vladimir.kovacevic

That did the trick! Thanks!

• 0 views • link

written 10.5 years ago by vladimir.kovacevic

Hi Vladimir!

It seems you had exactly the same problem I have now. When I run bcftools with the call command to convert a bcf file into vcf I get the following error message:

[vcf.c:842 bcf_hdr_read] invalid BCF2 magic string: only BCFv2.2 is supported

Could you please tell me how you've managed to solve this issue.

Thanks!

• 0 views • link

updated 6.8 years ago by Ram

• written 9.8 years ago by vasilipankratov •

You need to use samtools 1.3 if you're using bcftools 1.3. for samtools-0.1.18 and samtools-0.1.19, you should use the bcftools under the samtools build folder like "/share/apps/prog/samtools-0.1.19/bcftools/bcftools"

• 0 views • link

updated 6.8 years ago by Ram

• written 9.6 years ago by hlzi •