Browse Source

Add write support for DNA Strider format.

snapgene-write
Damien Goutte-Gattat 3 years ago
parent
commit
6b832d212a
  1. 90
      incenp/bio/seqio/XdnaIO.py
  2. 3
      incenp/bio/seqio/__init__.py

90
incenp/bio/seqio/XdnaIO.py

@ -22,9 +22,10 @@
from Bio import Alphabet
from Bio.Seq import Seq
from Bio.SeqIO.Interfaces import SequenceWriter
from Bio.SeqRecord import SeqRecord
from Bio.SeqFeature import SeqFeature, FeatureLocation
from struct import unpack
from Bio.SeqFeature import SeqFeature, FeatureLocation, ExactPosition
from struct import unpack, pack
from re import match
_seq_types = {
@ -154,3 +155,88 @@ def XdnaIterator(handle):
num_features -= 1
yield record
class XdnaWriter(SequenceWriter):
def write_file(self, records):
if len(records) != 1:
raise ValueError("A DNA Strider file can only contain one record")
record = records[0]
alptype = type(record.seq.alphabet)
if issubclass(alptype, Alphabet.DNAAlphabet):
seqtype = 1
elif issubclass(alptype, Alphabet.RNAAlphabet):
seqtype = 3
elif issubclass(alptype, Alphabet.ProteinAlphabet):
seqtype = 4
else:
seqtype = 0
if record.annotations.get('topology', 'linear') == 'circular':
topology = 1
else:
topology = 0
# Write header
self.handle.write(pack('>BBB25xII60xI11xB',
0, # version
seqtype, topology, len(record),
0, # negative length
len(record.description),
255 # end of header
))
# Actual sequence and comment
self.handle.write(str(record.seq))
self.handle.write(record.description)
self.handle.write(pack('>B', 0)) # Annotation section marker
self._write_pstring('0') # right-side overhang
self._write_pstring('0') # left-side overhand
# Write features
self.handle.write(pack('>B', len(record.features)))
for feature in record.features:
if type(feature.location.start) != ExactPosition or type(feature.location.end) != ExactPosition:
# Cannot store fuzzy locations, skip feature
continue
self._write_pstring(feature.qualifiers.get('label', [''])[0])
description = ''
for qname in feature.qualifiers:
if qname in ('label', 'translation'):
continue
for val in feature.qualifiers[qname]:
if len(description) > 0:
description = description + '\x0D'
description = description + '%s="%s"' % (qname, val)
self._write_pstring(description)
self._write_pstring(feature.type)
start = feature.location.start.position + 1 # 1-based coordinates
end = feature.location.end.position
strand = 1
if feature.location.strand == -1:
start, end = end, start
strand = 0
self._write_pstring(str(start))
self._write_pstring(str(end))
self.handle.write(pack('>BBBB', strand, 1, 0, 1))
self._write_pstring('127,127,127')
return 1
def _write_pstring(self, s):
self.handle.write(pack('>B', len(s)))
self.handle.write(s)

3
incenp/bio/seqio/__init__.py

@ -27,6 +27,9 @@ from . import SnapGeneIO
if not SeqIO._FormatToIterator.has_key('xdna'):
SeqIO._FormatToIterator['xdna'] = XdnaIO.XdnaIterator
if not SeqIO._FormatToWriter.has_key('xdna'):
SeqIO._FormatToWriter['xdna'] = XdnaIO.XdnaWriter
if not 'xdna' in SeqIO._BinaryFormats:
SeqIO._BinaryFormats.append('xdna')

Loading…
Cancel
Save