Question: Zero-Read Depth In Samtools Mpileup Causing Problems In Varscan2 With Somatic Module
0
gravatar for DoubleD
5.5 years ago by
DoubleD130
United States
DoubleD130 wrote:

hello. I searched far and wide for a solution to this problem, and saw a post about it maybe being fixed in VarScan 2.3.6, but I'm hainvg it nonetheless.

i have been compiling mpileups for use in Varscan2. I have been having some troubles, I think because my mpileups of normal&tumor have lines with zero-depth reads. I used a reference for alignment, with the commands

samtools mpileup -B -q 1 -f human_g1k_v37.fasta normal.bam tumor.bam > test.mpileup

but when running in Varscan with

java -jar VarScan.v2.3.6.jar somatic test.mpileup outputfile --mpileup 1

I get

"Parsing exception on line: 1 9999 N 0 * * 0"

or if I remove that line, then the next 0 read position fails

"1 10000 N 0 0"

or remove that and the next...

"1 29418 G 0 * * 0"

I'm wondering if there were something else I could do. The first three lines of the mpileup are as follows:

1 9999 N 0 * * 0
1 10000 N 0 0
1 10001 T 9 .^>.^2.^".^2.^%.^-.^+.^2. <@B>>@BBA 28 .^>.^2.^+.^> etc

Thank you, DD

ADD COMMENTlink modified 5.5 years ago by Biostar ♦♦ 20 • written 5.5 years ago by DoubleD130
1

This seems to be an ongoing issue with VarScan2 (there are other similar threads both here and on SEQanswers). Allegedly, most of these were supposed to be fixed in the version you're using, but...

If nothing else, you might just pass the mpileup through awk to remove the 0 coverage lines (something like samtools mpileup -q 1 -f human_g1k_v37.fasta normal.bam tumor.bam | awk '{if($4 != 0) print $0}' > test.mpileup).

ADD REPLYlink modified 5.5 years ago • written 5.5 years ago by Devon Ryan89k

Great suggestion, thank you. I'm not missing anything by filtering these spots, since there are zero reads there. If I filter those lines, then it shows an Exception at lines with 1 read, and filtering those, wtih 2 reads, etc. There's the default option that it won't call a variant below 8 reads in tumor or 6 in the normal, so perhaps for SNPs filtering everything below 6 would work.

I had several reads with 0 in the last column (the base quality)... so I filtered those too

ADD REPLYlink modified 5.5 years ago • written 5.5 years ago by DoubleD130
1

Dan Koboldt, the author of Varscan, is away at a conference right now, but I passed this along to him - he may be able to offer more insight or a bug fix when he returns.

ADD REPLYlink written 5.5 years ago by Chris Miller20k

Thanks for passing it along. I ended up using

cat /media/BAM/111_N_T.mpileup | awk '{if($4 >= 6) print $0}' | awk '{if($7 != 0) print $0}' | java -jar ~/Downloads/VarScan.v2.3.6.jar somatic -mpileup 1

as it was giving read errors up to 6 reads and then also even though I used -q 1 with samtools to make the pileup, there were mapping qualities of 0 so I had to eliminate the lines with a 0 at the end

ADD REPLYlink written 5.5 years ago by DoubleD130

Yikes, that's totally an implementation mistake in VarScan2 then, it really should not be doing that. You might try filing a bug report on the sourceforge page for varscan2 (there's a forum there as well, though it doesn't seem to be very active).

ADD REPLYlink written 5.5 years ago by Devon Ryan89k
1

FYI, TGI automatically monitors Biostar and SEQanswers to try and pick up questions about VarScan and our other tools.

ADD REPLYlink written 5.5 years ago by nnutter210
Please log in to add an answer.

Help
Access

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