Entering edit mode
6.0 years ago
lemon
•
0
I got my sequence data, and I want to call its snp and indel when reference genome exists. The command I used is as follows:
bcftools mpileup -Ou -f /annovar/Gadb/Ga.fa /annovar/Gadb/0613.bam | \
bcftools call -Ou -mv | \
bcftools filter -s LowQual -e '%QUAL<20 || DP>100' > 0613.flt.vcf
This is a whole genome, and it should contain indels but the results I got form this commond only contain snps , however, no indels.
Could anyone help me solve this problem?
Try running it without the LowQual filter?
Hello lemon ,
how have you checked that there are no indels in your vcf?
Visualize your bam file with e.g. igv to find a region where a clear insdel is. Retry the variant calling for only this region.
fin swimmer
thanks , when i try to use IGV, i met the same problem that indel were thought to be many snps as you can see in the following picture. it seems that actually a ''G'' was inserted
do u know how to deal with this problem? thank you very much.
Indeed this looks a bit strange. So the problem here isn't
bcftools
, because it rely on the alignment information. Which program and command have you used for alignment?I used this command to make alignment. thanks
Hello and thanks for posting it.
I don't know
biokanga
. I would trybwa mem
and see what happens then. If this works, one can go back to biokanga (if you still want/need to use it) and start finding out what's going on there.fin swimmer