-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathGEO_postprocess.pl
103 lines (55 loc) · 3.05 KB
/
GEO_postprocess.pl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
#!/usr/bin/perl
use strict;
#use Clone qw(clone);
unless (@ARGV == 2) {
&USAGE;
}
sub USAGE {
die '
perl ~/bin/perl/GEO_postprocess.pl <GEOres.fit.txt> <GPL2995.annot.gz>
GO-file has this structure:
Gene_ID <TAB> GOterm,GOterm,GOterm
Prod file has this structure:
Gene_ID <TAB> product
';
}
# check that DESeq ran okay
my $dir = shift;
#my $go = shift;
my $prod = shift;
unless (-s "$dir") {
die "\nSorry, the file $dir wasnt made\n";
}
# make nice scores
system "cat $dir | sed \'s/\"//g\' | tr \' \' \'\\t\' | sed \'s/ID=//g\' > $dir.scores.dat";
# make nice down
system "cat $dir.scores.dat | grep -v logFC | awk '\$6\<0.05' | awk '\$2\~\"-\"' | awk '\{print \$0\"\\tDOWN\"\}' > $dir.DOWN";
system "cat $dir.DOWN | cut -f1 | awk \'{print \$1}\' > $dir.down.list ";
system "gzcat $prod | cut -f1,2,3 | grep -w -f $dir.down.list > $dir.down.list.annot ";
#print "cat $dir.scores.dat | grep -v baseMean | awk '\$7\<0.05' | awk '\$3\~\"-\"' | awk '\{print \$0\"\\tDOWN\"\}' > $dir.DOWN\n";
#print "cat $dir.DOWN | cut -f1 | awk \'{print \$1}\' > $dir.down.list \n";
#print "cat $dir.scores.dat | awk '\$8\<0.05' | awk '\$3!\~\"-\"' | awk '\{print \$0\"\\tUP\"\}' > $dir.down\n";
#system "cat $dir/lib.dat.down | cut -f2 | awk -F\'\.\' \'{print \$1}\' > $dir/lib.dat.down.list ";
# make nice up
system "cat $dir.scores.dat | grep -v logCPM | awk '\$6\<0.05' | awk '\$2!\~\"-\"' | awk '\{print \$0\"\\tUP\"\}' > $dir.UP";
system "cat $dir.UP | cut -f1 | awk \'{print \$1}\' > $dir.up.list ";
#print "cat $dir.scores.dat | grep -v baseMean | awk '\$7\<0.05' | awk '\$3!\~\"-\"' | awk '\{print \$0\"\\tUP\"\}' > $dir.UP\n";
#print "cat $dir.UP | cut -f1 | awk \'{print \$1}\' > $dir.up.list \n";
system "gzcat $prod | cut -f1,2,3 | grep -w -f $dir.up.list > $dir.up.list.annot ";
# make nice non-sign
system "cat $dir.scores.dat | awk '\$3\>0.05' | awk '\{print \$0\"\\tNS\"\}' > $dir.NS";
#print "cat $dir.scores.dat | awk '\$7\>0.05' | awk '\{print \$0\"\\tNS\"\}' > $dir.NS\n";
#__END__
system "cat $dir.UP $dir.DOWN $dir.NS > $dir.txt ";
#print "cat $dir.UP $dir.DOWN $dir.NS > $dir.txt\n";
#print "\nRunning topGO down \nperl ~/bin/perl/topGO_starter.pl $dir.down.list $go BP $prod $dir.down &\n";
#system "perl ~/bin/perl/topGO_starter.pl $dir.down.list $go BP $prod $dir.down";
#print "\nRunning topGO up \nperl ~/bin/perl/topGO_starter.pl $dir.up.list $go BP $prod $dir.up &\n";
#system "perl ~/bin/perl/topGO_starter.pl $dir.up.list $go BP $prod $dir.up";
print "All done\n";
__END__
perl ~/bin/perl/foreach.pl "cat X | grep DOWN | awk '{print 2}' > X.down" *vs*/lib.dat | sed 's/ 2/ $2/' | tail -100 > down.sh
perl ~/bin/perl/foreach.pl "cat X | grep UP | awk '{print 2}' > X.up" *vs*/lib.dat | sed 's/ 2/ $2/' | tail -100 > up.sh
# do topGO for down and up
perl ~/bin/perl/topGO_starter.pl aPS_vs_aPS/lib.dat.down ../EMU_2013_06_09.gff.GOtop BP ../EMU_2013_06_09.prod aPS_vs_aPS/lib.dat.down;
perl ~/bin/perl/topGO_starter.pl aPS_vs_MCana/lib.dat.down ../EMU_2013_06_09.gff.GOtop BP ../EMU_2013_06_09.prod aPS_vs_MCana/lib.dat.down;