This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Varscan output variant file - is it the "true" file?

I want to run samtools and pipe the mpileup output into a VarScan call using the following conditions:

  • 30 mapping score
  • 30 phred quality score
  • output-snp format
  • min coverage = 50
  • min tumor coverage = 3
  • p-value = 0.01
  • min variant frequency = 4%
`samtools mpileup -B -q 30 -Q 30 -o normal-tumor.mpileup normal_sort.bam tumor_sort.bam | java -jar VarScan.v2.3.9.jar somatic normal-tumor.mpileup output --mpileup 1 --strand-filter 1 --output-snp --min-coverage 50 --min-coverage-tumor 3 --p-value 0.01 --min-var-freq 0.04`

Output:

[mpileup] 2 samples in 2 input files
Min coverage:   50x for Normal, 3x for Tumor
Min reads2:     2
Min strands2:   1
Min var freq:   0.04
Min freq for hom:       0.75
Normal purity:  1.0
Tumor purity:   1.0
Min avg qual:   15
P-value thresh: 0.01
Somatic p-value:        0.05
Reading input from /dev/stdin
Reading mpileup input...
27242828 positions in mpileup file
929562 had sufficient coverage for comparison
0 were called Reference
0 were mixed SNP-indel calls and filtered
0 were removed by the strand filter
929562 were called Germline
0 were called LOH
0 were called Somatic
0 were called Unknown
0 were called Variant

I'm expecting an output variant file but I'm only getting true and output.indel files. Is the true file the output variant file?

varscan
• 1,786 views • link
updated by ATpoint ATpoint • written by melissachua90 • melissachua90

Please spend some time learning how Unix pipes work. In a nutshell, if you pipe something then you do not name the output. So normal-tumor.mpileup should not appear in the first, nor the second command part. That is the whole idea of a pipe, not to save things to a file. VarScan2, despite being old, has a manual page, http://varscan.sourceforge.net/somatic-calling.html, I suggest you follow it.

• 0 views • link
written by ATpoint ATpoint

When normal-tumor.mpileup is not specified, I get a traceback error:

[mpileup] 1 samples in 1 input files
<mpileup> Set max per-file depth to 8000
Min coverage:   50x for Normal, 3x for Tumor
Min reads2:     2
Min strands2:   1
Min var freq:   0.04
Min freq for hom:       0.75
Normal purity:  1.0
Tumor purity:   1.0
Min avg qual:   15
P-value thresh: 0.01
Somatic p-value:        0.05
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...

Code:

samtools mpileup -B -q 30 -Q 30 -o normal_sort.bam tumor_sort.bam | java -jar Var
Scan.v2.3.9.jar somatic --mpileup 1 --strand-filter 1 --output-snp --min-coverage 50 --min-coverage-tumor 3 --p-value 0.01 --min-var-freq 0.04

Which basically means that VarScan does not recognize the input file. Discussed here: Using Samtools mpileup with Varscan https://sourceforge.net/p/varscan/discussion/1073559/thread/27563361/

• 0 views • link
written by melissachua90 • melissachua90

Googling for "varscan somatic pipe" would have also shown you my previous answer on this: How to get "normal.pileup" file for varscan somatic calling Again, please learn the principles of Unix and pipes. If a tool expects an input and it can be replaced by stdin then you often have to declare that with either - or something like /dev/stdin.

• 0 views • link
written by ATpoint ATpoint

Do you mean something like this:

samtools mpileup -B -q 30 -Q 30 -f hg19.fa normal_sort.bam tumor_sort.bam | java -jar VarScan.v2.3.9.jar somatic /dev/stdin output --mpileup 1 --strand-filter 1 --output-snp --min-coverage 50 --min-coverage-tumor 3 --p-value 0.01 --min-var-freq 0.04
• 0 views • link
written by melissachua90 • melissachua90

0 answers

No answers yet.

Log in to answer this question.

More posts like this