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

View of /FigKernelScripts/submit_ff_test.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.6 - (download) (as text) (annotate)
Thu Oct 8 18:54:18 2009 UTC (10 years, 5 months ago) by arodri7
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.5: +34 -10 lines
commit changes

#!/usr/bin/env /home/arodri7/FIGdisk/bin/run_perl

use Data::Dumper;
use Carp;
use FIG_Config;
use strict;
use FIG;
my $fig = new FIG;

my ($ffD, $genome, $outputD, $kmers) = @ARGV;
print STDERR "processing $genome\n";

my $assignF = "$outputD/$genome.assigned";
my $not_assignF = "$outputD/$genome.not_assigned";
#open(OUT,"| assign_using_ff -b $ffD > $assignF 2>&1") || die "could not open $assignF";
my $tmpF = "$FIG_Config::temp/fasta_figfam_qc.$$";
open (TMP, ">$tmpF");
foreach my $peg ($fig->pegs_of($genome))
{
    if (my $prot = $fig->get_translation($peg))
    {
	print TMP ">$peg\n$prot\n";
    }
}
close TMP;

if ($kmers)
{
    &FIG::run("(assign_using_ff_server_simple < $tmpF >> $assignF) 2>> $not_assignF");
#    open(OUT,"| (assign_using_ff_server_simple < $tmpF >> $assignF )>>& $not_assignF") || die "could not open $assignF";
#    close(OUT);
}
else
{
    open(OUT,"| assign_using_ff_parallel -s $ffD 8 $assignF < $tmpF 2>&1") || die "could not open $assignF";
    close(OUT);
}

if ($kmers)
{
    `cat $not_assignF >> $assignF`;
    `rm $not_assignF`;
}

open (FH,"<$assignF") || die "could not open $assignF";
open (OUT, ">$outputD/$genome.compare_log");
my ($total_count, $match_count, $diff_count, $unannotated_count, $annotated_count, $misannotated);
while (<FH>){
  my ($peg_id, $annotation);
  if ($_ =~ m/^(fig\|.*?)\s+was not placed into a FIGfam\n/){
    $peg_id = $1;
    $annotation = "NOT ANNOTATED";
    $unannotated_count++;
  }
  elsif ($_ =~ m/^(fig\|.*?)\t(.*)\n/){
    $peg_id = $1;
    $annotation = $2;
    $annotated_count++;
  }

  if ($peg_id){
      # determine if the seed annotation and the figfam annotation are the same
      my $seed_annotation = $fig->function_of($peg_id);

      my ($match_status);
      if ($seed_annotation eq $annotation){
	  $match_status = "MATCH";
	  $match_count++;
      }
      else{
	  if ($annotation ne "NOT ANNOTATED"){
	      $misannotated++;
	  }
	  $match_status = "DIFF";
	  $diff_count++;
      }
  
      print OUT "$match_status\t$peg_id\t$annotation\t$seed_annotation\n";
  }
}
print OUT "MATCH: $match_count\nDIFF: $diff_count\nANNOTATED: $annotated_count\nUNANNOTATED: $unannotated_count\n";
close FH;
close OUT;

open (FH, "$ffD/release_history/VERSION");
my $name;
while (my $line = <FH>){
    chomp ($line);
    if ($line =~ /^version/){
	($name) = ($line) =~ /version\t(.*)/;
	last;
    }
}
close FH;

$total_count = $annotated_count+$unannotated_count;
open (STATS, ">>$ffD/release_history/QC.dat");

print STATS "NAME\t$name\n";
print STATS "ORGANISM\t$genome\n";
print STATS "TOTAL\t$total_count\n";
print STATS "ANNOTATED\t$annotated_count\n";
print STATS "MISANNOTATED\t$misannotated\n";
print STATS "//\n";
close STATS;

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3