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

Diff of /FigKernelScripts/load_features.pl

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.11, Mon Aug 16 22:10:55 2004 UTC revision 1.20, Mon Dec 5 18:56:37 2005 UTC
# Line 1  Line 1 
1  # -*- perl -*-  # -*- perl -*-
2    #
3    # Copyright (c) 2003-2006 University of Chicago and Fellowship
4    # for Interpretations of Genomes. All Rights Reserved.
5    #
6    # This file is part of the SEED Toolkit.
7    #
8    # The SEED Toolkit is free software. You can redistribute
9    # it and/or modify it under the terms of the SEED Toolkit
10    # Public License.
11    #
12    # You should have received a copy of the SEED Toolkit Public License
13    # along with this program; if not write to the University of Chicago
14    # at info@ci.uchicago.edu or the Fellowship for Interpretation of
15    # Genomes at veronika@thefig.info or download a copy from
16    # http://www.theseed.org/LICENSE.TXT.
17    #
18    
19    
20  use strict;  use strict;
21  use FIG;  use FIG;
22  my $fig = new FIG;  my $fig = new FIG;
23    
24  use DBrtns;  use Tracer;
   
25    
26    Trace("Preparing to load features.") if T(2);
27    my ($mode, @genomes) = FIG::parse_genome_args(@ARGV);
28  my $temp_dir = "$FIG_Config::temp";  my $temp_dir = "$FIG_Config::temp";
29  my($organisms_dir) = "$FIG_Config::organisms";  my $organisms_dir = "$FIG_Config::organisms";
30    
31  my($genome,@types,$type,$id,$loc,@aliases,$aliases,$contig);  my($genome,@types,$type,$id,$loc,@aliases,$aliases,$contig);
32    
33  # usage: load_features [G1 G2 G3 ... ]  # usage: load_features [G1 G2 G3 ... ]
34    
35  open(REL,">$temp_dir/tmpfeat$$") || die "could not open $temp_dir/tmpfeat$$";  Open(\*REL, ">$temp_dir/tmpfeat$$");
36  open(ALIAS,"| sort -u > $temp_dir/tmpalias$$") || die "could not open $temp_dir/tmpalias$$";  Open(\*ALIAS, "| sort -T $temp_dir -u > $temp_dir/tmpalias$$");
   
 my $dbf = $fig->{_dbf};  
   
   
   
 my @genomes;  
 if (@ARGV == 0)  
 {  
     $dbf->drop_table( tbl => "features" );  
     $dbf->drop_table( tbl => "ext_alias" );  
   
     $dbf->create_table( tbl => 'ext_alias',  
                         flds => "id varchar(32), alias varchar(32), genome varchar(16)"  
                         );  
37    
     if ($FIG_Config::dbms eq "Pg")  
     {  
         $dbf->create_table( tbl  => "features",  
                             flds => "id varchar(32), idN INTEGER, type varchar(16),genome varchar(16),"  .  
                                     "location varchar(5000),"  .  
                                     "contig varchar(96), minloc INTEGER, maxloc INTEGER,"  .  
                                     "aliases TEXT"  
                             );  
     }  
     elsif ($FIG_Config::dbms eq "mysql")  
     {  
         $dbf->create_table( tbl  => "features",  
                             flds => "id varchar(32), idN INTEGER, type varchar(16),genome varchar(16),"  .  
                                     "location TEXT,"  .  
                                     "contig varchar(96), minloc INTEGER, maxloc INTEGER,"  .  
                                     "aliases TEXT"  
                             );  
     }  
38    
39      @genomes = $fig->genomes;  if ($mode eq 'all') {
40    
41      # Here we extract external aliases from the peg.synonyms table, when they can be inferred      # Here we extract external aliases from the peg.synonyms table, when they can be inferred
42      # accurately.      # accurately.
43      open(SYN,"<$FIG_Config::global/peg.synonyms") || die "could not open $FIG_Config::global/peg.synonyms";          Trace("Extracting external aliases from the peg.synonyms table.") if T(2);
44        open(\*SYN, "<$FIG_Config::global/peg.synonyms");
45      while (defined($_ = <SYN>))      while (defined($_ = <SYN>))
46      {      {
47          chop;          chop;
48          my($x,$y) = split(/\t/,$_);          my($x,$y) = split(/\t/,$_);
49          my @ids = map { $_ =~ /^([^,]+)/ } ($x,split(/;/,$y));                  my @ids = map { $_ =~ /^([^,]+),(\d+)/; [$1,$2] } ($x,split(/;/,$y));
50          my @fig = ();          my @fig = ();
51          my(@nonfig) = ();          my(@nonfig) = ();
52          foreach $_ (@ids)          foreach $_ (@ids)
53          {          {
54              if ($_ =~ /^fig\|/)                          if ($_->[0] =~ /^fig\|/)
55              {              {
56                  push(@fig,$_);                  push(@fig,$_);
57              }              }
# Line 74  Line 61 
61              }              }
62          }          }
63    
64          if (@fig == 1)                  foreach $x (@fig)
65          {          {
66              my $genome = &FIG::genome_of($fig[0]);                          my($peg,$peg_ln) = @$x;
67                            my $genome = &FIG::genome_of($peg);
68              foreach $_ (@nonfig)              foreach $_ (@nonfig)
69              {              {
70                  print ALIAS "$fig[0]\t$_\t$genome\n";                                  if ((@fig == 1) || ($peg_ln == $_->[1]))
71              }                                  {
72                                            print ALIAS "$peg\t$_->[0]\t$genome\n";
73                                            Trace("Alias record $peg, $_->[0] for $genome.") if T(4);
74          }          }
75      }      }
     close(SYN);  
76  }  }
 else  
 {  
     @genomes = @ARGV;  
     foreach $genome (@genomes)  
     {  
         $dbf->SQL("DELETE FROM features WHERE ( genome = \'$genome\' )");  
         $dbf->SQL("DELETE FROM ext_alias WHERE ( genome = \'$genome\' )");  
77      }      }
78        close(SYN);
79  }  }
80    
81  foreach $genome (@genomes)  foreach $genome (@genomes)
82  {  {
83            Trace("Processing $genome.") if T(3);
84      opendir(FEAT,"$organisms_dir/$genome/Features")      opendir(FEAT,"$organisms_dir/$genome/Features")
85          || die "could not open $genome/Features";          || die "could not open $genome/Features";
86      @types = grep { $_ =~ /^[a-zA-Z]+$/ } readdir(FEAT);      @types = grep { $_ =~ /^[a-zA-Z]+$/ } readdir(FEAT);
# Line 107  Line 91 
91          if ((-s "$organisms_dir/$genome/Features/$type/tbl") &&          if ((-s "$organisms_dir/$genome/Features/$type/tbl") &&
92              open(TBL,"<$organisms_dir/$genome/Features/$type/tbl"))              open(TBL,"<$organisms_dir/$genome/Features/$type/tbl"))
93          {          {
94  #           print STDERR "loading $genome/Features/$type/tbl\n";              Trace("Loading $genome/Features/$type/tbl") if T(4);
95              while (defined($_ = <TBL>))              my @tbl = <TBL>;
96                close(TBL);
97                my %seen;
98    
99                while ($_ = pop @tbl)
100              {              {
101                  chop;                  chop;
102                  ($id,$loc,@aliases) = split(/\t/,$_);                  ($id,$loc,@aliases) = split(/\t/,$_);
103                  if ($id)  
104                    if ($id && (! $seen{$id}))
105                  {                  {
106                        $seen{$id} = 1;
107                      my($minloc,$maxloc);                      my($minloc,$maxloc);
108                      if ($loc)                      if ($loc)
109                      {                      {
# Line 136  Line 126 
126                          my $alias;                          my $alias;
127                          foreach $alias (@aliases)                          foreach $alias (@aliases)
128                          {                          {
129                              if ($alias =~ /^(NP_|gi\||sp\|\tr\||kegg\||uni\|)/)                              if ($alias =~ /^([NXYZA]P_|gi\||sp\|\tr\||kegg\||uni\|)/)
130                              {                              {
131    
132                                  print ALIAS "$id\t$alias\t$genome\n";                                  print ALIAS "$id\t$alias\t$genome\tOVERRIDE\n";
133                                    Trace("$id override alias $alias for $genome") if T(4);
134                              }                              }
135                          }                          }
136                      }                      }
# Line 155  Line 146 
146                      }                      }
147                  }                  }
148              }              }
             close(TBL);  
149          }          }
150      }      }
151  }  }
152  close(REL);  close(REL);
153  close(ALIAS);  close(ALIAS);
154    Open(\*ALIASIN, "<$temp_dir/tmpalias$$");
155    Open(\*ALIASOUT, ">$temp_dir/tmpalias$$.1");
156    Trace("Parsing alias file.") if T(2);
157    $_ = <ALIASIN>;
158    while ($_ && ($_ =~ /^(\S+)/))
159    {
160        my @aliases = ();
161        my $curr = $1;
162        while ($_ && ($_ =~ /^(\S+)\t(\S+)(\t(\S+))?/) && ($1 eq $curr))
163        {
164                    push(@aliases,[$2,$3 ? 1 : 0]);
165                    $_ = <ALIASIN>;
166        }
167        my $x;
168        my $genome = &FIG::genome_of($curr);
169        foreach $x (@aliases)
170        {
171            if ($x->[1])
172            {
173                print ALIASOUT "$curr\t$x->[0]\t$genome\n";
174            }
175            else
176            {
177                my $i;
178                for ($i=0; ($i < @aliases) && ((! $aliases[$i]->[1]) || (! &same_class($x->[0],$aliases[$i]->[0]))); $i++) {}
179                if ($i == @aliases)
180                {
181                    print ALIASOUT "$curr\t$x->[0]\t$genome\n";
182                }
183            }
184        }
185    }
186    close(ALIASIN);
187    close(ALIASOUT);
188    unlink("$temp_dir/tmpalias$$");
189    
190  $dbf->load_table( tbl => "features",  $fig->reload_table($mode, 'features',
191                    file => "$temp_dir/tmpfeat$$" );                                     "id varchar(32), idN INTEGER, type varchar(16),genome varchar(16),"  .
192                                        "location TEXT,"  .
193                                        "contig varchar(96), minloc INTEGER, maxloc INTEGER,"  .
194                                        "aliases TEXT",
195                                            { features_id_ix => "id", features_org_ix => "genome",
196                                              features_type_ix => "type", features_beg_ix => "genome, contig, minloc" },
197                                            "$temp_dir/tmpfeat$$", \@genomes);
198    unlink("$temp_dir/tmpfeat$$");
199    
200  $dbf->load_table( tbl => "ext_alias",  $fig->reload_table($mode, 'ext_alias',
201                    file => "$temp_dir/tmpalias$$" );                                          "id varchar(32), alias varchar(32), genome varchar(16)",
202                                            { ext_alias_alias_ix => "alias", ext_alias_genome_ix => "genome",
203                                              ext_alias_ix_id => "id" },
204                                            "$temp_dir/tmpalias$$.1", \@genomes );
205    
206  if (@ARGV == 0)  unlink("$temp_dir/tmpalias$$.1");
207  {  Trace("Features loaded.") if T(2);
     $dbf->create_index( idx  => "ext_alias_alias_ix",  
                         tbl  => "ext_alias",  
                         type => "btree",  
                         flds => "alias" );  
   
     $dbf->create_index( idx  => "ext_alias_genome_ix",  
                         tbl  => "ext_alias",  
                         type => "btree",  
                         flds => "genome" );  
   
     $dbf->create_index( idx  => "ext_alias_id_ix",  
                         tbl  => "ext_alias",  
                         type => "btree",  
                         flds => "id" );  
   
     $dbf->create_index( idx  => "features_id_ix",  
                         tbl  => "features",  
                         type => "btree",  
                         flds => "id" );  
     $dbf->create_index( idx  => "features_org_ix",  
                         tbl  => "features",  
                         type => "btree",  
                         flds => "genome" );  
     $dbf->create_index( idx  => "features_type_ix",  
                         type => "btree",  
                         tbl  => "features",  
                         flds => "type" );  
     $dbf->create_index( idx  => "features_beg_ix",  
                         type => "btree",  
                         tbl  => "features",  
                         flds => "genome,contig,minloc" );  
208    
209      $dbf->vacuum_it("features")  sub same_class {
210        my($x,$y) = @_;
211    
212        my $class1 = &classA($x);
213        my $class2 = &classA($y);
214        return ($class1 && ($class1 eq $class2));
215    }
216    
217    sub classA {
218        my($alias) = @_;
219    
220        if ($alias =~ /^([^\|]+)\|/)
221        {
222                    return $1;
223        }
224        elsif ($alias =~ /^[NXYZA]P_[0-9\.]+$/)
225        {
226                    return "refseq";
227        }
228        else
229        {
230                    return "";
231        }
232  }  }
 unlink("$temp_dir/tmpfeat$$");  
 unlink("$temp_dir/tmpalias$$");  

Legend:
Removed from v.1.11  
changed lines
  Added in v.1.20

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3