This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Variant calling using samtools

I have cancer dataset containing 10 tumor and 10 control pairs. Each tumor or control dataset is 100 GB in size. I have a refrence sequence too. (mm9.fa). SO I need to do some variant calling using these available data. What I do is the following:-

samtools mpileup -g -f mm9.fa *.bam | bcftools view -bvcg - > var.raw.bcf (mpileup takes all the tumor and control pairs) bcftools view var.raw.bcf | vcfutils.pl varFilter -D100 > var.flt.vcf

But the entire process is agonizingly slow (>10 hrs and still nothing). What should I do? Can I find the variants individually and then merge them into one large bcf file? P.S: I am very new to this field and pardon my ignorance in some words written above. Thanks in advance.

alignment

2 answers

you have to try it with submission to a computer cluster. It will save you some time and computer access.

BBMap's variant-caller is roughly 80x faster than the pipeline you are using. Samtools has to be in your path (and I recommend samtools 1.4, which is much faster than older versions). The command for multiple samples would be:

callvariants.sh in=sample1.bam,sample2.bam,sample3.bam ref=mm9.fa ploidy=2 out=variants.vcf multisample

Thanks for the suggestion. I have 20 bam files. Can I write *.bam here? Its giving an error here. Or do I need to write each and every file name?

Sorry :) I'll modify it to allow *.bam but right now you have to list all of them. Also, please note I modified the command to add the term "multisample".

Thanks! Will do. This really helped.

Log in to answer this question.