Varscan fpfilter error

Question

[a01039012387]

Hello. I'm new to Varscan. I want to call the SNPs from the bam file using varscan in the order of somatic - processSomatic-fpfilter but I got an error message in the fpfilter process.

Here is my command.

perl fpfilter.pl --var-file [snp.Somatic.hc] --readcount-file [readcount] --output-file [fiinal]

and this is the error message.

Failed to interpret variant allele: V at fpfilter.pl line 333, <GEN2> line 1.

It seems to be a problem during the execution of the following command in fpfilter.

 **## iupactobase - Convert IUPAC ambiguity codes to variant bases ##

 sub iupactobase {
     my ( $allele1, $allele2 ) = @;

     return( $allele2 ) if( $allele2 eq "A" or $allele2 eq "C" or $allele2 eq "G" or $allele2 eq "T" );

     # Choose the most likely base-pair, or the default for triallelic variants
 if( $allele2 eq "M" ) {
      return( "C" ) if( $allele1 eq "A" );
      return( "A" ) if( $allele1 eq "C" );
      return( "A" );
  } elsif( $allele2 eq "R" ) {
      return( "G" ) if( $allele1 eq "A" );
      return( "A" ) if( $allele1 eq "G" );
      return( "A" );
   } elsif( $allele2 eq "W" ) {
       return( "T" ) if( $allele1 eq "A" );
       return( "A" ) if( $allele1 eq "T" );
       return( "A" );
   } elsif( $allele2 eq "S" ) {
      return( "C" ) if( $allele1 eq "G" );
      return( "G" ) if( $allele1 eq "C" );
      return( "C" );
   } elsif( $allele2 eq "Y" ) {
       return( "C" ) if( $allele1 eq "T" );
       return( "T" ) if( $allele1 eq "C" );
       return( "C" );
   } elsif( $allele2 eq "K" ) {
       return( "G" ) if( $allele1 eq "T" );
       return( "T" ) if( $allele1 eq "G" );
       return( "G" );
   }

   die "Failed to interpret variant allele: $allele2";
}*

But the input file does not seem to be a problem. Thus, I don't know what is a problem. Of course, there is no allele V in the input file. Below is the input file for fpfilter.

**readcount**

> 1 1225707 C 27 =:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00
> A:4:36.75:30.50:31.00:4:0:0.42:0.11:187.25:4:0.01:45.00:0.01
> C:23:56.30:37.22:36.30:4:19:0.51:0.01:402.78:18:0.38:58.39:0.38
> 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 1 6653482 G 31
> =: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:6:49.67:25.17:30.33:0:6:0.31:0.06:298.33:6:0.02:46.17:0.02
> G:25:58.76:29.52:36.68:22:3:0.52:0.00:564.20:25:0.64:71.88:0.64
> 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 1 9791175 T 32
> =:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00 A:14:57.79:32.21:36.43:14:0:0.64:0.03:368.86:14:0.32:55.07:0.32
> C:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00
> G:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00
> T:18:58.28:30.78:36.56:18:0:0.53:0.01:408.78:18:0.33:54.33:0.33
> N:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00

**snp.Somatic.hc**

> chrom position ref var normal_reads1 normal_reads2 normal_var_freq
> normal_gt tumor_reads1 tumor_reads2 tumor_var_freq tumor_gt
> somatic_status variant_p_value somatic_p_value tumor_reads1_plus
> tumor_reads1_minus tumor_reads2_plus tumor_reads2_minus
> normal_reads1_plus normal_reads1_minus normal_reads2_plus
> normal_reads2_minus 1 1225707 C A 30 0 0% C 23 4 14.81% M Somatic 1
> 0.044429255 4 19 4 0 5 25 0 0 1 6653482 G C 19 0 0% G 25 6 19.35% S Somatic 1 0.046334082 22 3 0 6 14 5 0 0

[ATpoint]

As you use the standalone fpfilter, this is probably an ancient version. Get the VarScan2.4.3 java executable from Github and run java -jar VarScan2 fpfilter. Edit: As you say you are new to varscan, if you want you can post your entire command from the aligned reads to the final file, and we can see if you used to most recent (and recommended) settings.

[a01039012387]

Oh, thanks for the reply.

Here is my entire commend for Varscan, and I already using varscan2.4.3.

> export GENOME=/home/bin/genome.fa
> samtools mpileup -f $GENOME -q 1 -B $path_new_bam_normal $path_new_bam_tumor > $path_mpileup
> java -jar VarScan.v2.4.3.source.jar somatic $path_mpileup --mpileup 1 --min-coverage 3 --min-var-freq 0.08 --somatic-p-value 0.05 --strand-filter 0 --output-snp $path_snp_indel --output-indel $path_snp_indel
> java -jar VarScan.v2.4.3.source.jar copynumber $path_mpileup $path_copynumber --mpileup 1
> java -jar VarScan.v2.4.3.source.jar processSomatic $path_snp_indel --min-tumor-freq 0.08
> bam-readcount -q 1 -b 20 -l $path_snp_hc_var -f /home/bin/genome.fa $path_new_bam_tumor > $path_readcount
> perl fpfilter.pl --var-file $path_snp_hc --readcount-file $path_readcount --output-file $path_final

and I run java -jar VarScan.v2.4.3.source.jar fpfilter /home/varscan/snp_indel/1.snp.Somatic.hc /home/varscan/snp_indel/1.readcount --output-file /home/varscan/snp_indel/1.final --dream3-settings 1 It works! but output file doesn't have any SNPs.

> 795 variants in input file 795 had a bam-readcount result 741 had
> reads1>=2 0 passed filters 795 failed filters
>         0 failed because no readcounts were returned
>         12 failed minimim variant count < 3
>         6 failed minimum variant freq < 0.05
>         0 failed minimum strandedness < 0.0
>         28 failed minimum reference readpos < 0.2
>         95 failed minimum variant readpos < 0.15
>         11 failed minimum reference dist3 < 0.2
>         7 failed minimum variant dist3 < 0.01
>         723 failed maximum reference MMQS > 30
>         791 failed maximum variant MMQS > 30
>         59 failed maximum MMQS diff (var - ref) > 50
>         31 failed maximum mapqual diff (ref - var) > 30
>         22 failed minimim ref mapqual < 20
>         108 failed minimim var mapqual < 30
>         0 failed minimim ref basequal < 15
>         0 failed minimim var basequal < 5
>         254 failed maximum RL diff (ref - var) > 0.05

[ATpoint]

Sorry for coming back so late. Well what shall I say, looks like all your initial variants are junk calls, at least based on the fpfilter. What kind of experiment was that?