9.6 bptutorial.pl


I've already shown you a little of the bptutorial.pl document. I ran and discussed a few of the short example programs in the preceding sections.

As you know, one of the easiest ways to get started with a programming system is to find some working and fairly generic programs in that system. You can read and run the programs, and then proceed to alter them using them as templates for your own programming development.

Bioperl comes with a directory of example programs, but the best place to begin looking for starting-off program code is right in the bptutorial.pl document itself. That .pl suffix on the name is the giveaway; the document is actually itself a program, cleverly designed so that you can read, and run, example programs that exercise the core parts of the Bioperl project.

The following explanation of the runnable programs that are part of bptutorial.pl appears at the end of the document (when you view it on the Web or as the output of perldoc bptutorial.pl ).

 V.2 Appendix: Tutorial demo scripts        The following scripts demonstrate many of the features of        bioperl. To run all the core demos, run:         > perl -w  bptutorial.pl 0        To run a subset of the scripts do         > perl -w  bptutorial.pl        and use the displayed help screen.        It may be best to start by just running one or two demos        at a time. For example, to run the basic sequence manipu-        lation demo, do:         > perl -w  bptutorial.pl 1        Some of the later demos require that you have an internet        connection and/or that you have an auxilliary bioperl        library and/or external cpan module and/or external pro-        gram installed.  They may also fail if you are not running        under Linux or Unix.  In all of these cases, the script        should fail "gracefully" simply saying the demo is being        skipped.  However if the script "crashes", simply run the        other demos individually (and perhaps send an email to        bioperl-l@bioperl.org detailing the problem :-). 

(Recall that the -w flag to Perl turns on warnings in almost the same manner as a use warnings; directive.)

To test my Bioperl installation, I started by running the basic sequence manipulation demo as suggested.

First, I thought I might copy the bptutorial.pl program file into my own working directory from the Bioperl distribution directory where I'd unpacked the source code. I wanted to put it in my own directory so as not to muddy up the Bioperl distribution directory with my own extraneous files. However, I discovered that the tutorial demo programs rely on a number of datafiles that are found in the t/data/ subdirectory of the Bioperl distribution. Running the program in my own directory gave me an error because the program evidently requires a FASTA -formatted file called dna1.fa :

 [tisdall]$ perl -w bptutorial.pl 1 Beginning sequence_manipulations and SeqIO example...  ------------- EXCEPTION  ------------- MSG: Could not open t/data/dna1.fa: No such file or directory STACK Bio::Root::IO::_initialize_io /usr/local/lib/perl5/site_perl/5.8.0/Bio/Root/IO. pm:264 STACK Bio::SeqIO::_initialize /usr/local/lib/perl5/site_perl/5.8.0/Bio/SeqIO.pm:437 STACK Bio::SeqIO::fasta::_initialize /usr/local/lib/perl5/site_perl/5.8.0/Bio/SeqIO/ fasta.pm:80 STACK Bio::SeqIO::new /usr/local/lib/perl5/site_perl/5.8.0/Bio/SeqIO.pm:355 STACK Bio::SeqIO::new /usr/local/lib/perl5/site_perl/5.8.0/Bio/SeqIO.pm:368 STACK main::_ _ANON_ _ bptutorial.pl:2758 STACK toplevel bptutorial.pl:3933 -------------------------------------- [tisdall]$ 

I decided it would be easier to run bptutorial.pl from the distribution directory. I entered that directory; on my Linux machine, I unpacked the Bioperl source code into /usr/local/src/bioperl-1.2.1. I then tried to run the demo:

 [tisdall]$ cd /usr/local/src/bioperl-1.2.1 [tisdall]$ perl -w bptutorial.pl 1 Beginning sequence_manipulations and SeqIO example...  First sequence in fasta format...  >Test1 AGCTTTTCATTCTGACTGCAACGGGCAATATGTCTCTGTGTGGATTAAAAAAAGAGTGTC TGATAGCAGCTTCTGAACTGGTTACCTGCCGTGAGTAAATTAAAATTTTATTGACTTAGG TCACTAAATACTTTAACCAATATAGGCATAGCGCACAGACAGATAAAAATTACAGAGTAC ACAACATCCATGAAACGCATTAGCACCACC Seq object display id is Test1 Sequence is  AGCTTTTCATTCTGACTGCAACGGGCAATATGTCTCTGTGTGGATTAAAAAAAGAGTGTCTGATAGCAGCTTCTGAACTGGTTAC CTGCCGTGAGTAAATTAAAATTTTATTGACTTAGGTCACTAAATACTTTAACCAATATAGGCATAGCGCACAGACAGATAAAAATT ACAGAGTACACAACATCCATGAAACGCATTAGCACCACC  Sequence from 5 to 10 is TTTCAT  Acc num is unknown  Moltype is dna  Primary id is Test1  Truncated Seq object sequence is TTTCAT  Reverse complemented sequence 5 to 10  is ATGAAA   Translated sequence 6 to 15 is LQRAICLCVD  Beginning 3-frame and alternate codon translation example...  ctgagaaaataa translated using method defaults   : LRK* ctgagaaaataa translated as a coding region (CDS): MRK Translating in all six frames:  frame: 0 forward: LRK*  frame: 0 reverse-complement: LFSQ  frame: 1 forward: *EN  frame: 1 reverse-complement: YFL  frame: 2 forward: EKI  frame: 2 reverse-complement: IFS Translating with all codon tables using method defaults: 1 : LRK* 2 : L*K* 3 : TRK* 4 : LRK* 5 : LSK* 6 : LRKQ 9 : LSN* 10 : LRK* 11 : LRK* 12 : SRK* 13 : LGK* 14 : LSNY 15 : LRK* 16 : LRK* 21 : LSN* [tisdall]$ 

