Question: Problem in using VarScan fpfilter
gravatar for Rahil
23 months ago by
Rahil170 wrote:

Hi, I ran bam-readcount to get coverage per nucleotide. I got the result after getting this warning `

“WARNING: In read DHDC08P1_0264:2:2204:2494:140080#CTTGTA: Couldn't find single-end mapping quality. Check to see if the SM tag is in BAM”

This was the command:

bam-readcount -f ../index/hg38.fa -q 1 -b 20 -l all.Germline.hc.var normal.sorted.bam > all.Germline.hc.readcount

` and this is a line of the output file:

chr1    931131  C   14  =:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  A:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  C:14:37.00:34.86:7.00:2:12:0.38:0.03:9.93:2:0.19:123.86:0.72    G:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  T:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  N:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  +CCCT:7:42.00:0.00:0.00:1:6:0.65:0.04:7.86:1:0.26:122.00:0.60

When I used the readcount file to filter out false positives from snp calling, I encountered this error:

Loading readcounts from all.Germline.hc.readcount...
Parsing variants from all.Germline.hc.var...
Exception thrown while filtering at chr1    931131: 3
java.lang.ArrayIndexOutOfBoundsException: 3
    at net.sf.varscan.FpFilter.<init>(
    at net.sf.varscan.VarScan.fpfilter(
    at net.sf.varscan.VarScan.main(

This was the command: java -jar VarScan.v2.3.9.jar fpfilter all.Germline.var all.Germline.readcount --min-ref-basequal 20 --min-var-basequal 20 --output-file all.Germline.hc.filtered.vcf --dream3-settings

Where is the problem?

Thanks very much for any help!

ADD COMMENTlink modified 22 months ago by Biostar ♦♦ 20 • written 23 months ago by Rahil170

all.Germline.hc.var must be a BED file. Also, you are not using the latest release of Varscan, which is I think 2.4.3. Better get the latest one from Github.

ADD REPLYlink modified 22 months ago • written 23 months ago by ATpoint40k

Thank you for your reply. As you suggested I got the latest release of VarScan( 2.4.3) but still get the same error. I think the problem is not related to my all.Germline.hc.var. As I found a BED file is a file with three columns (Chromosome Stop Start) and this is my all.Germline.hc.var file:

chr1    1046480 1046480
chr1    1182186 1182186
chr1    1217483 1217483

This is the command that I used to get the file:

awk 'BEGIN {OFS="\t"} {if (!/^#/) { isDel=(length($4) > length($5)) ? 1 : 0; print $1,($2+isDel),($2+isDel); }}' *.Germline.hc.vcf > all.Germline.hc.var

I am following this instruction link.


ADD REPLYlink modified 23 months ago • written 23 months ago by Rahil170
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1426 users visited in the last hour