Bio/AlignIO/psi.pm

   1 #
   2 # BioPerl module for Bio::AlignIO::psi
   3 #
   4 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   5 #
   6 # Cared for by Jason Stajich <jason@bioperl.org>
   7 #
   8 # Copyright Jason Stajich
   9 #
  10 # You may distribute this module under the same terms as perl itself
  11
  12 # POD documentation - main docs before the code
  13
  14 =head1 NAME
  15
  16 Bio::AlignIO::psi - Read/Write PSI-BLAST profile alignment files
  17
  18 =head1 SYNOPSIS
  19
  20 This module will parse PSI-BLAST output of the format seqid XXXX
  21
  22 =head1 DESCRIPTION
  23
  24 This is a parser for psi-blast blocks.
  25
  26 =head1 FEEDBACK
  27
  28 =head2 Mailing Lists
  29
  30 User feedback is an integral part of the evolution of this and other
  31 Bioperl modules. Send your comments and suggestions preferably to
  32 the Bioperl mailing list.  Your participation is much appreciated.
  33
  34   bioperl-l@bioperl.org                  - General discussion
  35   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  36
  37 =head2 Support
  38
  39 Please direct usage questions or support issues to the mailing list:
  40
  41 I<bioperl-l@bioperl.org>
  42
  43 rather than to the module maintainer directly. Many experienced and
  44 reponsive experts will be able look at the problem and quickly
  45 address it. Please include a thorough description of the problem
  46 with code and data examples if at all possible.
  47
  48 =head2 Reporting Bugs
  49
  50 Report bugs to the Bioperl bug tracking system to help us keep track
  51 of the bugs and their resolution. Bug reports can be submitted via
  52 the web:
  53
  54   https://github.com/bioperl/bioperl-live/issues
  55
  56 =head1 AUTHOR - Jason Stajich
  57
  58 Email jason@bioperl.org
  59
  60 =head1 APPENDIX
  61
  62 The rest of the documentation details each of the object methods.
  63 Internal methods are usually preceded with a _
  64
  65 =cut
  66
  67
  68 # Let the code begin...
  69
  70
  71 package Bio::AlignIO::psi;
  72 use vars qw($BlockLen $IdLength);
  73 use strict;
  74
  75 $BlockLen = 100;
  76 $IdLength = 13;
  77
  78 # Object preamble - inherits from Bio::Root::Root
  79
  80 use Bio::SimpleAlign;
  81 use Bio::LocatableSeq;
  82
  83 use base qw(Bio::AlignIO);
  84
  85 =head2 new
  86
  87  Title   : new
  88  Usage   : my $obj = Bio::AlignIO::psi->new();
  89  Function: Builds a new Bio::AlignIO::psi object
  90  Returns : Bio::AlignIO::psi
  91  Args    :
  92
  93 =cut
  94
  95 =head2 next_aln
  96
  97  Title   : next_aln
  98  Usage   : $aln = $stream->next_aln()
  99  Function: returns the next alignment in the stream
 100  Returns : Bio::Align::AlignI object
 101  Args    : NONE
 102
 103 See L<Bio::Align::AlignI>
 104
 105 =cut
 106
 107 sub next_aln {
 108     my ($self) = @_;
 109     my $aln;
 110     my %seqs;
 111     my @order;
 112     while( defined ($_ = $self->_readline ) ) {
 113         next if( /^\s+$/);
 114         if( !defined $aln ) {
 115             $aln = Bio::SimpleAlign->new();
 116         }
 117         my ($id,$s) = split;
 118         push @order, $id if( ! defined $seqs{$id});
 119         $seqs{$id} .= $s;
 120     }
 121     foreach my $id ( @order) {
 122     my $gaps = $seqs{$id} =~ tr/-/-/;
 123         my $seq = Bio::LocatableSeq->new(-seq => $seqs{$id},
 124                                         -id  => $id,
 125                                         -start => 1,
 126                                         -end   => length($seqs{$id}) - $gaps,
 127                                          -alphabet => $self->alphabet,
 128                     );
 129         $aln->add_seq($seq);
 130     }
 131     return $aln if defined $aln && $aln->num_sequences;
 132         return;
 133 }
 134
 135 =head2 write_aln
 136
 137  Title   : write_aln
 138  Usage   : $stream->write_aln(@aln)
 139  Function: writes the NCBI psi-format object (.aln) into the stream
 140  Returns : 1 for success and 0 for error
 141  Args    : Bio::Align::AlignI object
 142
 143 L<Bio::Align::AlignI>
 144
 145 =cut
 146
 147 sub write_aln {
 148         my ($self,$aln) = @_;
 149         unless( defined $aln && ref($aln) &&
 150                           $aln->isa('Bio::Align::AlignI') ) {
 151                 $self->warn("Must provide a valid Bio::Align::AlignI to write_aln");
 152                 return 0;
 153         }
 154         my $ct = 0;
 155         my @seqs = $aln->each_seq;
 156         my $len = 1;
 157         my $alnlen = $aln->length;
 158         my $idlen = $IdLength;
 159         my @ids = map { substr($_->display_id,0,$idlen) } @seqs;
 160         while( $len < ($alnlen + 1) ) {
 161                 my $start = $len;
 162                 my $end   = $len + $BlockLen;
 163                 $end = $alnlen if ( $end > $alnlen );
 164                 my $c = 0;
 165                 foreach my $seq ( @seqs ) {
 166                         $self->_print(sprintf("%-".$idlen."s %s\n",
 167                                                                                  $ids[$c++],
 168                                                                                  $seq->subseq($start,$end)));
 169                 }
 170                 $self->_print("\n");
 171                 $len += $BlockLen+1;
 172         }
 173         $self->flush if $self->_flush_on_write && defined $self->_fh;
 174         return 1;
 175 }
 176
 177 1;