9
|
1 #!/usr/bin/perl
|
|
2
|
|
3 use strict;
|
|
4
|
|
5 use Getopt::Long;
|
|
6
|
|
7 my $usage = qq~Usage:$0 <args> [<opts>]
|
|
8 where <args> are:
|
|
9 -h, --hapmap <Hapmap input file>
|
|
10 -m, --map <Map output file>
|
|
11 -g, --geno <Genotype output file>
|
|
12 -p, --path <Path for transpose executable>
|
|
13 ~;
|
|
14 $usage .= "\n";
|
|
15
|
|
16 my ($hapmap,$map,$geno,$path);
|
|
17
|
|
18
|
|
19 GetOptions(
|
|
20 "geno=s" => \$geno,
|
|
21 "map=s" => \$map,
|
|
22 "hapmap=s" => \$hapmap,
|
|
23 "path=s" => \$path,
|
|
24 );
|
|
25
|
|
26
|
|
27 die $usage
|
|
28 if ( !$geno || !$map || !$hapmap || !$path);
|
|
29
|
|
30 my $TRANSPOSE_EXE = "$path/transpose.awk";
|
|
31
|
|
32 my @snps;
|
|
33 my %chrom_pos;
|
|
34 my $num_line = 0;
|
|
35 open(my $O,">geno_transposed");
|
|
36 open(my $H,$hapmap);
|
|
37 while(<$H>)
|
|
38 {
|
|
39 $num_line++;
|
|
40 my $line = $_;
|
|
41 chomp($line);
|
|
42 $line =~s/\r//g;
|
|
43 $line =~s/\n//g;
|
|
44 my @infos = split(/\t/,$line);
|
|
45 if ($num_line == 1)
|
|
46 {
|
|
47 print $O "Ind_id";
|
|
48 for (my $i = 11; $i <= $#infos; $i++)
|
|
49 {
|
|
50 my $individual = $infos[$i];
|
|
51 print $O " " . $individual;
|
|
52 }
|
|
53 print $O "\n";
|
|
54 }
|
|
55 elsif ($num_line > 1)
|
|
56 {
|
|
57 my $snp = $infos[0];
|
|
58 my $variation = $infos[1];
|
|
59 my %scores;
|
|
60 if ($variation =~/(\w)\/(\w)/)
|
|
61 {
|
|
62 my $allele1 = $1;
|
|
63 my $allele2 = $2;
|
|
64 $scores{$allele1} = 0;
|
|
65 $scores{$allele2} = 1;
|
|
66 }
|
|
67 my $chrom = $infos[2];
|
|
68 my $pos = $infos[3];
|
|
69 $chrom_pos{$snp}{"chrom"} = $chrom;
|
|
70 $chrom_pos{$snp}{"pos"} = $pos;
|
|
71 push(@snps,$snp);
|
|
72 print $O "$snp";
|
|
73 for (my $i = 11; $i <= $#infos; $i++)
|
|
74 {
|
|
75 my $genotype = $infos[$i];
|
|
76 my @alleles = split("",$genotype);
|
|
77 if ($genotype ne "NN")
|
|
78 {
|
|
79 my $score = $scores{$alleles[0]} + $scores{$alleles[1]};
|
|
80 print $O " $score";
|
|
81 }
|
|
82 else
|
|
83 {
|
|
84 print $O " NA";
|
|
85 }
|
|
86 }
|
|
87 print $O "\n";
|
|
88 }
|
|
89 }
|
|
90 close($H);
|
|
91 close($O);
|
|
92
|
|
93 open(my $M,">$map");
|
|
94 print $M "SNP Chr Pos\n";
|
|
95 foreach my $snp(@snps)
|
|
96 {
|
|
97 print $M "$snp " . $chrom_pos{$snp}{"chrom"} . " ". $chrom_pos{$snp}{"pos"} . "\n";
|
|
98 }
|
|
99 close($M);
|
|
100
|
|
101 system("$TRANSPOSE_EXE geno_transposed >geno_transposed2");
|
|
102
|
|
103 open(my $F,">$geno");
|
|
104 open(my $G,"geno_transposed2");
|
|
105 while(<$G>)
|
|
106 {
|
|
107 my $line = $_;
|
|
108 $line =~s/ /\t/g;
|
|
109 print $F $line;
|
|
110 }
|
|
111 close($G);
|
|
112 close($F);
|
|
113
|
|
114 unlink("geno_transposed");
|
|
115 unlink("geno_transposed2");
|
|
116
|
|
117
|