File Coverage

blib/lib/BioX/Seq.pm
Criterion Covered Total %
statement 111 111 100.0
branch 45 46 97.8
condition 28 28 100.0
subroutine 19 19 100.0
pod 11 11 100.0
total 214 215 99.5


line stmt bran cond sub pod time code
1             package BioX::Seq 0.008006;
2              
3 3     3   231835 use 5.016;
  3         27  
4 3     3   16 use strict;
  3         4  
  3         55  
5 3     3   11 use warnings;
  3         5  
  3         198  
6              
7             use overload
8             '""' => \&_stringify,
9             '.=' => \&_concat,
10 3     3   3018 'bool' => sub{return(1)};
  3     101   2542  
  3         27  
  101         275  
11              
12             my %genetic_code = (
13            
14             # Standard codon table
15              
16             TTT => 'F' , TCT => 'S' , TAT => 'Y' , TGT => 'C' ,
17             TTC => 'F' , TCC => 'S' , TAC => 'Y' , TGC => 'C' ,
18             TTA => 'L' , TCA => 'S' , TAA => '*' , TGA => '*' ,
19             TTG => 'L' , TCG => 'S' , TAG => '*' , TGG => 'W' ,
20            
21             CTT => 'L' , CCT => 'P' , CAT => 'H' , CGT => 'R' ,
22             CTC => 'L' , CCC => 'P' , CAC => 'H' , CGC => 'R' ,
23             CTA => 'L' , CCA => 'P' , CAA => 'Q' , CGA => 'R' ,
24             CTG => 'L' , CCG => 'P' , CAG => 'Q' , CGG => 'R' ,
25              
26             ATT => 'I' , ACT => 'T' , AAT => 'N' , AGT => 'S' ,
27             ATC => 'I' , ACC => 'T' , AAC => 'N' , AGC => 'S' ,
28             ATA => 'I' , ACA => 'T' , AAA => 'K' , AGA => 'R' ,
29             ATG => 'M' , ACG => 'T' , AAG => 'K' , AGG => 'R' ,
30              
31             GTT => 'V' , GCT => 'A' , GAT => 'D' , GGT => 'G' ,
32             GTC => 'V' , GCC => 'A' , GAC => 'D' , GGC => 'G' ,
33             GTA => 'V' , GCA => 'A' , GAA => 'E' , GGA => 'G' ,
34             GTG => 'V' , GCG => 'A' , GAG => 'E' , GGG => 'G' ,
35            
36             # Extension containing all permutations using ambiguity codes
37            
38             AAU => 'N' , AAR => 'K' , AAY => 'N' , ACU => 'T' ,
39             ACM => 'T' , ACR => 'T' , ACW => 'T' , ACS => 'T' ,
40             ACY => 'T' , ACK => 'T' , ACV => 'T' , ACH => 'T' ,
41             ACD => 'T' , ACB => 'T' , ACN => 'T' , AGU => 'S' ,
42             AGR => 'R' , AGY => 'S' , ATU => 'I' , ATM => 'I' ,
43             ATW => 'I' , ATY => 'I' , ATH => 'I' , AUA => 'I' ,
44             AUC => 'I' , AUG => 'M' , AUT => 'I' , AUU => 'I' ,
45             AUM => 'I' , AUW => 'I' , AUY => 'I' , AUH => 'I' ,
46             CAU => 'H' , CAR => 'Q' , CAY => 'H' , CCU => 'P' ,
47             CCM => 'P' , CCR => 'P' , CCW => 'P' , CCS => 'P' ,
48             CCY => 'P' , CCK => 'P' , CCV => 'P' , CCH => 'P' ,
49             CCD => 'P' , CCB => 'P' , CCN => 'P' , CGU => 'R' ,
50             CGM => 'R' , CGR => 'R' , CGW => 'R' , CGS => 'R' ,
51             CGY => 'R' , CGK => 'R' , CGV => 'R' , CGH => 'R' ,
52             CGD => 'R' , CGB => 'R' , CGN => 'R' , CTU => 'L' ,
53             CTM => 'L' , CTR => 'L' , CTW => 'L' , CTS => 'L' ,
54             CTY => 'L' , CTK => 'L' , CTV => 'L' , CTH => 'L' ,
55             CTD => 'L' , CTB => 'L' , CTN => 'L' , CUA => 'L' ,
56             CUC => 'L' , CUG => 'L' , CUT => 'L' , CUU => 'L' ,
57             CUM => 'L' , CUR => 'L' , CUW => 'L' , CUS => 'L' ,
58             CUY => 'L' , CUK => 'L' , CUV => 'L' , CUH => 'L' ,
59             CUD => 'L' , CUB => 'L' , CUN => 'L' , GAU => 'D' ,
60             GAR => 'E' , GAY => 'D' , GCU => 'A' , GCM => 'A' ,
61             GCR => 'A' , GCW => 'A' , GCS => 'A' , GCY => 'A' ,
62             GCK => 'A' , GCV => 'A' , GCH => 'A' , GCD => 'A' ,
63             GCB => 'A' , GCN => 'A' , GGU => 'G' , GGM => 'G' ,
64             GGR => 'G' , GGW => 'G' , GGS => 'G' , GGY => 'G' ,
65             GGK => 'G' , GGV => 'G' , GGH => 'G' , GGD => 'G' ,
66             GGB => 'G' , GGN => 'G' , GTU => 'V' , GTM => 'V' ,
67             GTR => 'V' , GTW => 'V' , GTS => 'V' , GTY => 'V' ,
68             GTK => 'V' , GTV => 'V' , GTH => 'V' , GTD => 'V' ,
69             GTB => 'V' , GTN => 'V' , GUA => 'V' , GUC => 'V' ,
70             GUG => 'V' , GUT => 'V' , GUU => 'V' , GUM => 'V' ,
71             GUR => 'V' , GUW => 'V' , GUS => 'V' , GUY => 'V' ,
72             GUK => 'V' , GUV => 'V' , GUH => 'V' , GUD => 'V' ,
73             GUB => 'V' , GUN => 'V' , TAU => 'Y' , TAR => '*' ,
74             TAY => 'Y' , TCU => 'S' , TCM => 'S' , TCR => 'S' ,
75             TCW => 'S' , TCS => 'S' , TCY => 'S' , TCK => 'S' ,
76             TCV => 'S' , TCH => 'S' , TCD => 'S' , TCB => 'S' ,
77             TCN => 'S' , TGU => 'C' , TGY => 'C' , TTU => 'F' ,
78             TTR => 'L' , TTY => 'F' , TUA => 'L' , TUC => 'F' ,
79             TUG => 'L' , TUT => 'F' , TUU => 'F' , TUR => 'L' ,
80             TUY => 'F' , TRA => '*' , UAA => '*' , UAC => 'Y' ,
81             UAG => '*' , UAT => 'Y' , UAU => 'Y' , UAR => '*' ,
82             UAY => 'Y' , UCA => 'S' , UCC => 'S' , UCG => 'S' ,
83             UCT => 'S' , UCU => 'S' , UCM => 'S' , UCR => 'S' ,
84             UCW => 'S' , UCS => 'S' , UCY => 'S' , UCK => 'S' ,
85             UCV => 'S' , UCH => 'S' , UCD => 'S' , UCB => 'S' ,
86             UCN => 'S' , UGA => '*' , UGC => 'C' , UGG => 'W' ,
87             UGT => 'C' , UGU => 'C' , UGY => 'C' , UTA => 'L' ,
88             UTC => 'F' , UTG => 'L' , UTT => 'F' , UTU => 'F' ,
89             UTR => 'L' , UTY => 'F' , UUA => 'L' , UUC => 'F' ,
90             UUG => 'L' , UUT => 'F' , UUU => 'F' , UUR => 'L' ,
91             UUY => 'F' , URA => '*' , MGA => 'R' , MGG => 'R' ,
92             MGR => 'R' , YTA => 'L' , YTG => 'L' , YTR => 'L' ,
93             YUA => 'L' , YUG => 'L' , YUR => 'L' ,
94              
95             );
96              
97             sub new {
98              
99 133     133 1 465 my ($class, $seq, $id, $desc, $qual) = @_;
100              
101 133 100 100     493 if ( defined $seq && defined $qual
      100        
102             && (length($seq) != length($qual))) {
103 1         8 die "Sequence/quality length mismatch";
104             }
105              
106 132         220 my $self = bless {}, $class;
107 132   100     321 $self->{seq} = $seq // '';
108 132   100     270 $self->{id} = $id // undef;
109 132   100     292 $self->{desc} = $desc // undef;
110 132   100     301 $self->{qual} = $qual // undef;
111              
112 132         296 return $self;
113              
114             }
115              
116             sub seq : lvalue {
117              
118 29     29 1 701 my ($self,$new_val) = @_;
119 29 100       69 $self->{seq} = $new_val if (defined $new_val);
120 29         112 return $self->{seq};
121              
122             }
123              
124             sub id : lvalue {
125              
126 9     9 1 728 my ($self,$new_val) = @_;
127 9 100       22 $self->{id} = $new_val if (defined $new_val);
128 9         53 return $self->{id};
129              
130             }
131              
132             sub desc : lvalue {
133              
134 9     9 1 22 my ($self,$new_val) = @_;
135 9 100       18 $self->{desc} = $new_val if (defined $new_val);
136 9         44 return $self->{desc};
137              
138             }
139              
140             sub qual : lvalue {
141              
142 10     10 1 279 my ($self,$new_val) = @_;
143 10 100       20 $self->{qual} = $new_val if (defined $new_val);
144 10         44 return $self->{qual};
145              
146             }
147              
148             sub range {
149              
150 4     4 1 239 my ($self, $start, $end) = @_;
151 4 100 100     15 if ($start < 1 || $end > length($self->{seq})) {
152 2         18 warn "Range outside of sequence length\n";
153 2         9 return undef;
154             }
155 2         5 my $seq = substr $self->{seq}, $start-1, $end-$start+1;
156             my $qual = defined $self->{qual}
157 2 100       7 ? substr $self->{qual}, $start-1, $end-$start+1
158             : undef;
159             return __PACKAGE__->new(
160             $seq,
161             "$self->{id}_$start-$end",
162             $self->{desc},
163 2         9 $qual,
164             );
165              
166             }
167              
168             sub as_fasta {
169              
170 4     4 1 7 my ($self, $line_length) = @_;
171 4   100     22 my $l = $line_length // 60;
172 4 100       10 if (! defined $self->{id}) {
173 1         9 warn "Can't write FASTA with undefined ID\n";
174 1         6 return undef;
175             }
176 3         7 my $string = '>' . $self->{id};
177 3 100       7 $string .= ' ' . $self->{desc} if (defined $self->{desc});
178 3         6 $string .= "\n";
179 3         3 my $i = 0;
180 3         8 while ($i < length($self->{seq})) {
181 5         9 $string .= substr($self->{seq}, $i, $l) . "\n";
182 5         7 $i += $l;
183             }
184 3         11 return $string;
185              
186             }
187              
188             sub as_fastq {
189              
190 6     6 1 323 my ($self, $qual) = @_;
191 6 100       47 if (! defined $self->{id}) {
192 1         11 warn "Can't write FASTQ with undefined ID\n";
193 1         6 return undef;
194             }
195 5         16 my $string = '@' . $self->{id};
196 5 100       16 $string .= ' ' . $self->{desc} if (defined $self->{desc});
197 5         8 $string .= "\n";
198 5         7 $string .= $self->{seq};
199 5         6 $string .= "\n+\n";
200              
201             # check string lengths
202 5 100       15 if (defined $self->{qual}) {
    100          
203 3 100       7 if (length $self->{qual} != length $self->{seq}) {
204 1         7 die "Sequence/quality length mismatch";
205             }
206             }
207             # check that quality is defined somewhere
208             elsif (! defined $qual) {
209 1         7 die "Attempt to write FASTQ with undefined quality";
210             }
211              
212             # populate qual with constant quality if not defined
213             $string .= defined $qual
214             ? chr($qual+33) x length($self->{seq})
215 3 100       11 : $self->{qual};
216 3         4 $string .= "\n";
217 3         14 return $string;
218              
219             }
220              
221             sub as_input {
222              
223 2     2 1 5 my ($self, $arg2) = @_;
224             die "No input format found"
225 2 50       7 if (! defined $self->{_input_format});
226 2         9 my $method = 'as_' . $self->{_input_format};
227 2         18 $self->$method();
228              
229             }
230              
231             sub rev_com {
232              
233 7     7 1 258 my ($self) = @_;
234              
235 7         11 my $seq = $self->{seq};
236 7         12 $seq =~ tr/Xx/Nn/;
237 7 100       14 if (! _is_nucleic($seq) ) {
238 1         25 warn "Bad input sequence\n";
239 1         8 return undef;
240             }
241 6         14 $seq = reverse $seq;
242 6         9 $seq =~ tr
243             {ACGTMRWSYKVHDBNacgtmrwsykvhdbn-}
244             {TGCAKYWSRMBDHVNtgcakywsrmbdhvn-};
245              
246 6         7 my $qual = $self->{qual};
247 6 100       13 $qual = reverse $qual if (defined $qual);
248              
249             # If in void context, act in-place
250 6 100       11 if (! defined wantarray) {
251 1         2 $self->{seq} = $seq;
252 1         2 $self->{qual} = $qual;
253 1         2 return 1;
254             }
255              
256             # else return a new sequence object
257             return __PACKAGE__->new(
258             $seq,
259             $self->{id},
260             $self->{desc},
261 5         12 $qual,
262             );
263              
264             }
265              
266             sub translate {
267              
268 6     6 1 70 my ($self,$frame) = @_;
269              
270 6   100     18 $frame = $frame // 0;
271 6 100 100     31 die "$frame is not a valid frame (must be between 0 and 5)\n"
272             if ($frame < 0 || $frame > 5);
273              
274 4 100       10 my $seq = uc( $frame > 2 ? $self->rev_com->seq : $self->seq );
275 4         14 $seq = substr $seq, $frame%3;
276 4         7 $seq =~ tr/X/N/;
277 4 100       7 if (! _is_nucleic($seq) ) {
278 2         18 warn "Input doesn't look like DNA\n";
279 2         8 return undef;
280             }
281              
282 2   100     14 $seq = join('', map {$genetic_code{$_} // 'X'}
  5         29  
283             unpack 'A3' x int(length($seq)/3), $seq);
284              
285             # If in void context, act in-place
286 2 100       6 if (! defined wantarray) {
287 1         2 $self->{seq} = $seq;
288 1         2 $self->{qual} = undef;
289 1         3 return 1;
290             }
291              
292             # else return a new sequence object
293             return __PACKAGE__->new(
294             $seq,
295             $self->{id},
296             $self->{desc},
297 1         5 );
298              
299             }
300              
301             sub _stringify {
302              
303 125     125   878 my ($self) = @_;
304 125   100     5759 return $self->{seq} // '';
305              
306             }
307              
308             sub _concat {
309              
310 1     1   2 my ($self,$addition,$other) = @_;
311 1         3 $self->{seq} .= $addition;
312 1         2 return $self;
313              
314             }
315              
316             sub _is_nucleic {
317              
318 11     11   15 my ($seq) = @_;
319 11         140 return $seq !~ /[^ACGTUMRWSYKVHDBN-]/i;
320             }
321              
322              
323             1;
324              
325              
326             __END__