}
sub varFilter {
- my %opts = (d=>2, D=>10000, a=>2, W=>10, Q=>10, w=>10, p=>undef, 1=>1e-4, 2=>1e-100, 3=>0, 4=>1e-4, G=>0, S=>1000);
+ my %opts = (d=>2, D=>10000000, a=>2, W=>10, Q=>10, w=>10, p=>undef, 1=>1e-4, 2=>1e-100, 3=>0, 4=>1e-4, G=>0, S=>1000);
getopts('pd:D:W:Q:w:a:1:2:3:4:G:S:', \%opts);
die(qq/
Usage: vcfutils.pl varFilter [options] <in.vcf>
print; next;
}
next if ($t[4] eq '.'); # skip non-var sites
+ next if ($t[3] eq 'N'); # skip sites with unknown ref ('N')
# check if the site is a SNP
my $type = 1; # SNP
if (length($t[3]) > 1) {
$flt = 1 if ($flt == 0 && $mq >= 0 && $mq < $opts{Q});
$flt = 7 if ($flt == 0 && /PV4=([^,]+),([^,]+),([^,]+),([^,;\t]+)/
&& ($1<$opts{1} || $2<$opts{2} || $3<$opts{3} || $4<$opts{4}));
- $flt = 8 if ($flt == 0 && ((/MXGQ=(\d+)/ && $1 >= $opts{G}) || (/MXSP=(\d+)/ && $1 < $opts{S})));
+ $flt = 8 if ($flt == 0 && ((/MXGQ=(\d+)/ && $1 < $opts{G}) || (/MXSP=(\d+)/ && $1 >= $opts{S})));
my $score = $t[5] * 100 + $dp_alt;
my $rlen = length($t[3]) - 1; # $indel_score<0 for SNPs
} else { # SNP or MNP
for my $x (@staging) {
next if (($x->[0]&3) != 3 || $x->[4] + $x->[2] + $ow < $t[1]);
- $flt = 5;
+ if ($x->[4] + length($x->[7]) - 1 == $t[1] && substr($x->[7], -1, 1) eq substr($t[4], 0, 1)
+ && length($x->[7]) - length($x->[6]) == 1) {
+ $x->[1] = 5;
+ } else { $flt = 5; }
last;
}
# check MNP
$q = chr($q <= 126? $q : 126);
$seq .= $b;
$qual .= $q;
- } else { # an INDEL
+ } elsif ($t[4] ne '.') { # an INDEL
push(@gaps, [$t[1], length($t[3])]);
}
$last_pos = $t[1];