examples/tools/psw.pl

   1 #!/usr/bin/perl
   2
   3 # PROGRAM  : psw.pl
   4 # PURPOSE  : Simple driver for Bio::Tools::pSW
   5 # AUTHOR   : Ewan Birney birney@sanger.ac.uk
   6 # CREATED  : Tue Oct 27 1998
   7 #
   8 # INSTALLATION
   9 #
  10 # you almost certainly have to have installed bioperl
  11 # from the makefile system for this to work. This is
  12 # because this module use XS extensions (C source code
  13 # 'compiled into' perl)
  14 #
  15 # The lib system below is just so that I (ewan) can test it
  16 # on site...
  17 #
  18
  19 use lib "/nfs/disk100/pubseq/wise/PerlMod/";
  20
  21 #
  22 # This is a simple example script. We are going
  23 # to make 3 sequences directly from memory and
  24 # then align them once using blosum matrix and once
  25 # using a gonnet matrix. These matrices should
  26 # in the examples directory.
  27 #
  28
  29 use Bio::Tools::pSW;
  30
  31 # redundant, as Bio::Tools::pSW uses them, but useful to say
  32 # precisely what we are using ;)
  33
  34 use Bio::Seq;
  35 use Bio::SimpleAlign;
  36 use Bio::AlignIO;
  37
  38 # for legibility - write with newlines and then strip them!
  39
  40 $tseq = 'SKSESPKEPEQLRKLFIGGLSFETTDESLRSHFEQWGTLTDCVVMRDPNTKRSRGFGFVT
  41 YATVEEVDAAMNARPHKVDGRVVEPKRAVSREDSQRPGAHLTVKKIFVGGIKEDTEEHHL
  42 RDYFEQYGKIEVIEIMTDRGSGKKRGFAFVTFDDHDSVDKIVIQKYHTVNGHNCEVRKAL
  43 SKQEMASASSSQRGRSGSGNFGGGRGGGFGGNDNFGRGGNFSGRGGFGGSRGGGGYGGSG
  44 DGYNGFGNDGGYGGGGPGYSGGSRGYGSGGQGYGNQGSGYGGSGSYDSYNNGGGRGFGGG
  45 SGSNFGGGGSYNDFGNYNNQSSNFGPMKGGNFGGRSSGPYGGGGQYFAKPRNQGGYGGSS
  46 SSSSYGSGRRF';
  47
  48 $tseq =~ s/[^A-Z]//g;
  49
  50 $seq1 = Bio::Seq->new(-id=>'roa1_human',-seq=>$tseq);
  51
  52 $tseq = 'MVNSNQNQNGNSNGHDDDFPQDSITEPEHMRKLFIGGLDYRTTDENLKAHFEKWGNIVDV
  53 VVMKDPRTKRSRGFGFITYSHSSMIDEAQKSRPHKIDGRVVEPKRAVPRQDIDSPNAGAT
  54 VKKLFVGALKDDHDEQSIRDYFQHFGNIVDINIVIDKETGKKRGFAFVEFDDYDPVDKVV
  55 LQKQHQLNGKMVDVKKALPKQNDQQGGGGGRGGPGGRAGGNRGNMGGGNYGNQNGGGNWN
  56 NGGNNWGNNRGGNDNWGNNSFGGGGGGGGGYGGGNNSWGNNNPWDNGNGGGNFGGGGNNW
  57 NNGGNDFGGYQQNYGGGPQRGGGNFNNNRMQPYQGGGGFKAGGGNQGNYGGNNQGFNNGG
  58 NNRRY';
  59
  60 $tseq =~ s/[^A-Z]//g;
  61
  62 $seq2 = Bio::Seq->new(-id=>'roa1_drome',-seq=>$tseq);
  63
  64 $tseq = 'MHKSEAPNEPEQLRKLFIGGLSFETTDESLREHFEQWGTLTDCVVMRDPNSKRSRGFGFV
  65 TYLSTDEVDAAMTARPHKVDGRVVEPKRAVSREDSSRPGAHLTVKKIFVGGIKEDTEEDH
  66 LREYFEQYGKIEVIEIMTDRGSGKKRGFAFVTFEDHDSVDKIVIQKYHTVNNHNSQVRKA
  67 LSKQEMASVSGSQRERGGSGNYGSRGGFGNDNFGGRGGNFGGNRGGGGGFGNRGYGGDGY
  68 NGDGQLWWQPSLLGWNRGYGAGQGGGYGAGQGGGYGGGGQGGGYGGNGGYDGYNGGGSGF
  69 SGSGGNFGSSGGYNDFGNYNSQSSSNFGPMKGGNYGGGRNSGPYGGGYGGGSASSSSGYG
  70 GGRRF';
  71
  72 $tseq =~ s/[^A-Z]//g;
  73 $seq3 = Bio::Seq->new(-id=>'roa1_xenla',-seq=>$tseq);
  74
  75
  76 #
  77 # Now make an Alignment Factory with blosum62 as a matrix
  78 # gap -12 and ext -2
  79 #
  80
  81 $fac = Bio::Tools::pSW->new(-matrix => 'blosum62.bla',-gap => 12, -ext => 2);
  82
  83
  84 #
  85 # run seq1 vs seq2 and seq1 vs seq3 and write the output direct
  86 # to stdout using the 'pretty' method
  87 #
  88
  89 $fac->align_and_show($seq1,$seq2,STDOUT);
  90 print "Next alignment\n";
  91 $fac->align_and_show($seq1,$seq3,STDOUT);
  92
  93
  94 #
  95 # a different factory, using gonnet, and now make a simple align and
  96 # provide MSF format
  97 #
  98
  99 $fac = Bio::Tools::pSW->new(-matrix => 'gon250.bla',-gap => 12, -ext => 2);
 100
 101
 102 # switch on reporting this time and change the amount of memory it is allowed
 103
 104 print STDOUT "Doing the next calculation in limited memory, with a progress report\n";
 105
 106 $fac->report(1);
 107 $fac->kbyte(100);
 108
 109 $al = $fac->pairwise_alignment($seq1,$seq2);
 110
 111
 112 # write out a MSF file
 113 my $out = Bio::AlignIO->newFh('-fh'=> \*STDOUT,  '-format' => 'msf');
 114 my $status = print $out $al;
 115 #$al->write_MSF(\*STDOUT);
 116
 117
 118