File Coverage

Bio/AlignIO/xmfa.pm
Criterion Covered Total %
statement 66 75 88.0
branch 27 38 71.0
condition 11 28 39.2
subroutine 6 7 85.7
pod 3 3 100.0
total 113 151 74.8


line stmt bran cond sub pod time code
1             #
2             # BioPerl module for Bio::AlignIO::xmfa
3             #
4             # Copyright Chris Fields
5             #
6             # You may distribute this module under the same terms as perl itself
7             # POD documentation - main docs before the code
8              
9             =head1 NAME
10              
11             Bio::AlignIO::xmfa - XMFA MSA Sequence input/output stream
12              
13             =head1 SYNOPSIS
14              
15             Do not use this module directly. Use it via the L
16             class.
17              
18             =head1 DESCRIPTION
19              
20             This object can transform L objects from
21             XMFA flat file databases. For more information, see:
22              
23             http://asap.ahabs.wisc.edu/mauve-aligner/mauve-user-guide/mauve-output-file-formats.html
24              
25             This module is based on the AlignIO::fasta parser written by
26             Peter Schattner
27              
28             =head1 TODO
29              
30             Finish write_aln(), clean up code, allow LargeLocatableSeq (ie for
31             very large sequences a'la Mauve)
32              
33             =head1 FEEDBACK
34              
35             =head2 Support
36              
37             Please direct usage questions or support issues to the mailing list:
38              
39             I
40              
41             rather than to the module maintainer directly. Many experienced and
42             reponsive experts will be able look at the problem and quickly
43             address it. Please include a thorough description of the problem
44             with code and data examples if at all possible.
45              
46             =head2 Reporting Bugs
47              
48             Report bugs to the Bioperl bug tracking system to help us keep track
49             the bugs and their resolution. Bug reports can be submitted via the
50             web:
51              
52             https://github.com/bioperl/bioperl-live/issues
53              
54             =head1 AUTHORS
55              
56             Chris Fields
57              
58             =head1 APPENDIX
59              
60             The rest of the documentation details each of the object
61             methods. Internal methods are usually preceded with a _
62              
63             =cut
64              
65             # Let the code begin...
66              
67             package Bio::AlignIO::xmfa;
68 2     2   415 use strict;
  2         3  
  2         53  
69              
70 2     2   6 use base qw(Bio::AlignIO);
  2         2  
  2         1009  
71             our $WIDTH = 60;
72              
73             =head2 next_aln
74              
75             Title : next_aln
76             Usage : $aln = $stream->next_aln
77             Function: returns the next alignment in the stream.
78             Returns : Bio::Align::AlignI object - returns 0 on end of file
79             or on error
80             Args : -width => optional argument to specify the width sequence
81             will be written (60 chars by default)
82              
83             See L
84              
85             =cut
86              
87             sub next_aln {
88 3     3 1 10 my $self = shift;
89 3         13 my ($width) = $self->_rearrange([qw(WIDTH)],@_);
90 3   33     15 $self->width($width || $WIDTH);
91              
92 3         5 my ($name, $tempname, $seqchar);
93 3         15 my $aln = Bio::SimpleAlign->new();
94 3         3 my $seqs = 0;
95             # alignments
96 3         16 while (defined (my $entry = $self->_readline) ) {
97 75         62 chomp $entry;
98 75 100       144 if ( index($entry, '=') == 0 ) {
    100          
99 3 50 33     13 if (defined $name && $seqchar) {
100 3         6 my $seq = $self->_process_seq($name, $seqchar);
101 3         9 $aln->add_seq($seq);
102             }
103 3 50 33     25 if ($aln && $entry =~ m{score\s*=\s*(\d+)}) {
104 3         11 $aln->score($1);
105             }
106 3         3 $seqchar = '';
107 3         4 undef $name;
108 3         4 last;
109             } elsif ( $entry =~ m{^>.+$}xms) {
110 9 100       18 if ( defined $name ) {
111 6         14 my $seq = $self->_process_seq($name, $seqchar);
112 6         13 $aln->add_seq($seq);
113             }
114 9         6 $seqchar = '';
115 9         22 $name = $entry;
116             } else {
117 63         122 $seqchar .= $entry;
118             }
119             }
120            
121             # this catches last sequence if '=' is not present (Mauve)
122 3 50       6 if ( defined $name ) {
123 0         0 my $seq = $self->_process_seq($name, $seqchar);
124 0         0 $aln->add_seq($seq);
125             }
126 3 50       8 $aln->num_sequences ? return $aln : return;
127             }
128              
129             =head2 write_aln
130              
131             Title : write_aln
132             Usage : $stream->write_aln(@aln)
133             Function: writes the $aln object into the stream in xmfa format
134             Returns : 1 for success and 0 for error
135             Args : L object
136              
137             See L
138              
139             =cut
140              
141             sub write_aln {
142 2     2 1 9 my ($self,@aln) = @_;
143 2         6 my $width = $self->width;
144 2         5 my ($seq,$desc,$rseq,$name,$count,$length,$seqsub,$start,$end,$strand,$id);
145              
146 2         5 foreach my $aln (@aln) {
147 2 50 33     14 if( ! $aln || ! $aln->isa('Bio::Align::AlignI') ) {
148 0         0 $self->warn("Must provide a Bio::Align::AlignI object when calling write_aln");
149 0         0 next;
150             }
151             #if( $self->force_displayname_flat ) {
152             # $aln->set_displayname_flat(1);
153             #}
154 2         4 my $seqct = 1;
155 2         6 foreach $rseq ( $aln->each_seq() ) {
156 9   100     19 ($start, $end, $strand, $id) = ($rseq->start, $rseq->end, $rseq->strand || 0,
157             $rseq->display_id);
158 9 100       29 $strand = ($strand == 1) ? '+' :
    100          
159             ($strand == -1) ? '-' :
160             '';
161 9         44 $name = sprintf("%d:%d-%d %s %s",$seqct,$start,$end,$strand,$id);
162 9         14 $seq = $rseq->seq();
163 9   50     15 $desc = $rseq->description || '';
164 9 50       61 $self->_print (">$name $desc\n") or return ;
165 9         11 $count = 0;
166 9         8 $length = length($seq);
167 9 50 33     34 if(defined $seq && $length > 0) {
168 9         95 $seq =~ s/(.{1,$width})/$1\n/g;
169             } else {
170 0         0 $seq = "\n";
171             }
172 9 50       19 $self->_print($seq) || return 0;
173 9         11 $seqct++;
174             }
175 2         4 my $alndesc = '';
176 2 100       9 $alndesc = "score = ".$aln->score if ($aln->score);
177 2 50       8 $self->_print("= $alndesc\n") || return 0;
178            
179             }
180 2 50 33     6 $self->flush if $self->_flush_on_write && defined $self->_fh;
181 2         7 return 1;
182             }
183              
184             =head2 _get_len
185              
186             Title : _get_len
187             Usage :
188             Function: determine number of alphabetic chars
189             Returns : integer
190             Args : sequence string
191              
192             =cut
193              
194             sub _get_len {
195 0     0   0 my ($self,$seq) = @_;
196 0         0 $seq =~ s/[^A-Z]//gi;
197 0         0 return CORE::length($seq);
198             }
199              
200             =head2 width
201              
202             Title : width
203             Usage : $obj->width($newwidth)
204             $width = $obj->width;
205             Function: Get/set width of alignment
206             Returns : integer value of width
207             Args : on set, new value (a scalar or undef, optional)
208              
209              
210             =cut
211              
212             sub width{
213 5     5 1 7 my $self = shift;
214              
215 5 100       14 return $self->{'_width'} = shift if @_;
216 2   33     8 return $self->{'_width'} || $WIDTH;
217             }
218              
219             ####### PRIVATE #######
220              
221             sub _process_seq {
222 9     9   13 my ($self, $entry, $seq) = @_;
223 9         11 my ($start, $end, $strand, $seqname, $desc, $all);
224             # put away last name and sequence
225 9 50       52 if ( $entry =~ m{^>\s*\d+:(\d+)-(\d+)\s([+-]{1})(?:\s+(\S+)\s*(\S\.*)?)?} ) {
226 9         28 ($start, $end, $seqname, $desc) = ($1, $2, $4, $5);
227 9 100       27 $strand = ($3 eq '+') ? 1 : -1;
228             } else {
229 0         0 $self->throw("Line does not comform to XMFA format:\n$entry");
230             }
231 9   33     56 my $seqobj = Bio::LocatableSeq->new(
232             -nowarnonempty => 1,
233             -strand => $strand,
234             -seq => $seq,
235             -display_id => $seqname,
236             -description => $desc || $all,
237             -start => $start,
238             -end => $end,
239             -alphabet => $self->alphabet,
240             );
241 9         37 $self->debug("Reading $seqname\n");
242 9         12 return $seqobj;
243             }
244              
245             1;