diff --git a/variants/filter_svaba.pl b/variants/filter_svaba.pl index 0068be0f51812daf5f94cf10e683351f28122736..116f7da9956667fc1871165dc1f859e7e05deab5 100755 --- a/variants/filter_svaba.pl +++ b/variants/filter_svaba.pl @@ -114,7 +114,7 @@ W1:while (my $line = <IN>) { } print VCFOUT join("\t",$chrom,$pos,$id,$ref,$alt,$score,$filter,$newannot, $format,@gts),"\n"; - }elsif ($hash{SVTYPE} eq 'DEL' && $hash{SPAN} && $hash{SPAN} > 9999) { + } elsif ($hash{SVTYPE} eq 'DEL' && $hash{SPAN} && $hash{SPAN} > 9999) { my ($allele,$effect,$impact,$gene,$geneid,$feature, $featureid,$biotype,$rank,$codon,$aa,$pos_dna,$len_cdna, $cds_pos,$cds_len,$aapos,$aalen,$distance,$err) = split(/\|/,$keeptrx); @@ -230,6 +230,9 @@ close IN; foreach my $id (keys %svpairs) { + if ($id =~ m/815016443/) { + warn "debugging\n"; + } my $alt1 = $svpairs{$id}{1}{alt}; my $alt2 = $svpairs{$id}{2}{alt}; my $svtype; @@ -237,8 +240,10 @@ foreach my $id (keys %svpairs) { $svtype = 'DEL'; }elsif ($alt2 =~ m/^\w\[/ && $alt1 =~ m/^\]/) { $svtype = 'INS'; + }else { + $svtype = 'UNK'; } - if ($svtype eq 'INS' || ($svtype eq 'DEL' && $svpairs{$id}{1}{gene} !~ m/&/)) { + if ($svtype eq 'INS' || ($svtype eq 'DEL' && $svpairs{$id}{1}{gene} !~ m/&/ && $svpairs{$id}{1}{span} < 9999)) { if ($filter =~ m/LOWMAPQ|LowQual/i) { $filter = 'FailedQC'.$filter; } diff --git a/variants/svcalling.sh b/variants/svcalling.sh index e3f682cc252237e0aa2d1999c162aad407430771..41671c513a5ce7ea4ed7555a37518d9b5ab02bfc 100755 --- a/variants/svcalling.sh +++ b/variants/svcalling.sh @@ -123,7 +123,7 @@ then if [[ $filter == 1 ]] then zgrep '#CHROM' ${pair_id}.svaba.sv.vcf.gz > ${pair_id}.svaba.genefusion.txt - zcat ${pair_id}.svaba.sv.vcf.gz | $SNPEFF_HOME/scripts/vcfEffOnePerLine.pl |java -jar $SNPEFF_HOME/SnpSift.jar extractFields - CHROM POS ALT ID ANN[*].EFFECT ANN[*].GENE ANN[*].BIOTYPE FILTER FORMAT GEN[*] |grep -E 'gene_fusion|feature_fusion' | sort -u >> ${pair_id}.sgf.txt + zcat ${pair_id}.svaba.sv.vcf.gz | $SNPEFF_HOME/scripts/vcfEffOnePerLine.pl |java -jar $SNPEFF_HOME/SnpSift.jar extractFields - CHROM POS ALT ID ANN[*].EFFECT ANN[*].GENE ANN[*].BIOTYPE FILTER FORMAT GEN[*] |grep -E 'gene_fusion|feature_fusion' | sort -u > ${pair_id}.sgf.txt mv ${pair_id}.svaba.vcf.gz ${pair_id}.svaba.ori.vcf.gz perl $baseDir/filter_svaba.pl -t $tid -p $pair_id -i ${pair_id}.svaba.ori.vcf.gz -s ${pair_id}.svaba.sv.vcf.gz bgzip ${pair_id}.svaba.vcf diff --git a/variants/uniform_vcf_gt.pl b/variants/uniform_vcf_gt.pl index bef68c37b0c38a757dd087ae0d07e871e8b3fd76..661f7f21008897d707c753a3f9360eeacfcca108 100755 --- a/variants/uniform_vcf_gt.pl +++ b/variants/uniform_vcf_gt.pl @@ -27,15 +27,19 @@ while (my $line = <VCF>) { my ($key,$val) = split(/=/,$a); $hash{$key} = $val; } - if ($alt =~ m/^chr/) { - $chr2 =~ m/(chr\w+):(\d+)/; - $chr2=$1; + if ($alt =~ m/chr(\w+):(\d+)/i) { + $chr2='chr'.$1; $p2 = $2; $hash{CHR2} = $chr2; $hash{'END'} = $p2; $annot .= ";CHR2=$chr2;END=$p2"; + }elsif ($alt =~ m/CHR(\w+):(\d+)/i) { + $chr2='chr'.$1; + $p2 = $2; + $hash{CHR2} = 'chr'.$1; + $hash{END} = $2; + $annot .= ";CHR2=$chr2;END=$p2"; } - my @deschead = split(/:/,$format); my $newformat = 'GT:DP:AD:AO:RO'; my @newgts = (); @@ -61,6 +65,10 @@ while (my $line = <VCF>) { } elsif ($gtdata{AD} =~ m/^\d+$/){ $gtdata{AO} = $gtdata{AD}; $gtdata{RO} = $gtdata{DP} - $gtdata{AO}; + if ($gtdata{RO} < 0) { + $gtdata{DP} += $gtdata{AO}; + $gtdata{RO} = $gtdata{DP} - $gtdata{AO}; + } $gtdata{AD} = join(',',$gtdata{RO},$gtdata{AO}); } elsif (exists $gtdata{DV} && exists $gtdata{RV}) { $gtdata{AO} = $gtdata{DV} + $gtdata{RV};