/tools/human_genome_variation/linkToGProfile.pl

https://bitbucket.org/cistrome/cistrome-harvard/ · Perl · 55 lines · 40 code · 7 blank · 8 comment · 7 complexity · ac44394bc1c449d1e413adb287473faa MD5 · raw file 

    

  1. #!/usr/bin/env perl
    2. 

  2. 

  3. use strict;
    4. 

  4. use warnings;
    5. 

  5. 

  6. ###################################################
    7. 

  7. # linkToGProfile.pl
    8. 

  8. # Generates a link to gprofile for a list of gene IDs.
    9. 

  9. # g:Profiler a web-based toolset for functional profiling of gene lists from large-scale experiments (2007) NAR 35 W193-W200
    10. 

  10. ###################################################
    11. 

  11.  
    12. 

  12. if (!@ARGV or scalar @ARGV != 4) {
    13. 

  13.    print "usage: linkToGProfile.pl infile.tab 1basedCol idType outfile\n";
    14. 

  14.    exit 1;
    15. 

  15. }
    16. 

  16. 

  17. my $in = shift @ARGV;
    18. 

  18. my $col = shift @ARGV;
    19. 

  19. my $type = shift @ARGV;
    20. 

  20. my $out = shift @ARGV;
    21. 

  21. 

  22. if ($col < 1) { 
    23. 

  23.    print "ERROR the column number should be 1 based counting\n";
    24. 

  24.    exit 1;
    25. 

  25. }
    26. 

  26. my @gene;
    27. 

  27. open(FH, $in) or die "Couldn't open $in, $!\n";
    28. 

  28. while (<FH>) {
    29. 

  29.    chomp;
    30. 

  30.    my @f = split(/\t/);
    31. 

  31.    if (scalar @f < $col) {
    32. 

  32.       print "ERROR there is no column $col in $in\n";
    33. 

  33.       exit 1;
    34. 

  34.    }
    35. 

  35.    if ($f[$col-1]) { push(@gene, $f[$col-1]); }
    36. 

  36. }
    37. 

  37. close FH or die "Couldn't close $in, $!\n";
    38. 

  38.  
    39. 

  39. my $link = 'http://biit.cs.ut.ee/gprofiler/index.cgi?organism=hsapiens&query=GENELIST&r_chr=1&r_start=start&r_end=end&analytical=1&domain_size_type=annotated&term=&significant=1&sort_by_structure=1&user_thr=1.00&output=png&prefix=TYPE';
    40. 

  40. $link =~ s/TYPE/$type/;
    41. 

  41. my $g = join("+", @gene);
    42. 

  42. $link =~ s/GENELIST/$g/;
    43. 

  43. #print output
    44. 

  44. if (length $link > 2048) { 
    45. 

  45.    print "ERROR too many genes to fit in URL, please select a smaller set\n";
    46. 

  46.    exit;
    47. 

  47. }
    48. 

  48. open(FH, ">", $out) or die "Couldn't open $out, $!\n";
    49. 

  49. print FH "<html><head><title>g:Profiler link</title></head><body>\n",
    50. 

  50.       '<A TARGET=_BLANK HREF="', $link, '">click here to send list of identifiers to g:Profiler</A>', "\n",
    51. 

  51.       '</body></html>', "\n";
    52. 

  52. close FH or die "Couldn't close $out, $!\n";
    53. 

  53. 

  54. #also do link that prints text that could be pulled back into Galaxy?
    55. 

  55. exit;
    56.