FREEBAYES vcf file and Snpsift error
0
0
Entering edit mode
6.6 years ago

Hi!

I am back ! XD

Need a little help, now that I got happily my vcf file with freebaye I try to get my indels and snps, I found this program Snpsift pretty interesting and easy (if someone know another better would appreciate if tell...)

Well the fact when i try to get the indels I try with this code

cat /home/cri/Desktop/VCF/FREEBAYE.vcf | java -jar SnpSift.jar filter "( exists INDEL )" > indels.vcf


get this error message

Exception in thread "main" java.lang.RuntimeException: INFO field 'INDEL' not found in VCF header
at ca.mcgill.mcb.pcingola.snpSift.lang.expression.Field.getReturnType(Field.java:232)
at ca.mcgill.mcb.pcingola.snpSift.lang.expression.Field.eval(Field.java:56)
at ca.mcgill.mcb.pcingola.snpSift.lang.expression.ExpressionUnary.eval(ExpressionUnary.java:24)
at ca.mcgill.mcb.pcingola.snpSift.SnpSiftCmdFilter.evaluate(SnpSiftCmdFilter.java:133)
at ca.mcgill.mcb.pcingola.snpSift.SnpSiftCmdFilter.run(SnpSiftCmdFilter.java:299)
at ca.mcgill.mcb.pcingola.snpSift.SnpSiftCmdFilter.run(SnpSiftCmdFilter.java:269)
at ca.mcgill.mcb.pcingola.snpSift.SnpSift.run(SnpSift.java:372)
at ca.mcgill.mcb.pcingola.snpSift.SnpSift.main(SnpSift.java:70)


This is just because my vcf format don't allow me to get that or is something else Im missing?

Thanks!

vcf indel snp • 2.1k views
1
Entering edit mode

You can use awk '/^#/||length($4)!=length($5)' to extract indel lines. For freebayes, it is usually better to split complex variants with vt or something alike. By the way, it is freebayes, not freebaye.

1
Entering edit mode

related to the error, did your .vcf file have header? you can use head -n 29 /home/cri/Desktop/VCF/FREEBAYE.vcf cause from the error that is what the software says It can not find any header

0
Entering edit mode

mm seems like it hasnt right? why vcftools works fine? dont need header in that case?

head -n 29 freebayessample.vcf

## reference=/media/cri/CRIS_DATA/ULT/OTRO2/otro/hg19.fa

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> contig=<id=chr8,length=146364022> contig=ID=chr9,length=141213431 contig=ID=chr10,length=135534747 contig=ID=chr11,length=135006516 contig=ID=chr12,length=133851895 contig=ID=chr13,length=115169878 contig=ID=chr14,length=107349540 contig=ID=chr15,length=102531392 contig=ID=chr16,length=90354753 contig=ID=chr17,length=81195210 contig=ID=chr18,length=78077248 contig=ID=chr19,length=59128983 contig=ID=chr20,length=63025520 contig=ID=chr21,length=48129895 contig=ID=chr22,length=51304566 contig=ID=chrX,length=155270560 contig=ID=chrY,length=59373566

1
Entering edit mode

As you can see there is no header in file so when you use the command

-jar SnpSift.jar filter "( exists INDEL )"

there is no INDEL "which should be in header" so he can not extract or filter any thing,

why it works fine with vcftools; maybe the tool based on the coulm number, while snEff based on the column header name so for vcftools if there is no header there will not be any problem, but in case of snpEff there must be header or you will have this ugly error

0
Entering edit mode

thanks medhat

Anyhow later I will start over again all process to clear my mind and learn all the steps I have done....can get bit lost after the first time.... :)

now I just go forward until the end...generate my bed file atm

0
Entering edit mode

I experienced this problem using BBMap for alignment plus FreeBayes for variant calling. BBMap outputs v1.4 CIGAR format by default, which prevented FreeBayes from calling indels. Problem solved by adding the flag "sam=1.3" to BBMap.

Edit: It appears that OP's problem is a missing header, rather than indels not being called, so this post does not directly address the problem.

0
Entering edit mode

@cristina: Use reformat.sh from BBMap to convert your existing SAM v.1.4 tags to SAM v.1.3 in your BAM file. You do not need to re-do the alignments.

0
Entering edit mode

Thanks genomax, anyhow if I do that I need to redo all the steps afther alignment right?...... ^^"

At the end I found other tool that worked for me and called my snps and indels, if someone has the same problem as me, can use....

vcftools --vcf /home/cri/Desktop/VCF/vcftools/FREEBAYE.vcf --remove-indels --recode --recode-INFO-all --out /home/cri/Desktop/VCF/vcftools/outSNP.xls

vcftools --vcf /home/cri/Desktop/VCF/vcftools/FREEBAYE.vcf --keep-only-indels --recode --recode-INFO-all --out /home/cri/Desktop/VCF/vcftools/outindels.xls

thanks all for the tips and help...I need soon a big book to write all XD

0
Entering edit mode

thanks harold! I write the tip for next time I do aligment :)