That seemed to run without error. Now I wanted to see the Perl code that had run and produced that output.

Exploring a bit, I found that the POD documentation part of bptutorial.pl is roughly the first half of the file, and that the second half of the file contains the Perl code for the demos.

The author attributions and the libraries that are loaded are at the beginning of the Perl code section, somewhere near the middle of the bptutorial.pl file:

 #!/usr/bin/perl # PROGRAM  : bptutorial.pl # PURPOSE  : Demonstrate various uses of the bioperl package # AUTHOR   : Peter Schattner schattner@alum.mit.edu # CREATED  : Dec 15 2000 # REVISION : $Id: ch09.xml,v 1.2 2003/10/03 19:16:55 becki Exp chodacki $ use strict; use Bio::SimpleAlign; use Bio::AlignIO; use Bio::SeqIO; use Bio::Seq; 

That seems like a good bet for a list of the most important Bioperl modules. Actually, some other modules are loaded for individual demos; at various points, the following come into play:

 use Bio::SearchIO; use Bio::Root::IO; use Bio::MapIO; use Bio::TreeIO; use Bio::Perl 

Following those use Module directives comes the following:

 # subroutine references my ($access_remote_db, $index_local_db, $fetch_local_db,     $sequence_manipulations, $seqstats_and_seqwords,     $restriction_and_sigcleave, $other_seq_utilities, $run_remoteblast,     $run_standaloneblast,  $blast_parser, $bplite_parsing, $hmmer_parsing,     $run_clustalw_tcoffee, $run_psw_bl2seq, $simplealign,     $gene_prediction_parsing, $sequence_annotation, $largeseqs,     $run_tree, $run_map, $run_struct, $run_perl, $searchio_parsing,     $liveseqs, $demo_variations, $demo_xml, $display_help, $bpinspect1 ); # global variable file names.  Edit these if you want to try #out a tutorial script on a different file  Bio::Root::IO->catfile("t","data","ecolitst.fa"); # used in $sequence_manipulations my $dna_seq_file = Bio::Root::IO->catfile("t","data","dna1.fa");      # used in $other_seq_utilities and in $run_perl and $sequence_annotation my $amino_seq_file = Bio::Root::IO->catfile("t","data","cysprot1.fa");  # used in $blast_parser my $blast_report_file = Bio::Root::IO->catfile("t","data","blast.report");   # used in $bplite_parsing my $bp_parse_file1 = Bio::Root::IO->catfile("t","data","blast.report"); # used in $bplite_parsing my $bp_parse_file2 = Bio::Root::IO->catfile("t","data","psiblastreport.out"); # used in $bplite_parsing my $bp_parse_file3 = Bio::Root::IO->catfile("t","data","bl2seq.out");        # used in $run_clustalw_tcoffee my $unaligned_amino_file = Bio::Root::IO->catfile("t","data","cysprot1a.fa"); # used in $simplealign my $aligned_amino_file = Bio::Root::IO->catfile("t","data","testaln.pfam"); # other global variables my (@runlist, $n ); 

A look at the documentation ( perldoc Bio::Root::IO ) shows that the catfile method returns the pathname of a file, which is being saved in a scalar variable such as $dna_seq_file . This method is there for portability between operating systems, because operating systems each have their own syntax for specifying pathnames.

After that comes the code for the help screen, the output of which you've already seen; next comes the individual demo subroutines; and finally the code for the main subroutine, which you saw in the last section.

