# This module implements classes for nucleotide chains.
#
# Written by Konrad Hinsen
#
"""
Nucleic acid chains
"""
__docformat__ = 'restructuredtext'
from MMTK import Biopolymers, Bonds, ChemicalObjects, Collections, \
ConfigIO, Database, Universe, Utility
from Scientific.Geometry import Vector
from MMTK.Biopolymers import defineNucleicAcidResidue
#
# Nucleotides are special groups
#
[docs]class Nucleotide(Biopolymers.Residue):
"""
Nucleic acid residue
Nucleotides are a special kind of group. Like any other
group, they are defined in the chemical database. Each residue
has two or three subgroups ('sugar' and 'base', plus 'phosphate'
except for 5'-terminal residues) and is usually
connected to other residues to form a nucleotide chain. The database
contains three variants of each residue (5'-terminal, 3'-terminal,
non-terminal).
"""
def __init__(self, name = None, model = 'all'):
"""
:param name: the name of the nucleotide in the chemical database. This
is the full name of the residue plus the suffix
"_5ter" or "_3ter" for the terminal variants.
:type name: str
:param model: one of "all" (all-atom), "none" (no hydrogens),
"polar" (united-atom with only polar hydrogens),
"polar_charmm" (like "polar", but defining
polar hydrogens like in the CHARMM force field).
Currently the database has definitions only for "all".
:type model: str
"""
if name is not None:
blueprint = _residueBlueprint(name, model)
ChemicalObjects.Group.__init__(self, blueprint)
self.model = model
self._init()
[docs] def backbone(self):
"""
:returns: the sugar and phosphate groups
:rtype: :class:`~MMTK.ChemicalObjects.Group`
"""
bb = self.sugar
if hasattr(self, 'phosphate'):
bb = Collections.Collection([bb, self.phosphate])
return bb
[docs] def bases(self):
"""
:returns: the base group
:rtype: :class:`~MMTK.ChemicalObjects.Group`
"""
return self.base
def _residueBlueprint(name, model):
try:
blueprint = _residue_blueprints[(name, model)]
except KeyError:
if model == 'polar':
name = name + '_uni'
elif model == 'polar_charmm':
name = name + '_uni2'
elif model == 'none':
name = name + '_noh'
blueprint = Database.BlueprintGroup(name)
_residue_blueprints[(name, model)] = blueprint
return blueprint
_residue_blueprints = {}
#
# Nucleotide chains are molecules with added features.
#
[docs]class NucleotideChain(Biopolymers.ResidueChain):
"""
Nucleotide chain
Nucleotide chains consist of nucleotides that are linked together.
They are a special kind of molecule, i.e. all molecule operations
are available.
Nucleotide chains act as sequences of residues. If n is a NucleotideChain
object, then
* len(n) yields the number of nucleotides
* n[i] yields nucleotide number i
* n[i:j] yields the subchain from nucleotide number i up to but
excluding nucleotide number j
"""
def __init__(self, sequence, **properties):
"""
:param sequence: the nucleotide sequence. This can be a string
containing the one-letter codes, or a list
of two-letter codes (a "d" or "r" for the type of
sugar, and the one-letter base code), or a
:class:`~MMTK.PDB.PDBNucleotideChain` object.
If a PDBNucleotideChain object is supplied, the
atomic positions it contains are assigned to the
atoms of the newly generated nucleotide chain,
otherwise the positions of all atoms are undefined.
:keyword model: one of "all" (all-atom), "no_hydrogens" or "none"
(no hydrogens), "polar_hydrogens" or "polar"
(united-atom with only polar hydrogens),
"polar_charmm" (like "polar", but defining
polar hydrogens like in the CHARMM force field),
"polar_opls" (like "polar", but defining
polar hydrogens like in the latest OPLS force field).
Default is "all". Currently the database contains
definitions only for "all".
:type model: str
:keyword terminus_5: if True, the first residue is constructed
using the 5'-terminal variant, if False the
non-terminal version is used. Default is True.
:type terminus_5: bool
:keyword terminus_3: if True, the last residue is constructed
using the 3'-terminal variant, if False the
non-terminal version is used. Default is True.
:type terminus_3: bool
:keyword circular: if True, a bond is constructed
between the first and the last residue.
Default is False.
:type circular: bool
:keyword name: a name for the chain (a string)
:type name: str
"""
if sequence is not None:
hydrogens = self.binaryProperty(properties, 'hydrogens', 'all')
if hydrogens == 1:
hydrogens = 'all'
elif hydrogens == 0:
hydrogens = 'none'
term5 = self.binaryProperty(properties, 'terminus_5', True)
term3 = self.binaryProperty(properties, 'terminus_3', True)
circular = self.binaryProperty(properties, 'circular', False)
try:
model = properties['model'].lower()
except KeyError:
model = hydrogens
self.version_spec = {'hydrogens': hydrogens,
'terminus_5': term5,
'terminus_3': term3,
'model': model,
'circular': circular}
if isinstance(sequence[0], basestring):
conf = None
else:
conf = sequence
sequence = [r.name for r in sequence]
sequence = [Biopolymers._fullName(r) for r in sequence]
if term5:
if sequence[0].find('5ter') == -1:
sequence[0] += '_5ter'
if term3:
if sequence[-1].find('3ter') == -1:
sequence[-1] += '_3ter'
self.groups = []
n = 0
for residue in sequence:
n += 1
r = Nucleotide(residue, model)
r.name = r.symbol + '_' + `n`
r.sequence_number = n
r.parent = self
self.groups.append(r)
self._setupChain(circular, properties, conf)
is_nucleotide_chain = True
def __getslice__(self, first, last):
return NucleotideSubChain(self, self.groups[first:last])
[docs] def backbone(self):
"""
:returns: the sugar and phosphate groups of all nucleotides
:rtype: :class:`~MMTK.Collections.Collection`
"""
bb = Collections.Collection([])
for residue in self.groups:
try:
bb.addObject(residue.phosphate)
except AttributeError:
pass
bb.addObject(residue.sugar)
return bb
[docs] def bases(self):
"""
:returns: the base groups of all nucleotides
:rtype: :class:`~MMTK.Collections.Collection`
"""
return Collections.Collection([r.base for r in self.groups])
def _descriptionSpec(self):
kwargs = ','.join([name + '=' + `self.version_spec[name]`
for name in sorted(self.version_spec.keys())])
return "N", kwargs
def _typeName(self):
return ''.join([s.ljust(3) for s in self.sequence()])
def _graphics(self, conf, distance_fn, model, module, options):
if model != 'backbone':
return ChemicalObjects.Molecule._graphics(self, conf,
distance_fn, model,
module, options)
color = options.get('color', 'black')
material = module.EmissiveMaterial(color)
objects = []
for i in range(1, len(self.groups)-1):
a1 = self.groups[i].phosphate.P
a2 = self.groups[i+1].phosphate.P
p1 = a1.position(conf)
p2 = a2.position(conf)
if p1 is not None and p2 is not None:
bond_vector = 0.5*distance_fn(a1, a2, conf)
cut = bond_vector != 0.5*(p2-p1)
if not cut:
objects.append(module.Line(p1, p2, material = material))
else:
objects.append(module.Line(p1, p1+bond_vector,
material = material))
objects.append(module.Line(p2, p2-bond_vector,
material = material))
return objects
#
# Subchains are created by slicing chains or extracting a chain from
# a group of connected chains.
#
[docs]class NucleotideSubChain(NucleotideChain):
"""
A contiguous part of a nucleotide chain
NucleotideSubChain objects are the result of slicing operations on
NucleotideChain objects. They cannot be created directly.
NucleotideSubChain objects permit all operations of NucleotideChain
objects, but cannot be added to a universe.
"""
def __init__(self, chain, groups, name = ''):
self.groups = groups
self.atoms = []
self.bonds = []
for g in self.groups:
self.atoms = self.atoms + g.atoms
self.bonds = self.bonds + g.bonds
for i in range(len(self.groups)-1):
self.bonds.append(Bonds.Bond((self.groups[i].sugar.O_3,
self.groups[i+1].phosphate.P)))
self.bonds = Bonds.BondList(self.bonds)
self.name = name
self.parent = chain.parent
self.type = None
self.configurations = {}
self.part_of = chain
is_incomplete = True
def __repr__(self):
if self.name == '':
return 'SubChain of ' + repr(self.part_of)
else:
return ChemicalObjects.Molecule.__repr__(self)
__str__ = __repr__
#
# Type check functions
#
[docs]def isNucleotideChain(x):
"Returns True if x is a NucleotideChain."
return hasattr(x, 'is_nucleotide_chain')