-
Notifications
You must be signed in to change notification settings - Fork 0
/
transform_to_0and1_gatk.pl
executable file
·54 lines (35 loc) · 1.06 KB
/
transform_to_0and1_gatk.pl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
#! /usr/bin/perl -w
use warnings;
###transform 1/1, 0/1, 1/0 and 0/0 to 1, hz and 0.
my $input= $ARGV[0]; #e.g. "REALIGNMENT_PHAGE.raw1.gatk.vcf
my @temp; # array to store 0s and 1s
my @temp2; # array to store the original data so that the (GT:PL:DP:SP:GQ) data is retained
open(LIST,"<$input") || die "cannot find list file!\n";
while(<LIST>){
chomp ($_);
#To print the header
if ($_=~/#CHROM/){
print "$_\n";
}
#if ($_=~/^NC_002516.2||RGP42||PNUK73/){ #old
if ($_!~/^#/){#match any line that isn't a comment
@temp=split('\t',$_);
@temp2=@temp;
for ($i=9; $i < scalar(@temp); $i++){
if ($temp[$i]=~/^1:.*/){
$temp[$i] = 1;
}
if ($temp[$i]=~/^0:.*/){
$temp[$i] = 0;
}
if ($temp[$i]=~/^2:.*/){ # gatk tools puts a 2 when there are 3 different types of amino acid
$temp[$i] = 2;
}
if ($temp[$i]=~/\./){ # gatk puts a "." when it's not capable of making a call. kept just in case.
#gatk posa . quan no es capac de cridar en algunes mostres i en altres si. me'ls quedo per si de cas.
$temp[$i] = 999;
}
}
print join("\t", @temp),"\n";
}
}