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

View of /FigKernelScripts/svr_update_assignments.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (download) (as text) (annotate)
Tue Jan 19 18:14:49 2010 UTC (9 years, 10 months ago) by gdpusch
Branch: MAIN
CVS Tags: mgrast_dev_08112011, mgrast_dev_08022011, rast_rel_2014_0912, myrast_rel40, mgrast_dev_05262011, mgrast_dev_04082011, rast_rel_2010_0928, mgrast_version_3_2, mgrast_dev_12152011, mgrast_dev_06072011, rast_rel_2010_0526, rast_rel_2014_0729, mgrast_dev_02212011, rast_rel_2010_1206, mgrast_release_3_0, mgrast_dev_03252011, rast_rel_2010_0118, rast_rel_2011_0119, mgrast_release_3_0_4, mgrast_release_3_0_2, mgrast_release_3_0_3, mgrast_release_3_0_1, mgrast_dev_03312011, mgrast_release_3_1_2, mgrast_release_3_1_1, mgrast_release_3_1_0, mgrast_dev_04132011, mgrast_dev_04012011, rast_rel_2010_0827, myrast_33, rast_rel_2011_0928, mgrast_dev_04052011, mgrast_dev_02222011, mgrast_dev_10262011, HEAD
Changes since 1.2: +1 -1 lines
Fixed typo.

# -*- perl -*-

use SeedEnv;
use strict;
use warnings;

my $usage = "usage: svr_update_assignments [-suggest_only] [-no_blast] SeedDir AssignmentsFile";

my $do_blast = 1;
my $in_place_update = 1;
while (@ARGV && ($ARGV[0] =~ m/^-/)) {
    if    ($ARGV[0] =~ m/-suggest_only/) {
	$in_place_update = 0;
    }
    elsif ($ARGV[0] =~ m/-no_blast/) {
	$do_blast = 0;
    }
    else {
	die "Unknown argument $ARGV[0]\n$usage";
    }
    shift @ARGV;
}

my($dir,$assignmentsF,$fh);
(
 ($dir          = shift @ARGV) && (-d $dir) && open($fh,"<$dir/Features/peg/fasta") &&
 ($assignmentsF = shift @ARGV)
)
    || die $usage;

my %assignments;
if (-s $assignmentsF)
{
    %assignments = map { $_ =~ /^(\S+)\t(\S.*\S)/ ? ($1 => $2) : () } `cat $assignmentsF`;
}

my $annoObject = ANNOserver->new();
my $kmer_func_results_handle  = $annoObject->assign_function_to_prot( -input => $fh, -assignToAll => $do_blast, -kmer => 8);


if ($in_place_update) {
    open(KMEROV,">$dir/Kmer.overrides")
	|| die "Could not write-open Kmer overrides file \'$dir/Kmer.overrides\'";
}

while (my $tuple = $kmer_func_results_handle->get_next)
{
    my($peg, $function, $otu, $score, $non_overlapping,$overlapping) = @$tuple;
    my $old = $assignments{$peg};
    if ((not $old) || (($old ne $function) && ($non_overlapping >= 3)))
    {
	$assignments{$peg} = $function;
	print KMEROV "$peg\t$function\n" if $in_place_update;
    }
}
close(KMEROV) if $in_place_update;

my $assgn_fh = \*STDOUT;
if ($in_place_update) {
    open($assgn_fh, ">$assignmentsF")
	|| die "Could not write-open assignments file \'$assignmentsF\'";
}

foreach my $peg (sort { &SeedUtils::by_fig_id($a, $b) } (keys %assignments)) {
    print $assgn_fh ($peg, "\t", $assignments{$peg}, "\n");
}
close($assgn_fh) if $in_place_update;


#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#...Update Subsystems and Bindings...
#-----------------------------------------------------------------------
if ($in_place_update) {
    if (-d "$dir/Subsystems")
    {
	system "/bin/rm -rf $dir/Subsystems";
    }
    mkdir("$dir/Subsystems",0777) || die "could not make $dir/Subsystems";
    
    my @roles = map { [$assignments{$_},$_] } keys(%assignments);
    print STDERR &Dumper(\@roles);
    
    my $reconstruction = $annoObject->metabolic_reconstruction( -roles => \@roles);

    open(BINDINGS,">$dir/Subsystems/bindings") || die "could not open $dir/Subsystems/bindings";
    my %subsys;
    foreach my $tuple (@$reconstruction)
    {
	my($ss_and_var,$role,$peg) = @$tuple;
	if ($ss_and_var =~ /^(.*):([^:]+)$/)
	{
	    my $ss = $1;
	    my $var = $2;
	    $subsys{$ss} = $var;
	    print BINDINGS join("\t",($ss,$role,$peg)),"\n";
	}
    }
    close(BINDINGS);

    open(SS,">$dir/Subsystems/subsystems") || die "could not open $dir/Subsystems/subsystems";
    foreach $_ (sort keys(%subsys))
    {
	print SS join("\t",($_,$subsys{$_})),"\n";
    }
    close(SS);
}


						       

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3