manipulate_seq.pl

1	#!/usr/local/bin/perl -w
2	
3	# Use BioPerl to convert between sequence formats
4	# and do some basic sequence manipulation
5	# WI Biocomputing course - Bioinformatics for Biologists - October 2003
6	
7	use Bio::SeqIO;
8	use Bio::Tools::SeqStats;   # used to calculate MW
9	
10	$inFile = "BMP7.tfa";
11	
12	$in  = Bio::SeqIO->new('-file' => "$inFile" ,
13	                           '-format' => 'Fasta');
14	$out = Bio::SeqIO->new('-format' => 'Genbank');
15	
16	while ($seqobj = $in->next_seq())
17	{
18	   # Print sequence in output format
19	   $out->write_seq($seqobj);
20	
21	   print "\nRaw sequence:\n", $seqobj->seq(), "\n";
22	   print "\nSequence ID is ", $seqobj->display_id(), "\n";
23	   print "Sequence from 1 to 100:\n", $seqobj->subseq(1,100), "\n";
24	
25	   # Check if there's an accession number
26	   print "Sequence description: ", $seqobj->desc(), "\n";
27	
28	   print "Type of sequence: ", $type = $seqobj->alphabet(), "\n";
29	   if ($type eq "dna")
30	   {
31	      $rev_comp = $seqobj->revcom->seq();
32	      print "\nReverse complemented sequence:\n$rev_comp\n";
33	
34	      print "\nReverse complemented sequence from 1 to 100:\n",
35	          $seqobj->revcom->subseq(1, 100), "\n";
36	   }
37	   $seq_stats  =  Bio::Tools::SeqStats->new($seqobj);
38	   $weight = $seq_stats->get_mol_wt();
39	
40	   print "Molecular weight: $$weight[0]\n";
41	   # Need $$weight[0] because $weight is an array reference
42	}
43	
44