lib/Bio/AlignIO/metafasta.pm

   1 #
   2 # BioPerl module for Bio::AlignIO::metafasta
   3 #
   4 # Copyright Heikki Lehvaslaiho
   5 #
   6 # You may distribute this module under the same terms as perl itself
   7 #
   8 # POD documentation - main docs before the code
   9
  10 =head1 NAME
  11
  12 Bio::AlignIO::metafasta - Metafasta MSA Sequence input/output stream
  13
  14 =head1 SYNOPSIS
  15
  16 Do not use this module directly.  Use it via the L<Bio::AlignIO> class.
  17
  18 =head1 DESCRIPTION
  19
  20 This object can transform L<Bio::SimpleAlign> objects to and from
  21 metafasta flat file databases.
  22
  23 The format of a metafasta file is
  24
  25   >test/1-25
  26   ABCDEFHIJKLMNOPQRSTUVWXYZ
  27   &charge
  28   NBNAANCNJCNNNONNCNNUNNXNZ
  29   &chemical
  30   LBSAARCLJCLSMOIMCHHULRXRZ
  31
  32 where the sequence block is followed by one or several meta blocks.
  33 Each meta block starts with the ampersand character '&' in the first
  34 column and is immediately followed by the name of the meta data which
  35 continues until the new line. The meta data follows it. All
  36 characters, except new line, are important in meta data.
  37
  38 =head1 SEE ALSO
  39
  40 L<Bio::SeqIO::metafasta>
  41
  42 =head1 FEEDBACK
  43
  44 =head2 Support
  45
  46 Please direct usage questions or support issues to the mailing list:
  47
  48 I<bioperl-l@bioperl.org>
  49
  50 rather than to the module maintainer directly. Many experienced and
  51 reponsive experts will be able look at the problem and quickly
  52 address it. Please include a thorough description of the problem
  53 with code and data examples if at all possible.
  54
  55 =head2 Reporting Bugs
  56
  57 Report bugs to the Bioperl bug tracking system to help us keep track
  58 the bugs and their resolution.  Bug reports can be submitted via the
  59 web:
  60
  61   https://github.com/bioperl/bioperl-live/issues
  62
  63 =head1 AUTHOR - Heikki Lehvaslaiho
  64
  65 Email heikki-at-bioperl-dot-org
  66
  67 =head1 APPENDIX
  68
  69 The rest of the documentation details each of the object
  70 methods. Internal methods are usually preceded with a _
  71
  72 =cut
  73
  74 # Let the code begin...
  75
  76 package Bio::AlignIO::metafasta;
  77
  78 use vars qw($WIDTH);
  79 use strict;
  80
  81 use Bio::SimpleAlign;
  82 use Bio::Seq::Meta;
  83 use Bio::Seq::SeqFactory;
  84 use Bio::Seq::SeqFastaSpeedFactory;
  85
  86 use base qw(Bio::AlignIO);
  87
  88 BEGIN { $WIDTH = 60}
  89
  90 sub _initialize {
  91   my($self,@args) = @_;
  92   $self->SUPER::_initialize(@args);
  93   my ($width) = $self->_rearrange([qw(WIDTH)], @args);
  94   $width && $self->width($width);
  95 }
  96
  97 =head2 next_aln
  98
  99  Title   : next_aln
 100  Usage   : $aln = $stream->next_aln()
 101  Function: returns the next alignment in the stream.
 102  Returns : L<Bio::Align::AlignI> object - returns 0 on end of file
 103             or on error
 104  Args    : NONE
 105
 106 =cut
 107
 108 sub next_aln {
 109     my( $self ) = @_;
 110     my $seq;
 111     my $alphabet;
 112     local $/ = "\n>";
 113
 114     my $aln =  Bio::SimpleAlign->new();
 115
 116     while(defined (my $entry = $self->_readline)) {
 117         chomp($entry);
 118         if ($entry =~ m/\A\s*\Z/s) { # very first one
 119             return unless $entry = $self->_readline;
 120             chomp($entry);
 121         }
 122         $entry =~ s/^>//;
 123
 124         my ($top,$sequence) = split(/\n/,$entry,2);
 125         defined $sequence && $sequence =~ s/>//g;
 126
 127         my @metas;
 128         ($sequence, @metas) = split /\n&/, $sequence;
 129
 130         my ($id, $start, $end);
 131         if ( $top =~ /(\S+)\/(\d+)-(\d+)/ ) {
 132             $id = $1;
 133             $start = $2;
 134             $end = $3;
 135         }
 136         elsif ($top =~ /(\S+)/) {
 137             $id = $1;
 138             $start = 1;
 139             $end = length($sequence);
 140         }
 141
 142         defined $sequence && $sequence =~ s/\s//g; # Remove whitespace
 143
 144         $seq = Bio::Seq::Meta->new('-seq'        => $sequence,
 145                                    '-display_id' => $id,
 146                                    '-start'      => $start,
 147                                    '-end'        => $end,
 148                                    '-alphabet'   => $self->alphabet,
 149                                    );
 150
 151         foreach my $meta (@metas) {
 152             my ($name,$string) = split /\n/, $meta;
 153             $string =~ s/\n//g; # Remove newlines, spaces are important
 154             $seq->named_meta($name, $string);
 155         }
 156
 157         $aln->add_seq($seq);
 158
 159         # alignment needs seqs all the same length, pad with gaps
 160         my $alnlen = $aln->length;
 161         foreach my $seq ( $aln->each_seq ) {
 162                 if ( $seq->length < $alnlen ) {
 163                         my ($diff) = ($alnlen - $seq->length);
 164                         $seq->seq( $seq->seq() . "-" x $diff);
 165                 }
 166         }
 167     }
 168     return $aln if $aln->num_sequences;
 169         return;
 170 }
 171
 172 =head2 write_aln
 173
 174  Title   : write_aln
 175  Usage   : $stream->write_aln(@aln)
 176  Function: writes the $aln object into the stream in fasta format
 177  Returns : 1 for success and 0 for error
 178  Args    : L<Bio::Align::AlignI> object
 179
 180 =cut
 181
 182 sub write_aln {
 183     my ($self,@aln) = @_;
 184     my $width = $self->width;
 185
 186     foreach my $aln (@aln) {
 187         if( ! $aln || ! $aln->isa('Bio::Align::AlignI')  ) {
 188             $self->warn("Must provide a Bio::Align::AlignI object when calling write_aln");
 189             next;
 190         }
 191         foreach my $seq ( $aln->each_seq() ) {
 192             my $name = $aln->displayname($seq->get_nse);
 193
 194             my $str  = $seq->seq();
 195             if(length($str) > 0) {
 196                 $str =~ s/(.{1,$width})/$1\n/g;
 197             } else {
 198                 $str = "\n";
 199             }
 200             $self->_print (">",$name,"\n",$str) or return;
 201             if ($seq->isa('Bio::Seq::MetaI')) {
 202                 foreach my $meta ($seq->meta_names) {
 203                     my $str = $seq->named_meta($meta);
 204                     $str =~ s/(.{1,$width})/$1\n/g;
 205                     $self->_print ("&",$meta,"\n",$str);
 206                 }
 207             }
 208         }
 209     }
 210     $self->flush if $self->_flush_on_write && defined $self->_fh;
 211     return 1;
 212 }
 213
 214
 215 =head2 width
 216
 217  Title   : width
 218  Usage   : $obj->width($newval)
 219  Function: Get/Set the line width for METAFASTA output
 220  Returns : value of width
 221  Args    : newvalue (optional)
 222
 223
 224 =cut
 225
 226 sub width{
 227    my ($self,$value) = @_;
 228    if( defined $value) {
 229       $self->{'width'} = $value;
 230     }
 231     return $self->{'width'} || $WIDTH;
 232 }
 233
 234 1;