~ubuntu-branches/ubuntu/gutsy/bioperl/gutsy

« back to all changes in this revision

Viewing changes to scripts/contributed/expression_analysis/entropy.pl

  • Committer: Bazaar Package Importer
  • Author(s): Matt Hope
  • Date: 2004-04-18 14:24:11 UTC
  • mfrom: (1.2.1 upstream) (2.1.1 warty)
  • Revision ID: james.westby@ubuntu.com-20040418142411-gr92uexquw4w8liq
Tags: 1.4-1
* New upstream release
* Examples and working code are installed by default to usr/bin,
  this has been moved to usr/share/doc/bioperl/bin

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
#!/usr/local/bin/perl
2
 
 
3
 
# Entropy (2001)
4
 
 
5
 
open (ENT,'*.csv') or die "$!";     # file name
6
 
@allfeat = <ENT>;
7
 
close (ENT);
8
 
 
9
 
open (ENTOUT,'>entropy.csv');
10
 
 
11
 
foreach $feat (@allfeat) {
12
 
    $feat =~ s/\n//;
13
 
    ($name,$a,$b,$c,$d,$e,$f,$g,$h) = split(/,/,$feat);
14
 
    @feature = ($a,$b,$c,$d,$e,$f,$g,$h);
15
 
    @sortfeat = sort {$a <=> $b} @feature;
16
 
    $interval = ($sortfeat[7] - $sortfeat[0])/9;
17
 
    @compart = (0,0,0,0,0,0,0,0,0,0);
18
 
    for ($fe = 0;$fe <= 7;++$fe) {
19
 
        $cell = int (($feature[$fe]-$sortfeat[0])/$interval);
20
 
        $compart[$cell] += 0.125;
21
 
    }
22
 
    $entropy = -($compart[0]*log($compart[0]+0.00001)/log(2)
23
 
        +$compart[1]*log($compart[1]+0.00001)/log(2)
24
 
        +$compart[2]*log($compart[2]+0.00001)/log(2)
25
 
        +$compart[3]*log($compart[3]+0.00001)/log(2)
26
 
        +$compart[4]*log($compart[4]+0.00001)/log(2)
27
 
        +$compart[5]*log($compart[5]+0.00001)/log(2)
28
 
        +$compart[6]*log($compart[6]+0.00001)/log(2)
29
 
        +$compart[7]*log($compart[7]+0.00001)/log(2)
30
 
        +$compart[8]*log($compart[8]+0.00001)/log(2)
31
 
        +$compart[9]*log($compart[9]+0.00001)/log(2));
32
 
    push (@entropy,$entropy);
33
 
    print ENTOUT "$name,$a,$b,$c,$d,$e,$f,$g,$h,$entropy\n";
34
 
}
35
 
 
36
 
close (ENTOUT);
37
 
 
38
 
# Cut features genes with lowest 5% entropy
39
 
 
40
 
@sortent = sort{$a <=> $b} @entropy;
41
 
$entther = ($sortent[*] + $sortent[*])/2;   # *:lowest 5% number of features (genes)
42
 
 
43
 
open (CENT,'entropy.csv') or die "$!";
44
 
@allfeat2 = <CENT>;
45
 
close (CENT);
46
 
 
47
 
open (CENTOUT,'>relation.csv');
48
 
open (DROPOUT,'>dropout.csv');
49
 
 
50
 
foreach $feat2 (@allfeat2) {
51
 
    $feat2 =~ s/\n//;
52
 
    ($name,$a,$b,$c,$d,$e,$f,$g,$h,$entropy) = split(/,/,$feat2);
53
 
    # print "$entther\n";
54
 
    if ($entropy >= $entther) {
55
 
        print CENTOUT "$name,$a,$b,$c,$d,$e,$f,$g,$h\n";
56
 
        $entgene +=1;
57
 
    } else {
58
 
        print DROPOUT "$name,$a,$b,$c,$d,$e,$f,$g,$h\n";
59
 
        $dropgene +=1;
60
 
    }
61
 
}
62
 
 
63
 
print "Gene with sufficient entoropy:$entgene\n";
64
 
print "Dropout gene:$dropgene\n";
65
 
 
66
 
close (CENTOUT);
67
 
close (DROPOUT);