Mercurial > repos > dereeper > sniplay
annotate hapmap2mlmm/HapmapToMLMMFiles.pl @ 15:31c23d943c29 draft default tip
Uploaded
author | dereeper |
---|---|
date | Tue, 08 Jan 2019 08:47:56 -0500 |
parents | c6640c49fd01 |
children |
rev | line source |
---|---|
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 | |
10
c6640c49fd01
planemo upload commit 475f4d7d8442a0d75e103af326ae5881c4d2a4ac
dereeper
parents:
9
diff
changeset
|
101 system("gawk -E $TRANSPOSE_EXE geno_transposed >geno_transposed2"); |
9 | 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 |