Using Samtools mpileup with Varscan
2
0
Entering edit mode
9.3 years ago
Ron ★ 1.2k

Hi all,

I am using samtools mpile up and varscan to find variants.below is my command:

samtools mpileup \
  -f hg19.fa \
  -D 8000 \
  -S \
  -l stat3_jak1.bed \
  -b $input_dir/tcell_bam_list | \
    /bin/java -jar mutation_tools/VarScan.v2.3.7.jar mpileup2cns \
    --vcf-sample-list $input_dir/tcell_bam_list \
    --min-coverage 6 \
    --min-var-freq 0.08 \
    --p-value 0.01 \
    --output-vcf 1 > $output_dir/mutation_analysis_result.vcf

My bed file looks like(for two genes,STAT3 and JAK1)

chr17    40465342    40540513    uc002hzk.1    0    -    40468844    40500534    0    24    "2476,63,43,213,140,95,53,136,99,84,48,94,30,60,93,159,152,95,82,96,99,145,151,217,"    "0,3464,3857,8957,9679,9935,10248,11386,11638,12791,16085,16229,16422,18147,20348,20566,24110,24438,25406,25989,3223
chr1744,40465342954,"   40540513    uc002hzl.1    0    -    40467762    40500534    0    24    "2476,113,43,213,140,95,53,136,99,84,48,94,30,60,93,159,152,95,82,96,99,145,151,217,"    "0,3464,3857,8957,9679,9935,10248,11386,11638,12791,16085,16229,16422,18147,20348,20566,24110,24438,25406,25989,3223
chr1744,40465342954,"   40540513    uc002hzm.1    0    -    40467762    40500534    0    24    "2476,113,43,210,140,95,53,136,99,84,48,94,30,60,93,159,152,95,82,96,99,145,129,217,"    "0,3464,3857,8960,9679,9935,10248,11386,11638,12791,16085,16229,16422,18147,20348,20566,24110,24438,25406,25989,3223
chr1744,40465342954,"   40540513    uc010wgh.1    0    -    40467762    40497654    0    23    "2476,113,43,213,140,95,53,136,99,84,48,94,30,60,93,159,152,95,82,96,99,145,217,"    "0,3464,3857,8957,9679,9935,10248,11386,11638,12791,16085,16229,16422,18147,20348,20566,24110,24438,25406,25989,3223
chr1744,40465342    40540513    uc002hzn.1    0    -    40467762    40500534    0    24    "2476,113,43,213,140,95,53,136,99,84,48,94,30,60,93,159,152,95,82,96,99,145,129,217,"    "0,3464,3857,8957,9679,9935,10248,11386,11638,12791,16085,16229,16422,18147,20348,20566,24110,24438,25406,25989,3223
4,33244,35064,74954,"

tcell_bam_list contains the list of all bam files.

I get this following error:

Min coverage:    6
Min reads2:    2
Min var freq:    0.08
Min avg qual:    15
P-value thresh:    0.01
[mpileup] 24 samples in 24 input files
<mpileup> Set max per-file depth to 333
Input stream not ready, waiting for 5 seconds...
Input stream not ready, waiting for 5 seconds...
Input stream not ready, waiting for 5 seconds...
Input stream not ready, waiting for 5 seconds...
Input stream not ready, waiting for 5 seconds...
ERROR: Gave up waiting after 500 seconds...
next-gen RNA-Seq varscan samtools • 7.5k views
ADD COMMENT
0
Entering edit mode

I can't see any error; Varscan is just waiting for the first variants....

ADD REPLY
0
Entering edit mode

It gives up after a certain time.

ADD REPLY
0
Entering edit mode

Try -mpileup instead mpileup2cns

ADD REPLY
1
Entering edit mode
4.3 years ago
field654 ▴ 30

I had the same problem and figured it out.

That's just VarScan's fancy way of saying, can't find your input.

Very frustrating.

ADD COMMENT
0
Entering edit mode

This is it!! Registered on BioStars just to upvote you. Thanks!

ADD REPLY
0
Entering edit mode

Thank you, that's so encouraging. As a nonprofessional, it's especially enjoying to help.

ADD REPLY
0
Entering edit mode
9.3 years ago

Which VarScan command are you using? I know that there are problems with using pileup2snp/pileup2indel if you used mpileup (rather than pileup).

I don't think this is the error message that I would expect, but I do have a pipeline that includes variant calling with VarScan (if you are willing to modify some Perl code): https://docs.google.com/file/d/0B1xpw6_kQMKuaElwcTVXcnFmZ1E

ADD COMMENT

Login before adding your answer.

Traffic: 1857 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6