/
_dna.py
423 lines (364 loc) · 12.4 KB
/
_dna.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
# ----------------------------------------------------------------------------
# Copyright (c) 2013--, scikit-bio development team.
#
# Distributed under the terms of the Modified BSD License.
#
# The full license is in the file LICENSE.txt, distributed with this software.
# ----------------------------------------------------------------------------
import skbio
from skbio.util._decorator import classproperty, overrides
from skbio.util._decorator import stable
from ._nucleotide_mixin import NucleotideMixin, _motifs as _parent_motifs
from ._grammared_sequence import GrammaredSequence, DisableSubclassingMeta
class DNA(GrammaredSequence, NucleotideMixin,
metaclass=DisableSubclassingMeta):
r"""Store DNA sequence data and optional associated metadata.
Only characters in the IUPAC DNA character set [1]_ are supported.
Parameters
----------
sequence : str, Sequence, or 1D np.ndarray (np.uint8 or '\|S1')
Characters representing the DNA sequence itself.
metadata : dict, optional
Arbitrary metadata which applies to the entire sequence.
positional_metadata : Pandas DataFrame consumable, optional
Arbitrary per-character metadata. For example, quality data from
sequencing reads. Must be able to be passed directly to the Pandas
DataFrame constructor.
interval_metadata : IntervalMetadata
Arbitrary interval metadata which applies to intervals within
a sequence to store interval features (such as genes on the
DNA sequence).
lowercase : bool or str, optional
If ``True``, lowercase sequence characters will be converted to
uppercase characters in order to be valid IUPAC DNA characters. If
``False``, no characters will be converted. If a str, it will be
treated as a key into the positional metadata of the object. All
lowercase characters will be converted to uppercase, and a ``True``
value will be stored in a boolean array in the positional metadata
under the key.
validate : bool, optional
If ``True``, validation will be performed to ensure that all sequence
characters are in the IUPAC DNA character set. If ``False``, validation
will not be performed. Turning off validation will improve runtime
performance. If invalid characters are present, however, there is
**no guarantee that operations performed on the resulting object will
work or behave as expected.** Only turn off validation if you are
certain that the sequence characters are valid. To store sequence data
that is not IUPAC-compliant, use ``Sequence``.
See Also
--------
RNA
GrammaredSequence
Notes
-----
Subclassing is disabled for DNA, because subclassing makes
it possible to change the alphabet, and certain methods rely on the
IUPAC alphabet. If a custom sequence alphabet is needed, inherit directly
from ``GrammaredSequence``.
References
----------
.. [1] Nomenclature for incompletely specified bases in nucleic acid
sequences: recommendations 1984.
Nucleic Acids Res. May 10, 1985; 13(9): 3021-3030.
A Cornish-Bowden
Examples
--------
>>> from skbio import DNA
>>> DNA('ACCGAAT')
DNA
--------------------------
Stats:
length: 7
has gaps: False
has degenerates: False
has definites: True
GC-content: 42.86%
--------------------------
0 ACCGAAT
Convert lowercase characters to uppercase:
>>> DNA('AcCGaaT', lowercase=True)
DNA
--------------------------
Stats:
length: 7
has gaps: False
has degenerates: False
has definites: True
GC-content: 42.86%
--------------------------
0 ACCGAAT
"""
@classproperty
@overrides(NucleotideMixin)
def complement_map(cls):
comp_map = {
'A': 'T', 'T': 'A', 'G': 'C', 'C': 'G', 'Y': 'R', 'R': 'Y',
'S': 'S', 'W': 'W', 'K': 'M', 'M': 'K', 'B': 'V', 'D': 'H',
'H': 'D', 'V': 'B', 'N': 'N'
}
comp_map.update({c: c for c in cls.gap_chars})
return comp_map
@classproperty
@overrides(GrammaredSequence)
def definite_chars(cls):
return set("ACGT")
@classproperty
@overrides(GrammaredSequence)
def degenerate_map(cls):
return {
"R": set("AG"), "Y": set("CT"), "M": set("AC"), "K": set("TG"),
"W": set("AT"), "S": set("GC"), "B": set("CGT"), "D": set("AGT"),
"H": set("ACT"), "V": set("ACG"), "N": set("ACGT")
}
@classproperty
@overrides(GrammaredSequence)
def default_gap_char(cls):
return '-'
@classproperty
@overrides(GrammaredSequence)
def gap_chars(cls):
return set('-.')
@property
def _motifs(self):
return _motifs
@stable(as_of="0.4.0")
def transcribe(self):
"""Transcribe DNA into RNA.
DNA sequence is assumed to be the coding strand. Thymine (T) is
replaced with uracil (U) in the transcribed sequence.
Returns
-------
RNA
Transcribed sequence.
See Also
--------
translate
translate_six_frames
Notes
-----
DNA sequence's metadata, positional, and interval
metadata are included in the transcribed RNA sequence.
Examples
--------
Transcribe DNA into RNA:
>>> from skbio import DNA
>>> dna = DNA('TAACGTTA')
>>> dna
DNA
--------------------------
Stats:
length: 8
has gaps: False
has degenerates: False
has definites: True
GC-content: 25.00%
--------------------------
0 TAACGTTA
>>> dna.transcribe()
RNA
--------------------------
Stats:
length: 8
has gaps: False
has degenerates: False
has definites: True
GC-content: 25.00%
--------------------------
0 UAACGUUA
"""
seq = self._string.replace(b'T', b'U')
metadata = None
if self.has_metadata():
metadata = self.metadata
positional_metadata = None
if self.has_positional_metadata():
positional_metadata = self.positional_metadata
interval_metadata = None
if self.has_interval_metadata():
interval_metadata = self.interval_metadata
# turn off validation because `seq` is guaranteed to be valid
return skbio.RNA(seq, metadata=metadata,
positional_metadata=positional_metadata,
interval_metadata=interval_metadata,
validate=False)
@stable(as_of="0.4.0")
def translate(self, *args, **kwargs):
"""Translate DNA sequence into protein sequence.
DNA sequence is assumed to be the coding strand. DNA sequence is first
transcribed into RNA and then translated into protein.
Parameters
----------
args : tuple
Positional arguments accepted by ``RNA.translate``.
kwargs : dict
Keyword arguments accepted by ``RNA.translate``.
Returns
-------
Protein
Translated sequence.
See Also
--------
RNA.reverse_transcribe
RNA.translate
translate_six_frames
transcribe
Notes
-----
DNA sequence's metadata are included in the translated protein
sequence. Positional metadata are not included.
Examples
--------
Translate DNA into protein using NCBI's standard genetic code (table ID
1, the default genetic code in scikit-bio):
>>> from skbio import DNA
>>> dna = DNA('ATGCCACTTTAA')
>>> dna.translate()
Protein
--------------------------
Stats:
length: 4
has gaps: False
has degenerates: False
has definites: True
has stops: True
--------------------------
0 MPL*
Translate the same DNA sequence using a different NCBI genetic code
(table ID 3, the yeast mitochondrial code) and specify that translation
must terminate at the first stop codon:
>>> dna.translate(3, stop='require')
Protein
--------------------------
Stats:
length: 3
has gaps: False
has degenerates: False
has definites: True
has stops: False
--------------------------
0 MPT
"""
return self.transcribe().translate(*args, **kwargs)
@stable(as_of="0.4.0")
def translate_six_frames(self, *args, **kwargs):
"""Translate DNA into protein using six possible reading frames.
DNA sequence is assumed to be the coding strand. DNA sequence is first
transcribed into RNA and then translated into protein. The six possible
reading frames are:
* 1 (forward)
* 2 (forward)
* 3 (forward)
* -1 (reverse)
* -2 (reverse)
* -3 (reverse)
Translated sequences are yielded in this order.
Parameters
----------
args : tuple
Positional arguments accepted by ``RNA.translate_six_frames``.
kwargs : dict
Keyword arguments accepted by ``RNA.translate_six_frames``.
Yields
------
Protein
Translated sequence in the current reading frame.
See Also
--------
RNA.translate_six_frames
translate
transcribe
Notes
-----
This method is faster than (and equivalent to) performing six
independent translations using, for example:
``(seq.translate(reading_frame=rf)
for rf in GeneticCode.reading_frames)``
DNA sequence's metadata are included in each translated protein
sequence. Positional metadata are not included.
Examples
--------
Translate DNA into protein using the six possible reading frames and
NCBI's standard genetic code (table ID 1, the default genetic code in
scikit-bio):
>>> from skbio import DNA
>>> dna = DNA('ATGCCACTTTAA')
>>> for protein in dna.translate_six_frames():
... protein
... print('')
Protein
--------------------------
Stats:
length: 4
has gaps: False
has degenerates: False
has definites: True
has stops: True
--------------------------
0 MPL*
<BLANKLINE>
Protein
--------------------------
Stats:
length: 3
has gaps: False
has degenerates: False
has definites: True
has stops: False
--------------------------
0 CHF
<BLANKLINE>
Protein
--------------------------
Stats:
length: 3
has gaps: False
has degenerates: False
has definites: True
has stops: False
--------------------------
0 ATL
<BLANKLINE>
Protein
--------------------------
Stats:
length: 4
has gaps: False
has degenerates: False
has definites: True
has stops: False
--------------------------
0 LKWH
<BLANKLINE>
Protein
--------------------------
Stats:
length: 3
has gaps: False
has degenerates: False
has definites: True
has stops: True
--------------------------
0 *SG
<BLANKLINE>
Protein
--------------------------
Stats:
length: 3
has gaps: False
has degenerates: False
has definites: True
has stops: False
--------------------------
0 KVA
<BLANKLINE>
"""
return self.transcribe().translate_six_frames(*args, **kwargs)
@overrides(GrammaredSequence)
def _repr_stats(self):
"""Define custom statistics to display in the sequence's repr."""
stats = super(DNA, self)._repr_stats()
stats.append(('GC-content', '{:.2%}'.format(self.gc_content())))
return stats
_motifs = _parent_motifs.copy()
# Leave this at the bottom
_motifs.interpolate(DNA, "find_motifs")