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

Diff of /FigKernelScripts/dp.pl

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1, Thu Jun 11 19:35:46 2009 UTC revision 1.3, Thu Jun 18 14:42:48 2009 UTC
# Line 1  Line 1 
1    #!/usr/bin/env /scratch/Ross/FIGdisk/env/cee/bin/perl
2    
3    BEGIN {
4        @INC = qw(
5                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib
6                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/FigKernelPackages
7                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/WebApplication
8                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/FortyEight
9                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/PPO
10                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/RAST
11                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/MGRAST
12                  /disks/space0/scratch/Ross/FIGdisk/dist/releases/dev/cee/lib/SeedViewer
13                  /scratch/Ross/FIGdisk/dist/dev/cee/lib
14                  /scratch/Ross/FIGdisk/dist/dev/cee/lib/FigKernelPackages
15                  /scratch/Ross/FIGdisk/env/cee/lib/perl5/5.8.8/x86_64-linux
16                  /scratch/Ross/FIGdisk/env/cee/lib/perl5/5.8.8
17                  /scratch/Ross/FIGdisk/env/cee/lib/perl5/site_perl/5.8.8/x86_64-linux
18                  /scratch/Ross/FIGdisk/env/cee/lib/perl5/site_perl/5.8.8
19                  /scratch/Ross/FIGdisk/env/cee/lib/perl5/site_perl
20                  .
21                  /scratch/Ross/FIGdisk/config
22    
23    );
24    }
25    use Data::Dumper;
26    use Carp;
27    use FIG_Config;
28    $ENV{'BLASTMAT'} = "/scratch/Ross/FIGdisk/BLASTMAT";
29    $ENV{'FIG_HOME'} = "/scratch/Ross/FIGdisk";
30    # end of tool_hdr
31    ########################################################################
32  #  Use  #  Use
33  #  #
34  #       perl dp.pl -k 4 sample1 sample2  #       perl dp.pl -k 4 sample1 sample2
# Line 7  Line 38 
38  #  #
39  #  Here the "dot product" will return a value from 0 to 1, where the higher the value  #  Here the "dot product" will return a value from 0 to 1, where the higher the value
40  #  the more similar the samples.  #  the more similar the samples.
41    #
42    #  The command will work if either file contains a single line containing a fabric.
43    #
44    
45  use Carp;  use Carp;
46  use Data::Dumper;  use Data::Dumper;
47  use strict;  use strict;
48    
49  my($N,$sample1,$sample2);  my($N,$sample1,$sample2);
50  my $usage = "usage: dp [-k n] Sample1 Sample2";  my $usage = "usage: dp [-k n] [-j] Sample1 Sample2";
51    
52    my $jaccard = 0;
53  $N = 4;  $N = 4;
54  while ( @ARGV && $ARGV[0] =~ s/^-// )  while ( @ARGV && $ARGV[0] =~ s/^-// )
55  {  {
56      $_ = shift;      $_ = shift;
57      if ( s/^k// ) { $N    = $_ || shift; next }      if ( s/^k// ) { $N    = $_ || shift; next }
58        elsif   ( s/^j// ) { $jaccard = 1 }
59      else      else
60      {      {
61          print STDERR "Bad flag '$_'\n", $usage;          print STDERR "Bad flag '$_'\n", $usage;
# Line 37  Line 73 
73  my $v1       = &load_sample($sample1,$role_num,$N);  my $v1       = &load_sample($sample1,$role_num,$N);
74  my $v2       = &load_sample($sample2,$role_num,$N);  my $v2       = &load_sample($sample2,$role_num,$N);
75    
76  my $dot_prod = &dot_prod($v1,$v2);  # note we compute the dot prod of two unit vectors created  my $dot_prod = &similarity($v1,$v2,$jaccard);  # note we compute the dot prod of two unit vectors created
77                                      # by normalizing v1 and v2                                      # by normalizing v1 and v2
78    
79  print "$dot_prod\n";  print "$dot_prod\n";
# Line 67  Line 103 
103      {      {
104          chomp;          chomp;
105          my @flds = split(/\t/,$_);          my @flds = split(/\t/,$_);
106            if (@flds > 8000)
107            {
108                @$role_counts = @flds;
109            }
110            else
111            {
112          my $func = $flds[$N-1];          my $func = $flds[$N-1];
113          if ($func)          if ($func)
114          {          {
# Line 81  Line 123 
123              }              }
124          }          }
125      }      }
126        }
127        close(IN);
128      return $role_counts;      return $role_counts;
129  }  }
130    
# Line 109  Line 153 
153      return $uv;      return $uv;
154  }  }
155    
156  sub dot_prod {  sub similarity {
157      my($v1,$v2) = @_;      my($v1,$v2,$jaccard) = @_;
158    
159      $v1 = &unit_vector($v1);      $v1 = &unit_vector($v1);
160      $v2 = &unit_vector($v2);      $v2 = &unit_vector($v2);
161    
162        if ($jaccard)
163        {
164            my $dp = &dot_prod($v1,$v2);
165            my @v3 = @$v1;
166            my $i;
167            for ($i=0; ($i < @$v2); $i++)
168            {
169                if ($v3[$i] && $v2->[$i])
170                {
171                    $v3[$i] += $v2->[$i];
172                }
173                elsif ($v2->[$i])
174                {
175                    $v3[$i] = $v2->[$i];
176                }
177            }
178            my $sq = &dot_prod(\@v3,\@v3);
179            return $dp / $sq;
180        }
181        else
182        {
183            return &dot_prod($v1,$v2);
184        }
185    }
186    
187    sub dot_prod {
188        my($v1,$v2) = @_;
189    
190    
191      my $tot = 0;      my $tot = 0;
192      my $i;      my $i;
193      for ($i=0; ($i < @$v1); $i++)      for ($i=0; ($i < @$v1); $i++)
194      {      {
195          if (defined($v1->[$i]) && defined($v2->[$i]))          if ($v1->[$i] && $v2->[$i])
196          {          {
197              $tot += $v1->[$i] * $v2->[$i];              $tot += $v1->[$i] * $v2->[$i];
198          }          }

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.3

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3