At the very end of the bptutorial.pl file, I found the part of the code that launches all the demos:

 ## "main" program follows #no strict 'refs';     @runlist = @ARGV;     # display help if no option     if (scalar(@runlist)=  =0) {&$display_help;};     # argument = 0 means run tests 1 thru 22     if ($runlist[0] =  = 0) {@runlist = (1..22); };     foreach $n  (@runlist) {         if ($n =  =100) {my $object = $runlist[1]; &$bpinspect1($object); last;}         if ($n =  =1) {&$sequence_manipulations; next;}         if ($n =  =2) {&$seqstats_and_seqwords; next;}         if ($n =  =3) {&$restriction_and_sigcleave; next;}         if ($n =  =4) {&$other_seq_utilities; next;}         if ($n =  =5) {&$run_perl; next;}         if ($n =  =6) {&$searchio_parsing; next;}         if ($n =  =7) {&$bplite_parsing; next;}         if ($n =  =8) {&$hmmer_parsing; next;}         if ($n =  =9) {&$simplealign ; next;}         if ($n =  =10) {&$gene_prediction_parsing; next;}         if ($n =  =11) {&$access_remote_db; next;}         if ($n =  =12) {&$index_local_db; next;}         if ($n =  =13) {&$fetch_local_db; next;}         if ($n =  =14) {&$sequence_annotation; next;}         if ($n =  =15) {&$largeseqs; next;}         if ($n =  =16) {&$liveseqs; next;}         if ($n =  =17) {&$run_struct; next;}         if ($n =  =18) {&$demo_variations; next;}         if ($n =  =19) {&$demo_xml; next;}         if ($n =  =20) {&$run_tree; next;}         if ($n =  =21) {&$run_map; next;}         if ($n =  =22) {&$run_remoteblast; next;}         if ($n =  =23) {&$run_standaloneblast; next;}         if ($n =  =24) {&$run_clustalw_tcoffee; next;}         if ($n =  =25) {&$run_psw_bl2seq; next;}          &$display_help;     } ## End of "main" 

So, I searched for sequence_manipulation and found the following code for this, the first of the bptutorial.pl demos:

 ################################################# # sequence_manipulations  (  ): # $sequence_manipulations = sub {     my ($infile, $in, $out, $seqobj);     $infile = $dna_seq_file;     print "\nBeginning sequence_manipulations and SeqIO example... \n";     # III.3.1 Transforming sequence files (SeqIO)     $in  = Bio::SeqIO->new('-file' => $infile ,                            '-format' => 'Fasta');     $seqobj = $in->next_seq(  );     # perl "tied filehandle" syntax is available to SeqIO,     # allowing you to use the standard <> and print operations     # to read and write sequence objects, eg:     #$out = Bio::SeqIO->newFh('-format' => 'EMBL');     $out = Bio::SeqIO->newFh('-format' => 'fasta');     print "First sequence in fasta format... \n";     print $out $seqobj;     # III.4 Manipulating individual sequences     # The following methods return strings     print "Seq object display id is ",     $seqobj->display_id(  ), "\n"; # the human read-able id of the sequence     print "Sequence is ",     $seqobj->seq(  )," \n";        # string of sequence     print "Sequence from 5 to 10 is ",     $seqobj->subseq(5,10)," \n"; # part of the sequence as a string     print "Acc num is ",     $seqobj->accession_number(  ), " \n"; # when there, the accession number     print "Moltype is ",     $seqobj->alphabet(  ), " \n";    # one of 'dna','rna','protein'     print "Primary id is ", $seqobj->primary_seq->primary_id(  )," \n";     # a unique id for this sequence irregardless     #print "Primary id is ", $seqobj->primary_id(  ), " \n";     # a unique id for this sequence irregardless     # of its display_id or accession number     # The following methods return an array of  Bio::SeqFeature objects     $seqobj->top_SeqFeatures; # The 'top level' sequence features     $seqobj->all_SeqFeatures; # All sequence features, including sub     # seq features     # The following methods returns new sequence objects,     # but do not transfer features across     # truncation from 5 to 10 as new object     print "Truncated Seq object sequence is ",     $seqobj->trunc(5,10)->seq(  ), " \n";     # reverse complements sequence     print "Reverse complemented sequence 5 to 10  is ",     $seqobj->trunc(5,10)->revcom->seq, "  \n";     # translation of the sequence     print "Translated sequence 6 to 15 is ",     $seqobj->translate->subseq(6,15), " \n";     my $c = shift;     $c = 'ctgagaaaataa';     print "\nBeginning 3-frame and alternate codon translation example... \n";     my $seq = new Bio::PrimarySeq('-SEQ' => $c,                                   '-ID' => 'no.One');     print "$c translated using method defaults   : ",     $seq->translate->seq, "\n";     # Bio::Seq uses same sequence methods as PrimarySeq     my $seq2 = new Bio::Seq('-SEQ' => $c, '-ID' => 'no.Two');     print "$c translated as a coding region (CDS): ",     $seq2->translate(undef, undef, undef, undef, 1)->seq, "\n";     print "\nTranslating in all six frames:\n";     my @frames = (0, 1, 2);     foreach my $frame (@frames) {         print  " frame: ", $frame, " forward: ",         $seq->translate(undef, undef, $frame)->seq, "\n";         print  " frame: ", $frame, " reverse-complement: ",         $seq->revcom->translate(undef, undef, $frame)->seq, "\n";     }     print "Translating with all codon tables using method defaults:\n";     my @codontables = qw( 1 2 3 4 5 6 9 10 11 12 13 14 15 16 21 );     foreach my $ct (@codontables) {         print $ct, " : ",         $seq->translate(undef, undef, undef, $ct)->seq, "\n";     }     return 1; } ; ################################################# 


Mastering Perl for Bioinformatics
Mastering Perl for Bioinformatics
ISBN: 0596003072
EAN: 2147483647
Year: 2003
Pages: 156

flylib.com © 2008-2017.
If you may any questions please contact us: flylib@qtcs.net