[Bio] / FigKernelScripts / update_lineage.pl Repository:
ViewVC logotype

Annotation of /FigKernelScripts/update_lineage.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (view) (download) (as text)

1 : wilke 1.1 use FIG;
2 :     use strict;
3 :     use warnings;
4 :    
5 : wilke 1.3 use vars qw( $opt_d $opt_g $opt_t $opt_f $opt_s $opt_u $opt_v);
6 : wilke 1.1 use Getopt::Std;
7 :    
8 :    
9 : wilke 1.3 getopts('d:go:f:s:uvt:');
10 : wilke 1.1
11 : wilke 1.2 my $success;
12 : wilke 1.3 my $usage = "get_lineage (-g | -f FILE ) [-s skip_n_first_lines] -d DIR\n";
13 : wilke 1.1
14 : wilke 1.2 unless ($opt_d and -d $opt_d) {
15 :     print $usage;
16 :     exit;
17 :     }
18 :    
19 :     my $org_dir = $opt_d;
20 : wilke 1.1
21 : wilke 1.3 my @list; # organism list
22 :    
23 : wilke 1.1 my $fig = new FIG;
24 :    
25 :     my $counter = 0;
26 :    
27 : wilke 1.3 if ( $opt_g ){
28 :    
29 :     my @genomes;
30 : wilke 1.2
31 : wilke 1.3 if ($opt_t){
32 :     push @genomes , $opt_t ;
33 :     }
34 :     else{
35 :     @genomes = $fig->genomes;
36 :     }
37 :    
38 :     print scalar @genomes , " genomes are currently in SEED\n";
39 :    
40 :     foreach my $genome (@genomes) {
41 :    
42 :     my $overview=request_ncbi( $org_dir , $genome);
43 :    
44 :     }
45 :     }
46 : wilke 1.1
47 : wilke 1.3 if ($opt_f and -f $opt_f ){
48 :    
49 :     open(FILE , $opt_f) or die "Can't open $opt_f for reading!\n";
50 :    
51 :     if ( $opt_s ){
52 :     my $skip = $opt_s;
53 :     while ( $skip ){
54 :     my $line = <FILE>;
55 :     $skip-- ;
56 :     }
57 :     }
58 :    
59 :    
60 :     while (my $line = <FILE>){
61 :     chomp $line;
62 :     my @fields = split "\t" , $line;
63 :    
64 :     my $overview;
65 :     if (scalar @fields == 3){
66 :     $overview = {
67 :     taxonomy_id => $fields[1],
68 :     seed_id => $fields[0],
69 :     lineage => $fields[2],
70 :     };
71 :     }
72 :     elsif (scalar @fields == 5) {
73 :     $overview = {
74 :     taxonomy_id => $fields[1],
75 :     seed_id => $fields[0],
76 :     lineage => $fields[2],
77 :     seed_org => $fields[3],
78 :     ncbi_org => $fields[4],
79 :     };
80 :     }
81 :     else{
82 :     # print STDERR "Wrong table format,\n";
83 :     # print STDERR scalar @fields , " fields in file,\n";
84 :     # print STDERR "either 3 or 5 columns required\n";
85 :     # print STDERR $line , "\n";
86 :     my $tax_id = get_taxonomy_id($org_dir , $fields[0]);
87 :    
88 :     $overview = $fig->get_organism_info_from_ncbi($tax_id);
89 :     #$overview->{ taxonomy_id } = $fields[1] unless $overview->{ taxonomy_id };
90 :     $overview->{ seed_org } = $fig->genus_species($fields[0] ) ;
91 :     $overview->{ ncbi_org } = $overview->{ scientific_name };
92 :     unless ( $overview->{ lineage } ){
93 :     print STDERR "Still a problem with\t" , $line , "\n";
94 :     next;
95 :     }
96 :     }
97 :     push @list , $overview;
98 : wilke 1.1 }
99 :    
100 : wilke 1.3 foreach my $entry ( @list ){
101 :     # my $overview = request_ncbi( $org_dir , $entry->{ seed_id } );
102 :     }
103 :    
104 :     }
105 :    
106 :     if ( $opt_v ){
107 :     _statistic( $fig , $org_dir , \@list )
108 :     }
109 :    
110 :    
111 :     if ($opt_u){
112 :    
113 :     foreach my $entry ( @list ){
114 :     update_organism_dir( $fig , $org_dir , $entry->{ seed_id } , $entry);
115 :     }
116 :    
117 :     }
118 :    
119 :     sub _statistic{
120 :     my ( $fig , $dir , $org , $list) = @_;
121 :     my $statistic;
122 : wilke 1.2
123 : wilke 1.3 my $missing;
124 :     foreach my $entry (@$list){
125 :    
126 :    
127 :    
128 :     # get seed words
129 :    
130 :     my @seed_words = split " " , $entry->{ seed_org };
131 :    
132 :     foreach my $word ( split " " , $entry->{ ncbi_org } ) {
133 :     print join " " , @seed_words , "\n";
134 :     for ( my $i=0 ; $i< scalar @seed_words; $i++ ) {
135 :     delete $seed_words[$i] if ( $word eq $seed_words[$i])
136 :     }
137 :     print join " " , @seed_words , "\n";
138 :     }
139 :    
140 :    
141 :    
142 :     }
143 : wilke 1.2
144 : wilke 1.1 }
145 :    
146 :     sub update_organism_dir{
147 :     my ( $fig , $dir , $org , $overview) = @_;
148 :    
149 :     return 0 unless (-d "$dir/$org");
150 : wilke 1.3
151 :     return 0 unless ( $org );
152 :     # update TAXONOMY file
153 :     # 1) write to log
154 :     open(LOG ,">>$dir/$org/.change_log" ) || die "could not open $dir/$org/.change_log";
155 :     open(TMP,"$dir/$org/TAXONOMY") || die "could not open $dir/$org/TAXONOMY";
156 :    
157 :     my $line = <TMP>;
158 :     chomp $line;
159 :     print LOG $line , "\t" , time , "\n";
160 :    
161 :     close(TMP);
162 :     close(LOG);
163 : wilke 1.1
164 : wilke 1.3 # 2) update entries
165 :     open(TMP,">$dir/$org/TAXONOMY") || die "could not open $dir/$org/TAXONOMY";
166 :     print TMP $overview->{lineage}."; ".$fig->genus_species($org)."\n";
167 :     close(TMP);
168 : wilke 1.1
169 :     open(TMP,">$dir/$org/TAXONOMY_ID") || die "could not open $dir/$org/TAXONOMY";
170 :     print TMP $overview->{taxonomy_id}."\n";
171 :     close(TMP);
172 :    
173 :     chmod(0777,"$dir/$org/TAXONOMY_ID");
174 :    
175 :     #update database entry
176 :    
177 :     my $old_lineage = $fig->taxonomy_of($org);
178 :    
179 :     return 1;
180 :     }
181 : wilke 1.3
182 :    
183 :     sub get_taxonomy_id{
184 :     my ($dir , $org) = @_;
185 :    
186 :     my $file = "$dir/$org/TAXONOMY_ID";
187 :    
188 :     my $tax_id;
189 :    
190 :     if (-f $file){
191 :     open(TAX, $file) || die "could not open $file";
192 :     $tax_id = <TAX>;
193 :     chomp $tax_id;
194 :     close(TAX);
195 :     }
196 :     else{
197 :     ( $tax_id ) = $org =~ /(\d+)\.\d+/ ;
198 :     unless ( $tax_id ){
199 :     print STDERR "No taxonomy id for $org\n";
200 :     }
201 :     }
202 :    
203 :     return $tax_id;
204 :     }
205 :    
206 :    
207 :     sub request_ncbi{
208 :     my ($org_dir , $genome) = @_;
209 :    
210 :     my $tax_id = get_taxonomy_id($org_dir , $genome);
211 :     #print STDERR "Tax id is $tax_id\n";
212 :     my ( $overview ) = $fig->get_organism_info_from_ncbi($tax_id);
213 :     my $seed_lineage = $overview->{lineage} . "; ".$fig->genus_species($genome);
214 :    
215 :     print $genome , "\t";
216 :     print "$tax_id\t$seed_lineage" ;
217 :     if ($fig->genus_species($genome) eq $overview->{scientific_name} ){
218 :     print "\t\t\n" ;
219 :     }
220 :     else{
221 :     print "\t".$fig->genus_species($genome)."\t".$overview->{scientific_name}."\n" ;
222 :     print STDERR $genome , "\t$tax_id\t$seed_lineage" ;
223 :     print STDERR "\t".$fig->genus_species($genome)."\t".$overview->{scientific_name}."\n" ;
224 :     }
225 :    
226 :     #$success= update_organism_dir( $fig , $dir, $genome , $overview);
227 :     #{ exit } unless $success;
228 :    
229 :     $counter ++;
230 :     #exit if ($counter > 10);
231 :     sleep 5 ;
232 :    
233 :     return $overview;
234 :     }

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3