[Bio] / FigKernelPackages / ExpressionDir.pm Repository:
ViewVC logotype

Annotation of /FigKernelPackages/ExpressionDir.pm

Parent Directory Parent Directory | Revision Log Revision Log


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

1 : olson 1.1 package ExpressionDir;
2 :    
3 :     use Data::Dumper;
4 :     use strict;
5 :     use SeedAware;
6 :     use File::Copy;
7 : olson 1.2 use File::Temp 'tempfile';
8 : olson 1.1 use File::Spec::Functions;
9 :     use base 'Class::Accessor';
10 :     use Carp;
11 :     use Fcntl ':seek';
12 :    
13 : olson 1.4 __PACKAGE__->mk_accessors(qw(genome_dir expr_dir genome_id));
14 : olson 1.1
15 : olson 1.5 our @probe_parsers = qw(parse_probe_format_1lq
16 :     parse_probe_format_1
17 :     parse_probe_format_2
18 :     parse_probe_format_3);
19 :    
20 : olson 1.1 =head3 new
21 :    
22 : olson 1.4 my $edir = ExpressionDir->new($expr_dir);
23 :    
24 :     Create a new ExpressionDir object from an existing expression dir.
25 : olson 1.1
26 :    
27 :     =cut
28 :    
29 :     sub new
30 :     {
31 : olson 1.4 my($class, $expr_dir) = @_;
32 :    
33 :     my $gfile = catfile($expr_dir, "GENOME_ID");
34 :     open(GF, "<", $gfile) or die "Cannot open $expr_dir/GENOME_ID: $!";
35 :     my $genome_id = <GF>;
36 :     chomp $genome_id;
37 :     close(GF);
38 :    
39 :     my $self = {
40 :     genome_dir => catfile($expr_dir, $genome_id),
41 :     genome_id => $genome_id,
42 :     expr_dir => $expr_dir,
43 :     };
44 :     return bless $self, $class;
45 :     }
46 :    
47 :     =head3 create
48 : olson 1.1
49 : olson 1.4 my $edir = ExpressionDir->create($expr_dir, $genome_id, $genome_src)
50 :    
51 :     Create a new expression directory from the given genome id and genome data source.
52 :     The data source is either a FIG or FIGV object, or a SAPserver object
53 :     that points at a server from which the data can be extracted.
54 :    
55 :     =cut
56 :    
57 :     sub create
58 :     {
59 :     my($class, $expr_dir, $genome_id, $genome_src) = @_;
60 :    
61 :     if (! -d $expr_dir)
62 : olson 1.1 {
63 : olson 1.4 mkdir($expr_dir);
64 : olson 1.1 }
65 : olson 1.4 my $genome_dir = catfile($expr_dir, $genome_id);
66 :     if (! -d $genome_dir)
67 : olson 1.1 {
68 : olson 1.4 mkdir($genome_dir);
69 : olson 1.1 }
70 : olson 1.4
71 :     open(GF, ">", catfile($expr_dir, "GENOME_ID"));
72 :     print GF "$genome_id\n";
73 :     close(GF);
74 : olson 1.1
75 :     my $self = {
76 :     genome_dir => $genome_dir,
77 : olson 1.4 genome_id => $genome_id,
78 :     expr_dir => $expr_dir,
79 : olson 1.1 };
80 : olson 1.4 bless $self, $class;
81 :    
82 :     if (ref($genome_src) =~ /^FIG/)
83 :     {
84 :     $self->create_from_fig($genome_src);
85 :     }
86 :     elsif (ref($genome_src =~ /^SAP/))
87 :     {
88 :     $self->create_from_sap($genome_src);
89 :     }
90 :     else
91 :     {
92 :     confess "Unknown genome source\n";
93 :     }
94 :     return $self;
95 :     }
96 :    
97 :     sub create_from_fig
98 :     {
99 :     my($self, $fig) = @_;
100 :    
101 :     my $gdir = $fig->organism_directory($self->genome_id);
102 : olson 1.5
103 :     if (! -d $gdir)
104 :     {
105 :     confess "Genome directory $gdir not found";
106 :     }
107 :    
108 : olson 1.4 copy(catfile($gdir, "contigs"), catfile($self->genome_dir, "contigs"));
109 :     mkdir(catfile($self->genome_dir, "Features"));
110 :     my @pegs;
111 :     for my $ftype (qw(peg rna))
112 :     {
113 :     my $ofdir = catfile($gdir, "Features", $ftype);
114 :     my $nfdir = catfile($self->genome_dir, "Features", $ftype);
115 :    
116 :     if (open(OT, "<", "$ofdir/tbl"))
117 :     {
118 :     mkdir($nfdir);
119 :     open(NT, ">", "$nfdir/tbl") or confess "Cannot write $nfdir/tbl: $!";
120 :     while (<OT>)
121 :     {
122 :     my($id) = /^(\S+)\t/;
123 :     if (!$fig->is_deleted_fid($id))
124 :     {
125 :     print NT $_;
126 :     push(@pegs, $id) if $ftype eq 'peg';
127 :     }
128 :     }
129 :     close(OT);
130 :     close(NT);
131 :     }
132 :     copy("$ofdir/fasta", "$nfdir/fasta");
133 :     }
134 :    
135 :     open(SS, ">", catfile($self->genome_dir, "subsystem.data"));
136 :     my %phash = $fig->subsystems_for_pegs_complete(\@pegs);
137 :     for my $peg (keys %phash)
138 :     {
139 :     my $list = $phash{$peg};
140 :     next unless $list;
141 :    
142 :     for my $ent (@$list)
143 :     {
144 :     print SS join("\t", $peg, @$ent), "\n";
145 :     }
146 :     }
147 :     close(SS);
148 :    
149 :     open(AF, ">", catfile($self->genome_dir, "assigned_functions"));
150 :     my $fns = $fig->function_of_bulk(\@pegs);
151 :     for my $peg (@pegs)
152 :     {
153 :     print AF join("\t", $peg, $fns->{$peg}), "\n";
154 :     }
155 :     close(AF);
156 :     }
157 :    
158 :     sub create_from_sap
159 :     {
160 :     my($self, $sap) = @_;
161 :     confess "create_from_sap not yet implemented";
162 : olson 1.1 }
163 :    
164 : olson 1.5 sub parse_probe_format_1lq
165 :     {
166 :     my($self, $in_file, $out_file) = @_;
167 :    
168 :     my($fh);
169 :    
170 :     if ($in_file !~ /\.1lq$/)
171 :     {
172 :     return undef;
173 :     }
174 :    
175 :     open($fh, "<", $in_file) or confess "Cannot open $in_file for reading: $!";
176 :    
177 :     my $out;
178 :     open($out, ">", $out_file) or confess "Cannot open $out for writing: $!";
179 :    
180 :     # Skip 3 header lines.
181 :     $_ = <$fh> for 1..3;
182 :    
183 :     while (defined($_ = <$fh>))
184 :     {
185 :     if ($_ =~ /(\d+)\s+(\d+)\s+([ACGT]+)\s+(-?\d+)\s/)
186 :     {
187 :     if (length($3) < 15)
188 :     {
189 :     close($fh);
190 :     close($out);
191 :     confess "Bad length at line $. of $in_file";
192 :     return undef;
193 :     }
194 :     next if ($4 =~ /\d+3$/); #mismatch probe
195 :     my($x,$y,$seq) = ($1,$2,$3);
196 :     $seq = scalar reverse $seq;
197 :     print $out "$x\_$y\t$seq\n";
198 :     }
199 :     else
200 :     {
201 :     #
202 :     # We expect some lines not to match.
203 :     #
204 :     }
205 :     }
206 :    
207 :     close($fh);
208 :     close($out);
209 :     return 1;
210 :     }
211 :    
212 : olson 1.1 sub parse_probe_format_1
213 :     {
214 :     my($self, $in_file, $out_file) = @_;
215 :    
216 :     my($fh);
217 :    
218 :     open($fh, "<", $in_file) or confess "Cannot open $in_file for reading: $!";
219 :     my $l = <$fh>;
220 :     chomp $l;
221 :     $l =~ s/\r//;
222 :     my @hdrs = split(/\t/, $l);
223 :     my %hdrs;
224 :     $hdrs{$hdrs[$_]} = $_ for 0..$#hdrs;
225 :    
226 :     my $x_col = $hdrs{"Probe X"};
227 :     my $y_col = $hdrs{"Probe Y"};
228 :     my $seq_col = $hdrs{"Probe Sequence"};
229 :     if (!(defined($x_col) && defined($y_col) && defined($seq_col)))
230 :     {
231 :     close($fh);
232 :     return undef;
233 :     }
234 :    
235 :     my $out;
236 :     open($out, ">", $out_file) or confess "Cannot open $out for writing: $!";
237 :    
238 :     while (<$fh>)
239 :     {
240 :     chomp;
241 :     s/\r//g;
242 :     my @flds = split(/\t/,$_);
243 :     my($x,$y,$seq);
244 :     $x = $flds[$x_col];
245 :     $y = $flds[$y_col];
246 :     $seq = $flds[$seq_col];
247 :     my $id = "$x\_$y";
248 :     print $out "$id\t$seq\n";
249 :     }
250 :     close($fh);
251 :     close($out);
252 :     return 1;
253 :     }
254 :    
255 :     sub parse_probe_format_2
256 :     {
257 :     my($self, $in_file, $out_file) = @_;
258 :    
259 :     my($fh);
260 :    
261 :     local $/ = "\n>";
262 :    
263 :     open($fh, "<", $in_file) or confess "Cannot open $in_file for reading: $!";
264 :     my $l = <$fh>;
265 :     chomp $l;
266 :     $l =~ s/\r//;
267 :    
268 :     if ($l !~ /^>?\S+:(\d+):(\d+);\s+Interrogation_Position=\d+;\s+Antisense;\n([ACGT]+)/s)
269 :     {
270 :     close($fh);
271 :     return undef;
272 :     }
273 :     seek($fh, 0, SEEK_SET);
274 :    
275 :     my $out;
276 :     open($out, ">", $out_file) or confess "Cannot open $out for writing: $!";
277 :    
278 :     while (<$fh>)
279 :     {
280 :     chomp;
281 :    
282 :     if ($_ =~ /^>?\S+:(\d+):(\d+);\s+Interrogation_Position=\d+;\s+Antisense;\n([ACGT]+)/s)
283 :     {
284 :     if (length($3) < 15)
285 :     {
286 :     close($fh);
287 :     confess "Bad length at line $. of $in_file";
288 :     }
289 :     print $out "$1\_$2\t$3\n";
290 :     }
291 :     else
292 :     {
293 :     confess "Bad input at line $. of $in_file";
294 :     }
295 :     }
296 :     close($out);
297 :     close($fh);
298 :     return 1;
299 :     }
300 :    
301 :     #
302 :     # Our "native" format, used for passing through pre-parsed data.
303 :     #
304 :     sub parse_probe_format_3
305 :     {
306 :     my($self, $in_file, $out_file) = @_;
307 :    
308 :     my($fh);
309 :    
310 :     open($fh, "<", $in_file) or confess "Cannot open $in_file for reading: $!";
311 :     my $l = <$fh>;
312 :     chomp $l;
313 :     $l =~ s/\r//;
314 :    
315 :     if ($l !~ /^\d+_\d+\t[ACGT]+$/)
316 :     {
317 :     close($fh);
318 :     return undef;
319 :     }
320 :     seek($fh, 0, SEEK_SET);
321 :    
322 :     my $out;
323 :     open($out, ">", $out_file) or confess "Cannot open $out for writing: $!";
324 :    
325 :     while (<$fh>)
326 :     {
327 :     if ($_ =~ /^\d+_\d+\t[ACGT]+$/)
328 :     {
329 :     print $out $_;
330 :     }
331 :     else
332 :     {
333 :     confess "Bad input at line $. of $in_file";
334 :     }
335 :     }
336 :     close($out);
337 :     close($fh);
338 :     return 1;
339 :     }
340 :    
341 :     sub compute_probe_to_peg
342 :     {
343 :     my($self, $probes) = @_;
344 :    
345 : olson 1.5 my($probe_suffix) = $probes =~ /(\.[^.]+)$/;
346 :    
347 :     my $my_probes = catfile($self->expr_dir, "probes.in$probe_suffix");
348 : olson 1.1
349 :     copy($probes, $my_probes) or confess "Cannot copy $probes to $my_probes: $!";
350 :    
351 : olson 1.5 my $probes_fasta = catfile($self->expr_dir, "probes");
352 : olson 1.1
353 :     #
354 :     # Attempt to translate probe file.
355 :     #
356 :     my $success;
357 : olson 1.5 for my $meth (@probe_parsers)
358 : olson 1.1 {
359 :     if ($self->$meth($my_probes, $probes_fasta))
360 :     {
361 :     print STDERR "Translated $probes to $probes_fasta using $meth\n";
362 :     $success = 1;
363 :     last;
364 :     }
365 :     else
366 :     {
367 :     print STDERR "Failed to translate $probes to $probes_fasta using $meth\n";
368 :     }
369 :     }
370 :     if (!$success)
371 :     {
372 :     confess "Could not translate $probes\n";
373 :     }
374 :    
375 :     my $peg_probe_table = catfile($self->expr_dir, 'peg.probe.table');
376 :     my $probe_occ_table = catfile($self->expr_dir, 'probe.occ.table');
377 :    
378 : olson 1.2 my $feature_dir = catfile($self->genome_dir, "Features");
379 :     my @tbls;
380 :     for my $ftype (qw(peg rna))
381 :     {
382 :     my $tfile = catfile($feature_dir, $ftype, 'tbl');
383 :     if (-f $tfile)
384 :     {
385 :     push(@tbls, $tfile);
386 :     }
387 :     }
388 :     if (@tbls == 0)
389 :     {
390 :     confess "Could not find any tbl files in $feature_dir";
391 :     }
392 :    
393 : olson 1.4 $self->run([executable_for("make_probes_to_genes"),
394 :     $probes_fasta,
395 :     catfile($self->genome_dir, 'contigs'),
396 :     $tbls[0],
397 :     $peg_probe_table,
398 :     $probe_occ_table,
399 :     @tbls[1..$#tbls],
400 :     ],
401 :     { stderr => catfile($self->expr_dir, 'problems') });
402 : olson 1.1
403 :    
404 : olson 1.4 $self->run([executable_for("remove_multiple_occurring_probes"),
405 :     $peg_probe_table,
406 :     $probe_occ_table],
407 :     { stdout => catfile($self->expr_dir, 'peg.probe.table.no.multiple') } );
408 : olson 1.1
409 :     $self->make_missing_probes($peg_probe_table, $probes_fasta,
410 :     catfile($self->expr_dir, 'probe.no.match'));
411 :     }
412 : olson 1.2
413 : olson 1.1 sub make_missing_probes
414 :     {
415 :     my($self, $probe_table, $probes, $output) = @_;
416 :     open(MATCH,"<", $probe_table) or die "Cannot open $probe_table: $!";
417 :     open(PROBES,"<", $probes) or die "Cannot open $probes: $!";
418 :     open(OUTPUT, ">", $output) or die "Cannot open $output: $!";
419 :     my %locations;
420 :     while(<MATCH>)
421 :     {
422 :     chomp;
423 :     my($peg,$loc)=split "\t";
424 :     $locations{$loc} = $peg;
425 :     }
426 :    
427 :     while(<PROBES>)
428 :     {
429 :     chomp;
430 :     my($loc,$seq) = split "\t";
431 :     print OUTPUT $loc, "\n" if ! exists $locations{$loc};
432 :     }
433 :     close(MATCH);
434 :     close(PROBES);
435 :     close(OUTPUT);
436 :     }
437 :    
438 : olson 1.2 #
439 :     # we don't copy the experiment files in here because
440 :     # they may be very large. This may change.
441 :     #
442 :     # We do copy the cdf.
443 :     #
444 :     sub compute_rma_normalized
445 :     {
446 :     my($self, $cdf_file, $expt_dir) = @_;
447 :    
448 :     my $my_cdf = catfile($self->expr_dir, "expr.cdf");
449 :     copy($cdf_file, $my_cdf) or confess "Cannot copy $cdf_file to $my_cdf: $!";
450 :    
451 :     #
452 :     # We need to build the R library for this cdf.
453 :     #
454 :     my($fh, $tempfile) = tempfile();
455 :     #m = make.cdf.package("S_aureus.cdf", cdf.path="..",packagename="foo",package.path="/tmp")
456 :    
457 :     my $cdf_path = $self->expr_dir;
458 :     my $libdir = catfile($self->expr_dir, "r_lib");
459 :     -d $libdir or mkdir $libdir;
460 :     my $pkgdir = catfile($self->expr_dir, "r_pkg");
461 :     -d $pkgdir or mkdir $pkgdir;
462 :    
463 :     print $fh "library(makecdfenv);\n";
464 :     print $fh qq(make.cdf.package("expr.cdf", cdf.path="$cdf_path", packagename="datacdf", package.path="$pkgdir", species="genome name");\n);
465 :     close($fh);
466 :     system("Rscript", $tempfile);
467 :     system("R", "CMD", "INSTALL", "-l", $libdir, "$pkgdir/datacdf");
468 :    
469 :     local($ENV{R_LIBS}) = $libdir;
470 : olson 1.4 $self->run([executable_for("RunRMA"),
471 :     "data",
472 :     catfile($self->expr_dir, "peg.probe.table"),
473 :     catfile($self->expr_dir, "probe.no.match"),
474 :     $expt_dir,
475 :     $self->expr_dir]);
476 :    
477 : olson 1.2 my $output = catfile($self->expr_dir, "rma_normalized.tab");
478 :     if (! -f $output)
479 :     {
480 :     confess("Output file $output was not generated");
481 :     }
482 :     }
483 :    
484 :     sub compute_atomic_regulons
485 :     {
486 : olson 1.4 my($self, $pearson_cutoff) = @_;
487 :    
488 :     $pearson_cutoff ||= 0.7;
489 : olson 1.2
490 :     my $coreg_clusters = catfile($self->expr_dir, "coregulated.clusters");
491 :     my $coreg_subsys = catfile($self->expr_dir, "coregulated.subsys");
492 :     my $merged_clusters = catfile($self->expr_dir, "merged.clusters");
493 : olson 1.4 my $probes_always_on = catfile($self->expr_dir, "probes.always.on");
494 :     my $pegs_always_on = catfile($self->expr_dir, "pegs.always.on");
495 :    
496 :     $self->run([executable_for("call_coregulated_clusters_on_chromosome"), $self->expr_dir],
497 :     { stdout => $coreg_clusters });
498 :     $self->run([executable_for("make_coreg_conjectures_based_on_subsys"), $self->expr_dir],
499 :     { stdout => $coreg_subsys });
500 :    
501 :     $self->run([executable_for("filter_and_merge_gene_sets"), $self->expr_dir, $coreg_clusters, $coreg_subsys],
502 :     { stdout => $merged_clusters });
503 :     $self->run([executable_for("get_ON_probes"), $self->expr_dir, $probes_always_on, $pegs_always_on]);
504 :    
505 : olson 1.5 if (-s $probes_always_on == 0)
506 :     {
507 :     confess "No always-on probes were found";
508 :     }
509 :    
510 : olson 1.4 $self->run([executable_for("Pipeline"), $pegs_always_on, $merged_clusters, $self->expr_dir],
511 :     { stdout => catfile($self->expr_dir, "comments.by.Pipeline.R") });
512 :    
513 :     $self->run([executable_for("SplitGeneSets"), $merged_clusters, $pearson_cutoff, $self->expr_dir],
514 :     { stdout => catfile($self->expr_dir, "split.clusters") });
515 :    
516 : olson 1.5 $self->run([executable_for("compute_atomic_regulons_for_dir"), $self->expr_dir]);
517 : olson 1.4 }
518 :    
519 :     sub run
520 :     {
521 :     my($self, $cmd, $redirect) = @_;
522 :    
523 :     print "Run @$cmd\n";
524 :     my $rc = system_with_redirect($cmd, $redirect);
525 :     if ($rc != 0)
526 :     {
527 :     confess "Command failed: @$cmd\n";
528 :     }
529 : olson 1.2 }
530 :    
531 : olson 1.3 sub get_experiment_names
532 :     {
533 :     my($self) = @_;
534 :     my $f = catfile($self->expr_dir, "experiment.names");
535 :     my $fh;
536 :     open($fh, "<", $f) or confess "Could not open $f: $!";
537 :     my @out = map { chomp; my($num, $name) = split(/\t/); $name } <$fh>;
538 :     close($fh);
539 :     return @out;
540 :     }
541 :    
542 :     sub get_experiment_on_off_calls
543 :     {
544 :     my($self, $expt) = @_;
545 :    
546 :     my $f= catfile($self->expr_dir, "final_on_off_calls.txt");
547 :     my $fh;
548 :     open($fh, "<", $f) or confess "Could not open $f: $!";
549 :     my $names = <$fh>;
550 :     chomp $names;
551 :     my @names = split(/\t/, $names);
552 :     my $idx = 0;
553 :     my $expt_idx;
554 :     foreach my $n (@names)
555 :     {
556 :     if ($n eq $expt)
557 :     {
558 :     $expt_idx = $idx;
559 :     last;
560 :     }
561 :     $idx++;
562 :     }
563 :     if (!defined($expt_idx))
564 :     {
565 :     confess("Could not find experiment $expt in $f");
566 :     }
567 :    
568 :     my $calls = {};
569 :     while (<$fh>)
570 :     {
571 :     chomp;
572 :     my($peg, @calls) = split(/\t/);
573 :     #
574 :     # +1 because the row[0] element is the peg, and our index is
575 :     # zero-based.
576 :     #
577 :     $calls->{$peg} = $calls[$expt_idx + 1];
578 :     }
579 :    
580 :     close($fh);
581 :     return($calls);
582 :    
583 :     }
584 :    
585 :     =head3 save_model_gene_activity
586 :    
587 :     $e->save_model_gene_activity($data)
588 :    
589 :     Save the results of a modeling run for a given experiment.
590 :    
591 :     $data is of the form { experiment_id => $data_hash }
592 :    
593 :     =cut
594 :    
595 :     sub save_model_gene_activity
596 :     {
597 :     my($self, $data) = @_;
598 :     }
599 : olson 1.2
600 : olson 1.4 sub all_features
601 :     {
602 :     my($self, $type) = @_;
603 :    
604 :     my @ftypes;
605 :     my $fdir = catfile($self->genome_dir, "Features");
606 :     if (defined($type))
607 :     {
608 :     @ftypes = ($type);
609 :     }
610 :     else
611 :     {
612 :     opendir(D, $fdir);
613 :     @ftypes = grep { -f catfile($fdir, $_) && /^\./ } readdir(D);
614 :     closedir(D);
615 :     }
616 :     my @out;
617 :     for my $ftype (@ftypes)
618 :     {
619 :     if (open(TBL, "<", catfile($fdir, $ftype, "tbl")))
620 :     {
621 :     push(@out, map { /^(\S+)/; $1 } <TBL>);
622 :     close(TBL);
623 :     }
624 :     }
625 :     return @out;
626 :     }
627 :    
628 :     sub fid_locations
629 :     {
630 :     my($self, $fids) = @_;
631 :    
632 :     my %fids;
633 :     $fids{$_}++ for @$fids;
634 :    
635 :     my $genome_id = $self->genome_id;
636 :    
637 :     my $fdir = catfile($self->genome_dir, "Features");
638 :     opendir(D, $fdir);
639 :     my @ftypes = grep { -d catfile($fdir, $_) && ! /^\./ } readdir(D);
640 :     closedir(D);
641 :    
642 :     my $out = {};
643 :     for my $ftype (@ftypes)
644 :     {
645 :     if (open(TBL, "<", catfile($fdir, $ftype, "tbl")))
646 :     {
647 :     while (<TBL>)
648 :     {
649 :     my($id, $locs) = /^(\S+)\t(\S+)\t/;
650 :    
651 :     if ($fids{$id})
652 :     {
653 :     $out->{$id} = "$genome_id:" . SeedUtils::boundary_loc($locs);
654 :     }
655 :     }
656 :     close(TBL);
657 :     }
658 :     }
659 :     return $out;
660 :     }
661 :    
662 :     sub ids_in_subsystems
663 :     {
664 :     my($self) = @_;
665 :    
666 :     open(SS, "<", catfile($self->genome_dir, "subsystem.data"));
667 :     my $res = {};
668 :     while (<SS>)
669 :     {
670 :     chomp;
671 :     my($peg, $ss, $role, $variant) = split(/\t/);
672 :     $ss =~ s/\s+/_/g;
673 :     push(@{$res->{$ss}->{$role}}, $peg);
674 :     }
675 :     close(SS);
676 :     return $res;
677 :     }
678 :    
679 :     sub ids_to_subsystems
680 :     {
681 :     my($self, $ids) = @_;
682 :    
683 :     my %ids;
684 :     $ids{$_} = 1 for @$ids;
685 :    
686 :     open(SS, "<", catfile($self->genome_dir, "subsystem.data"));
687 :     my $res = {};
688 :     while (<SS>)
689 :     {
690 :     chomp;
691 :     my($peg, $ss, $role, $variant) = split(/\t/);
692 :     if ($ids{$peg})
693 :     {
694 :     push(@{$res->{$peg}}, $ss);
695 :     }
696 :     }
697 :     close(SS);
698 :     return $res;
699 :     }
700 :    
701 :     sub ids_to_functions
702 :     {
703 :     my($self, $ids) = @_;
704 :     open(AF, "<", catfile($self->genome_dir, "assigned_functions"));
705 :     my %ids;
706 :     $ids{$_} = 1 for @$ids;
707 :     my $res = {};
708 :    
709 :     while (<AF>)
710 :     {
711 :     chomp;
712 :     my($id, $fn) = split(/\t/);
713 :     $res->{$id} = $fn if $ids{$id};
714 :     }
715 :     close(AF);
716 :     return $res;
717 :     }
718 :    
719 : olson 1.1 1;
720 : olson 1.2
721 :    

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3