Skip to content

P
pdbfixer
Commit 4ea8b0bf by peastman
Options

Merge branch 'peastman-ids'

parents ffa808a3 e41a190b
Showing with 31 additions and 40 deletions
pdbfixer/pdbfixer.py
...@@ -6,7 +6,7 @@ Simbios, the NIH National Center for Physics-Based Simulation of ...@@ -6,7 +6,7 @@ Simbios, the NIH National Center for Physics-Based Simulation of
Biological Structures at Stanford, funded under the NIH Roadmap for Biological Structures at Stanford, funded under the NIH Roadmap for
Medical Research, grant U54 GM072970. See https://simtk.org. Medical Research, grant U54 GM072970. See https://simtk.org.
Portions copyright (c) 2013-2014 Stanford University and the Authors. Portions copyright (c) 2013-2015 Stanford University and the Authors.
Authors: Peter Eastman Authors: Peter Eastman
Contributors: Contributors:
...@@ -221,7 +221,6 @@ class PDBFixer(object): ...@@ -221,7 +221,6 @@ class PDBFixer(object):
self.pdb = app.PDBFile(structure) self.pdb = app.PDBFile(structure)
self.topology = self.pdb.topology self.topology = self.pdb.topology
self.positions = self.pdb.positions self.positions = self.pdb.positions
self.structureChains = list(self.structure.iter_chains())
# Load the templates. # Load the templates.
...@@ -462,7 +461,6 @@ class PDBFixer(object): ...@@ -462,7 +461,6 @@ class PDBFixer(object):
modeller.delete(allChains[i] for i in chainIndices) modeller.delete(allChains[i] for i in chainIndices)
self.topology = modeller.topology self.topology = modeller.topology
self.positions = modeller.positions self.positions = modeller.positions
self.structureChains = [self.structureChains[i] for i in range(len(self.structureChains)) if i not in chainIndices]
return return
...@@ -481,17 +479,17 @@ class PDBFixer(object): ...@@ -481,17 +479,17 @@ class PDBFixer(object):
>>> missing_residues = fixer.missingResidues >>> missing_residues = fixer.missingResidues
""" """
chains = [c for c in self.structureChains if any(atom.record_name == 'ATOM' for atom in c.iter_atoms())] chains = [c for c in self.topology.chains() if len(list(c.residues())) > 0]
chainWithGaps = {} chainWithGaps = {}
# Find the sequence of each chain, with gaps for missing residues. # Find the sequence of each chain, with gaps for missing residues.
for chain in chains: for chain in chains:
minResidue = min(r.number for r in chain.iter_residues()) minResidue = min(int(r.id) for r in chain.residues())
maxResidue = max(r.number for r in chain.iter_residues()) maxResidue = max(int(r.id) for r in chain.residues())
residues = [None]*(maxResidue-minResidue+1) residues = [None]*(maxResidue-minResidue+1)
for r in chain.iter_residues(): for r in chain.residues():
residues[r.number-minResidue] = r.get_name() residues[int(r.id)-minResidue] = r.name
chainWithGaps[chain] = residues chainWithGaps[chain] = residues
# Try to find the chain that matches each sequence. # Try to find the chain that matches each sequence.
...@@ -500,7 +498,7 @@ class PDBFixer(object): ...@@ -500,7 +498,7 @@ class PDBFixer(object):
chainOffset = {} chainOffset = {}
for sequence in self.structure.sequences: for sequence in self.structure.sequences:
for chain in chains: for chain in chains:
if chain.chain_id != sequence.chain_id: if chain.id != sequence.chain_id:
continue continue
if chain in chainSequence: if chain in chainSequence:
continue continue
...@@ -515,15 +513,15 @@ class PDBFixer(object): ...@@ -515,15 +513,15 @@ class PDBFixer(object):
# Now build the list of residues to add. # Now build the list of residues to add.
self.missingResidues = {} self.missingResidues = {}
for structChain, topChain in zip(self.structureChains, self.topology.chains()): for chain in self.topology.chains():
if structChain in chainSequence: if chain in chainSequence:
offset = chainOffset[structChain] offset = chainOffset[chain]
sequence = chainSequence[structChain].residues sequence = chainSequence[chain].residues
gappedSequence = chainWithGaps[structChain] gappedSequence = chainWithGaps[chain]
index = 0 index = 0
for i in range(len(sequence)): for i in range(len(sequence)):
if i < offset or i >= len(gappedSequence)+offset or gappedSequence[i-offset] is None: if i < offset or i >= len(gappedSequence)+offset or gappedSequence[i-offset] is None:
key = (topChain.index, index) key = (chain.index, index)
if key not in self.missingResidues: if key not in self.missingResidues:
self.missingResidues[key] = [] self.missingResidues[key] = []
residueName = sequence[i] residueName = sequence[i]
...@@ -556,16 +554,16 @@ class PDBFixer(object): ...@@ -556,16 +554,16 @@ class PDBFixer(object):
# Now add ones based on MODRES records. # Now add ones based on MODRES records.
modres = dict(((m.chain_id, m.number, m.residue_name), m.standard_name) for m in self.structure.modified_residues) modres = dict(((m.chain_id, str(m.number), m.residue_name), m.standard_name) for m in self.structure.modified_residues)
for structChain, topChain in zip(self.structureChains, self.topology.chains()): for chain in self.topology.chains():
for structResidue, topResidue in zip(structChain.iter_residues(), topChain.residues()): for residue in chain.residues():
key = (structChain.chain_id, structResidue.number, structResidue.name) key = (chain.id, residue.id, residue.name)
if key in modres: if key in modres:
replacement = modres[key] replacement = modres[key]
if replacement == 'DU': if replacement == 'DU':
replacement = 'DT' replacement = 'DT'
if replacement in self.templates: if replacement in self.templates:
nonstandard[topResidue] = replacement nonstandard[residue] = replacement
self.nonstandardResidues = [(r, nonstandard[r]) for r in sorted(nonstandard, key=lambda r: r.index)] self.nonstandardResidues = [(r, nonstandard[r]) for r in sorted(nonstandard, key=lambda r: r.index)]
def replaceNonstandardResidues(self): def replaceNonstandardResidues(self):
......
pdbfixer/ui.py
import simtk.openmm.app as app import simtk.openmm.app as app
from simtk.openmm.app.internal.pdbstructure import PdbStructure
import simtk.unit as unit import simtk.unit as unit
from pdbfixer import PDBFixer, substitutions, proteinResidues, dnaResidues, rnaResidues from pdbfixer import PDBFixer, substitutions, proteinResidues, dnaResidues, rnaResidues
import uiserver import uiserver
...@@ -74,7 +73,7 @@ def addResiduesPageCallback(parameters, handler): ...@@ -74,7 +73,7 @@ def addResiduesPageCallback(parameters, handler):
for i, key in enumerate(keys): for i, key in enumerate(keys):
if 'add'+str(i) not in parameters: if 'add'+str(i) not in parameters:
del fixer.missingResidues[key] del fixer.missingResidues[key]
displayMissingAtomsPage() displayConvertResiduesPage()
def convertResiduesPageCallback(parameters, handler): def convertResiduesPageCallback(parameters, handler):
for i in range(len(fixer.nonstandardResidues)): for i in range(len(fixer.nonstandardResidues)):
...@@ -135,7 +134,7 @@ def displayDeleteChainsPage(): ...@@ -135,7 +134,7 @@ def displayDeleteChainsPage():
content = "DNA" content = "DNA"
else: else:
content = ', '.join(set(residues)) content = ', '.join(set(residues))
table += ' <tr><td>%d</td><td>%d</td><td>%s</td><td><input type="checkbox" name="include%d" checked></td></tr>\n' % (chain.index+1, len(residues), content, i) table += ' <tr><td>%s</td><td>%d</td><td>%s</td><td><input type="checkbox" name="include%d" checked></td></tr>\n' % (chain.id, len(residues), content, i)
uiserver.setContent(header+loadHtmlFile("removeChains.html") % (numChains, table)) uiserver.setContent(header+loadHtmlFile("removeChains.html") % (numChains, table))
def displayAddResiduesPage(): def displayAddResiduesPage():
...@@ -145,14 +144,16 @@ def displayAddResiduesPage(): ...@@ -145,14 +144,16 @@ def displayAddResiduesPage():
displayConvertResiduesPage() displayConvertResiduesPage()
return return
table = "" table = ""
chains = list(fixer.topology.chains())
for i, key in enumerate(sorted(fixer.missingResidues)): for i, key in enumerate(sorted(fixer.missingResidues)):
residues = fixer.missingResidues[key] residues = fixer.missingResidues[key]
chain = fixer.structureChains[key[0]] chain = chains[key[0]]
if key[1] < len(chain.residues): chainResidues = list(chain.residues())
offset = chain.residues[key[1]].number-len(residues)-1 if key[1] < len(chainResidues):
offset = int(chainResidues[key[1]].id)-len(residues)-1
else: else:
offset = chain.residues[-1].number offset = int(chainResidues[-1].id)
table += ' <tr><td>%d</td><td>%d to %d</td><td>%s</td><td><input type="checkbox" name="add%d" checked></td></tr>\n' % (key[0]+1, offset+1, offset+len(residues), ', '.join(residues), i) table += ' <tr><td>%s</td><td>%d to %d</td><td>%s</td><td><input type="checkbox" name="add%d" checked></td></tr>\n' % (chain.id, offset+1, offset+len(residues), ', '.join(residues), i)
uiserver.setContent(header+loadHtmlFile("addResidues.html") % table) uiserver.setContent(header+loadHtmlFile("addResidues.html") % table)
def displayConvertResiduesPage(): def displayConvertResiduesPage():
...@@ -161,10 +162,6 @@ def displayConvertResiduesPage(): ...@@ -161,10 +162,6 @@ def displayConvertResiduesPage():
if len(fixer.nonstandardResidues) == 0: if len(fixer.nonstandardResidues) == 0:
displayMissingAtomsPage() displayMissingAtomsPage()
return return
indexInChain = {}
for structChain, topChain in zip(fixer.structureChains, fixer.topology.chains()):
for structResidue, topResidue in zip(structChain.iter_residues(), topChain.residues()):
indexInChain[topResidue] = structResidue.number
table = '' table = ''
nucleotides = ['DA', 'DC', 'DG', 'DT', 'A', 'C', 'G', 'T'] nucleotides = ['DA', 'DC', 'DG', 'DT', 'A', 'C', 'G', 'T']
for i in range(len(fixer.nonstandardResidues)): for i in range(len(fixer.nonstandardResidues)):
...@@ -179,7 +176,7 @@ def displayConvertResiduesPage(): ...@@ -179,7 +176,7 @@ def displayConvertResiduesPage():
if res == replaceWith: if res == replaceWith:
selected = ' selected' selected = ' selected'
options += '<option value="%s"%s>%s</option>' % (res, selected, res) options += '<option value="%s"%s>%s</option>' % (res, selected, res)
table += ' <tr><td>%d</td><td>%s %d</td><td><select name="residue%d">%s</select></td><td><input type="checkbox" name="convert%d" checked></td></tr>\n' % (residue.chain.index+1, residue.name, indexInChain[residue], i, options, i) table += ' <tr><td>%s</td><td>%s %s</td><td><select name="residue%d">%s</select></td><td><input type="checkbox" name="convert%d" checked></td></tr>\n' % (residue.chain.id, residue.name, residue.id, i, options, i)
uiserver.setContent(header+loadHtmlFile("convertResidues.html") % table) uiserver.setContent(header+loadHtmlFile("convertResidues.html") % table)
def displayMissingAtomsPage(): def displayMissingAtomsPage():
...@@ -191,10 +188,6 @@ def displayMissingAtomsPage(): ...@@ -191,10 +188,6 @@ def displayMissingAtomsPage():
fixer.addMissingAtoms() fixer.addMissingAtoms()
displayAddHydrogensPage() displayAddHydrogensPage()
return return
indexInChain = {}
for structChain, topChain in zip(fixer.structureChains, fixer.topology.chains()):
for structResidue, topResidue in zip(structChain.iter_residues(), topChain.residues()):
indexInChain[topResidue] = structResidue.number
table = "" table = ""
for residue in allResidues: for residue in allResidues:
atoms = [] atoms = []
...@@ -202,7 +195,7 @@ def displayMissingAtomsPage(): ...@@ -202,7 +195,7 @@ def displayMissingAtomsPage():
atoms.extend(atom.name for atom in fixer.missingAtoms[residue]) atoms.extend(atom.name for atom in fixer.missingAtoms[residue])
if residue in fixer.missingTerminals: if residue in fixer.missingTerminals:
atoms.extend(atom for atom in fixer.missingTerminals[residue]) atoms.extend(atom for atom in fixer.missingTerminals[residue])
table += ' <tr><td>%d</td><td>%s %d</td><td>%s</td></tr>\n' % (residue.chain.index+1, residue.name, indexInChain[residue], ', '.join(atoms)) table += ' <tr><td>%s</td><td>%s %s</td><td>%s</td></tr>\n' % (residue.chain.id, residue.name, residue.id, ', '.join(atoms))
uiserver.setContent(header+loadHtmlFile("addHeavyAtoms.html") % table) uiserver.setContent(header+loadHtmlFile("addHeavyAtoms.html") % table)
def displayAddHydrogensPage(): def displayAddHydrogensPage():
......
tests/test_removechains.py
...@@ -9,7 +9,7 @@ def remove_chain_ids_and_verify(pdbid, chain_ids_to_remove, expected_chain_ids_r ...@@ -9,7 +9,7 @@ def remove_chain_ids_and_verify(pdbid, chain_ids_to_remove, expected_chain_ids_r
# Remove specified chains. # Remove specified chains.
fixer.removeChains(chainIds=chain_ids_to_remove) fixer.removeChains(chainIds=chain_ids_to_remove)
# Check to make sure asserted chains remain. # Check to make sure asserted chains remain.
chain_ids_remaining = [c.chain_id for c in fixer.structureChains] chain_ids_remaining = [c.id for c in fixer.topology.chains()]
assert_items_equal(chain_ids_remaining, expected_chain_ids_remaining) assert_items_equal(chain_ids_remaining, expected_chain_ids_remaining)
def test_removechain_ids(): def test_removechain_ids():
...@@ -25,7 +25,7 @@ def remove_chain_indices_and_verify(pdbid, chain_indices_to_remove, expected_cha ...@@ -25,7 +25,7 @@ def remove_chain_indices_and_verify(pdbid, chain_indices_to_remove, expected_cha
# Remove specified chains. # Remove specified chains.
fixer.removeChains(chainIndices=chain_indices_to_remove) fixer.removeChains(chainIndices=chain_indices_to_remove)
# Check to make sure asserted chains remain. # Check to make sure asserted chains remain.
chain_ids_remaining = [c.chain_id for c in fixer.structureChains] chain_ids_remaining = [c.id for c in fixer.topology.chains()]
assert_items_equal(chain_ids_remaining, expected_chain_ids_remaining) assert_items_equal(chain_ids_remaining, expected_chain_ids_remaining)
def test_removechain_indices(): def test_removechain_indices():
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment