forester/archive/perl/bootstrapCounter.pl

   1 #!/usr/bin/perl -w
   2 #
   3 # bootstrapCounter.pl
   4 # -------------------
   5 #
   6 # Copyright (C) 2001 Washington University School of Medicine
   7 # and Howard Hughes Medical Institute
   8 # All rights reserved
   9 #
  10 # Author: Christian M. Zmasek
  11 #         zmasek@genetics.wustl.edu
  12 #         http://www.genetics.wustl.edu/eddy/people/zmasek/
  13 #
  14 # Created: 04/04/01
  15 #
  16 # Last modified 08/16/01
  17 #
  18 #
  19 # Objective. Determines the distribution of top orthology bootstrap values
  20 #            of a Xrio.pl output file.
  21 #
  22 # Usage. "bootstrapCounter.pl <infile = Xrio.pl-output> <outfile>"
  23 #
  24 # Important. The result of this is meaningful ONLY if the thresholds
  25 #            for output of the RIO analysis are set to zero (L=0 R=0).
  26 #
  27 # Format for infile:
  28 # ...
  29 #
  30 # # ############################################################################
  31 # # Annotation: B0511.6 CE17345   helicase (ST.LOUIS) TR:O61815 protein_id:AAC17654.1
  32 # # HMM       : ABC_tran
  33 # # score     : -59.6
  34 # # E-value   : 1.1
  35 # # Query has not been aligned (score lower than gathering cutoff).
  36 # # ############################################################################
  37 #
  38 #
  39 # # ############################################################################
  40 # # Annotation: B0511.7 CE17346    (ST.LOUIS) TR:O61817 protein_id:AAC17655.1
  41 # # HMM       : FHA
  42 # # score     : 71.6
  43 # # E-value   : 1.7e-17
  44 # RIO - Resampled Inference of Orthologs
  45 # Version: 1.000
  46 # ------------------------------------------------------------------------------
  47 # Alignment file: /tmp/Xriopl9846081980/Full-FHA
  48 # Alignment     : FHA domain
  49 # HMM           : FHA
  50 # Query file    : /tmp/Xriopl9846081980/__queryfile__
  51 # ==============================================================================
  52 #
  53 # Query         : CE17346.FHA_CAEEL/45-114
  54 #
  55 # Number (in %) of observed orthologies (o) and super orthologies (s) to query
  56 # in bootstrapped trees, evolutionary distance to query:
  57 #
  58 # Sequence              Description                                                                  # o[%] s[%]  distance
  59 # --------              -----------                                                                  ---- ----  --------
  60 # YC67_MYCTU/308-372    -                                                                              20   14  1.577840
  61 # FRAH_ANASP/204-277    FRAH PROTEIN.                                                                  17   16  1.532670
  62 # ABA2_NICPL/557-633    ZEAXANTHIN EPOXIDASE PRECURSOR (EC 1.14.-.-).                                  14   11  1.885700
  63 # ABA2_LYCES/563-639    ZEAXANTHIN EPOXIDASE PRECURSOR (EC 1.14.-.-).                                  14   11  2.140000
  64 #
  65 #
  66 #
  67 # Distance values (based on ML branch length values on consensus tree)
  68 # --------------------------------------------------------------------
  69 # Given the thresholds for distance calculations:
  70 # No sequence is considered orthologous to query.
  71 #
  72 # ...
  73
  74
  75
  76 use strict;
  77
  78 my $VERSION            = 0.200;
  79
  80 my $infile             = "";
  81 my $outfile            = "";
  82 my $return_line        = "";
  83 my $results            = 0;
  84 my $o_bootstraps       = 0;
  85 my $s_bootstraps       = 0;
  86 my @o_bootstraps_array = ();
  87 my @s_bootstraps_array = ();
  88 my $total              = 0;
  89 my $i                  = 0;
  90
  91
  92 if ( @ARGV != 2 ) {
  93     &errorInCommandLine();
  94     exit ( -1 );
  95 }
  96
  97 $infile  = $ARGV[ 0 ];
  98 $outfile = $ARGV[ 1 ];
  99
 100 if ( -e $outfile ) {
 101     die "\n$0: <<$outfile>> already exists.\n";
 102 }
 103 unless ( ( -s $infile ) && ( -f $infile ) && ( -T $infile ) ) {
 104     die "\n$0: <<$infile>> does not exist, is empty, or is not a plain textfile.\n";
 105 }
 106
 107
 108 open( IN, "$infile" ) || die "\n$0: Cannot open file <<$infile>>: $!\n";
 109
 110 $results = 0;
 111 for ( $i = 0; $i <= 100; ++$i ) {
 112     $s_bootstraps_array[ $i ] = $o_bootstraps_array[ $i ] = 0;
 113 }
 114
 115 while ( $return_line = <IN> ) {
 116
 117     if ( $return_line =~ /^\s*--------\s+/ ) {
 118         $results = 1;
 119     }
 120     elsif ( $return_line =~ /^\s*Distance\s+values\s+/i ) {
 121         $results = 0;
 122     }
 123     elsif ( $results == 1 && $return_line =~ /^\s*!NO\s+ORTHOLOGS/ ) {
 124         $o_bootstraps_array[ 0 ]++;
 125         $s_bootstraps_array[ 0 ]++;
 126         $total++;
 127         $results = 0;
 128     }
 129     elsif ( $results == 1 && $return_line =~ /(\S+)\s+(\S+)\s+\S+\s*$/ ) {
 130         $o_bootstraps = $1;
 131         $s_bootstraps = $2;
 132         $results = 0;
 133         if ( $o_bootstraps > 100 || $s_bootstraps > 100
 134         || $o_bootstraps < 0 ) {
 135             print "o bootstraps: $o_bootstraps\n";
 136             print "s bootstraps: $s_bootstraps\n";
 137             die "\n\n$0: Error: Boostrap value(s) out of range.\n\n";
 138         }
 139
 140         $total++;
 141         $o_bootstraps_array[ $o_bootstraps ]++;
 142         $s_bootstraps_array[ $s_bootstraps ]++;
 143
 144     }
 145 }
 146
 147 close( IN );
 148
 149
 150 open( OUT, ">$outfile" ) || die "\n$0: Cannot create file \"$outfile\": $!\n";
 151
 152 print OUT "bootstrapCounter.pl version: $VERSION\n\n";
 153 print OUT "Distribution of top bootstrap values\n\n";
 154 print OUT "Input file : $infile\n";
 155 print OUT "Output file: $outfile\n";
 156 print OUT "Date       : ".`date`."\n";
 157 print OUT "Total: $total\n\n";
 158 print OUT "top-orthology-bootstraps vs. count:\n";
 159 for ( $i = 0; $i < @o_bootstraps_array; ++$i ) {
 160     print OUT "$i $o_bootstraps_array[ $i ]\n";
 161 }
 162 print OUT "\n\ntop-super-orthology-bootstraps vs. count:\n";
 163 for ( $i = 0; $i < @s_bootstraps_array; ++$i ) {
 164     print OUT "$i $s_bootstraps_array[ $i ]\n";
 165 }
 166 close( OUT );
 167
 168 print( "\nDone.\n\n" );
 169
 170 exit( 0 );
 171
 172
 173
 174 sub errorInCommandLine {
 175     print "\n";
 176     print " bootstrapCounter.pl version: $VERSION\n";
 177     print " Usage: \"bootstrapCounter.pl <infile = Xrio.pl-output> <outfile>\"\n";
 178     print " Important: The result of this is meaningful ONLY if the thresholds\n";
 179     print " for output of the RIO analysis are set to zero (L=0 R=0).\n";
 180     print "\n";
 181     exit( -1 );
 182 }
 183
 184