0
|
1 #!/usr/bin/perl
|
|
2 # computes difference distances between motif models
|
|
3
|
|
4 if( $#ARGV != 1 ){
|
|
5 print "Usage: ./Expression_Pattern_Distance.pl model1 model2\n";
|
|
6 exit;
|
|
7 }
|
|
8
|
|
9 @model1 = &read_model( $ARGV[0] );
|
|
10 @model2 = &read_model( $ARGV[1] );
|
|
11 if( $#model1 != $#model2 ){
|
|
12 print "Error: model dimensions do not match\n";
|
|
13 exit;
|
|
14 }
|
|
15
|
|
16
|
|
17 #
|
|
18 # second distance: sum of absolute difference, divided by absolute sum
|
|
19 # bounded between 0 and 1; equals 1 if all signs are different
|
|
20 #
|
|
21 $dist = 0;
|
|
22 $denom = 0;
|
|
23 for($i=0; $i< $classes; $i++){
|
|
24 for($j=0; $j< $dims; $j++){
|
|
25 if( $model1[$i][$j] ne "nan" && $model2[$i][$j] ne "nan" ){
|
|
26 $dist+= abs($model1[$i][$j] - $model2[$i][$j]);
|
|
27 $denom+= abs($model1[$i][$j]) + abs($model2[$i][$j]);
|
|
28 }
|
|
29 }
|
|
30 }
|
|
31 $dist/=$denom;
|
|
32
|
|
33 printf("%0.3f\n", $dist);
|
|
34
|
|
35
|
|
36 exit;
|
|
37
|
|
38
|
|
39
|
|
40
|
|
41 ################################
|
|
42
|
|
43 # function to read a motif model from file $_[0]
|
|
44 sub read_model{
|
|
45 $dims = 0;
|
|
46 $classes = 0;
|
|
47 open(IN,"$_[0]") || die "Error: can't open file $_[0]\n";
|
|
48 @modpars = ();
|
|
49 # burn the first two lines
|
|
50 $line = <IN>;
|
|
51 $line = <IN>;
|
|
52 while( $line = <IN> ){
|
|
53 chomp($line);
|
|
54 @parts = split(' ',$line);
|
|
55 if($parts[0] !~ /#/ ){
|
|
56 @parts2 = ();
|
|
57 for($i=1; $i<= $#parts; $i++){
|
|
58 if( $parts[$i] ne "NA" ){
|
|
59 push(@parts2, $parts[$i]);
|
|
60 }
|
|
61 else{
|
|
62 push(@parts2, 0);
|
|
63 }
|
|
64 }
|
|
65 push(@modpars, [@parts2] );
|
|
66 $classes++;
|
|
67 $dims = $#parts;
|
|
68 }
|
|
69 }
|
|
70 close(IN);
|
|
71 return @modpars;
|
|
72 }
|
|
73
|