freebayes with ogap

[DNAjunk]

Hi !

I am trying to use freeBayes for variant calling and run into problems with the example provided in the freeBayes package.

When I run the following pipeline given in the package, I get the warnings stated below for most of the reads in my bamfiles

#!/bin/bash

outdir=OUTPUT
mkdir -p $outdir
reference=REFERENCE.fasta
bamlist=BAMFILENAMES.txt
cnvmap=CNVMAP.bed
region=REGION

time /share/home/erik/bin/bamtools merge -region $region \
$(for file in $(cat $bamlist); do echo " -in "$file; done) \
| time /share/home/erik/bin/ogap -z -R 25 -C 20 -Q 20 -S 0 -f $reference \
| time /share/home/erik/bin/bamleftalign -f $reference \
| time /share/home/erik/bin/samtools calmd -EAru - $reference 2>/dev/null \
| time /share/home/erik/bin/freebayes \
--min-alternate-count 2 \
--min-alternate-qsum 40 \
--pvar 0.0001 \
--use-mapping-quality \
--posterior-integration-limits 1,3 \
--genotype-variant-threshold 4 \
--site-selection-max-iterations 3 \
--genotyping-max-iterations 25 \
--max-complex-gap 3 \
--cnv-map $cnvmap \
--stdin \
--region $region \
-f $reference \
| gzip >$outdir/$region.vcf.gz



the output consists of many warnings (exceptions, basically for every read in my bamfiles) of the form

exception when realigning HWI-D00432:140126:C3KE2ACXX:1:2301:14872:89688 at position chr1:236849712 GTCCGCCCGCACGTGCTCGGGGCGGGGCGGCTGGGGCCCGACGGGCGGGGCCGGGGGCGGTGGCGGCGGGGGCGTGGCCGGGCGGCCACGTGACGGGCGGCGCGGCTATTAAGCCGCGCGG CCGGGGGCGGTGGCGGCGGGGG
CGCACGTGCTCGGGGCGGGGCGGCTGGGGCCCGACGGGCGGGGCCGGGGGCGGTGGCGGCGGGGGCGTGGCCGGGCGGCCACGTGACGGGCGGCGCGGCTATTAAGCCGCGCGGCAGCTGCTCGCAGCCGGAGCTGGTGCTT
CGGTGGCGGCGGGGGCGTGGCCGGGCGGCCACGTGACGGGCGG


I could not figure out why this happens and how to solve the problem.
It seems that the message is generated by ogap.cpp .
Any help is appreciated.
Thanks!

[hugorody]

Hey,

you should try these steps:

#index the reference fasta
samtools faidx ref.fasta

#convert all SAM files to BAM files by using samtools and the reference fasta indexed (*.fasta.fai)
samtools view -bt ref.fasta.fai alignment.sam > alignment.bam

#sort the BAM files
samtools sort alignment.bam $alignment.sorted

#index BAM files
samtools index $indexar

#apply BAQ
samtools calmd -Abr alignment.sorted.bam ref.fasta > alignment.baq.bam

#index the *.baq.bam files
samtools index alignment.baq.bam

#SNP call using the *.baq.bam files
freebayes --min-base-quality 20 -p 1 --fasta-reference ref.fasta alignment.baq.bam > snpcall.vcf

[DNAjunk]

Thanks for the advices!
I also got the information from the author of the program that one should use freebayes without the ogap program according to the freebayes manual.