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

Annotation of /FigKernelScripts/find_essential_reactions.pl

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : formsma 1.1 # _*_ Perl _*_
2 :     #
3 :     #
4 :     #
5 :     use strict;
6 :     use FIGV;
7 :     use Scenario;
8 :    
9 :     my ($genome_id,$orgdir) = @ARGV;
10 :    
11 :     my $debug = 0;
12 :    
13 :     #We want to basically generate a bunch of MPS files
14 :     # that are missing one reaction and determine what ones are nessisary for flux.
15 :     #
16 :    
17 :    
18 :     unless (scalar (@ARGV) >= 1)
19 :     {
20 :     print STDERR "usage: find_essential_reactions <genome_id>\n";
21 :     exit(-1);
22 :     }
23 :     my $fig;
24 :     if($orgdir)
25 :     {
26 :     $fig = new FIGV($orgdir);
27 :     }
28 :     else
29 :     {
30 :     $fig = new FIGV;
31 :     }
32 :     chomp $genome_id;
33 :    
34 :     my $filebase = $fig->model_directory($genome_id)."/Analysis/";
35 :    
36 :     my %reactions_tested;
37 :     my %flux_values;
38 :     while(1)
39 :     {
40 :     #open file for reading and writing.
41 :     open(NEW,">".$filebase."test.mps");
42 :     open(MPS,"<".$filebase."$genome_id.mps");
43 :     my $current_reaction = "";
44 :     while(<MPS>)
45 :     {
46 :     chomp;
47 :     if(!($_ =~ /(R\d{5}[\ |_r])/) || defined $reactions_tested{$1}
48 :     || ($1 ne $current_reaction && $current_reaction ne ""))
49 :     {
50 :     print NEW $_."\n";
51 :     }
52 :     elsif($current_reaction eq "" && defined $1)
53 :     {
54 :     print "Mached $1 as reaction\n";
55 :     $current_reaction = $1;
56 :     }
57 :     }
58 :     close(MPS);
59 :     close(NEW);
60 :     if($current_reaction eq "") {last;}
61 :     print "Testing reaction $current_reaction\n";
62 :     $reactions_tested{$current_reaction} = 1;
63 :     system(("/vol/biotools/bin/glpsol","--max",$filebase."test.mps","-o",$filebase."test.output"));
64 :    
65 :     #determine if we have flux!
66 :     my $flux;
67 :     open(OUT,$filebase."test.output");
68 :     while(<FLUX>)
69 :     {
70 :     if($_ =~ /^Objective/)
71 :     {
72 :     my @temp = split " ",$_;
73 :     chomp @temp;
74 :     if($temp[3] > 0)
75 :     {
76 :     $flux_values{$current_reaction} = $temp[3];
77 :     }
78 :     }
79 :     }
80 :     close(OUT);
81 :    
82 :     }
83 :     open(OUT,">".$filebase."essential_reactions.txt");
84 :     foreach (sort keys %flux_values)
85 :     {
86 :     if($flux_values{$_} == 0)
87 :     {
88 :     print OUT "$_ needed for biomass\n";
89 :     }
90 :     }
91 :     close(OUT);
92 :     system(("rm",$filebase."test.mps"));

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3