From b80eff3f8fdb13614bc4f1311afde6589ca09bbc Mon Sep 17 00:00:00 2001
From: Gabriel Studer <gabriel.studer@unibas.ch>
Date: Tue, 4 Oct 2016 16:25:45 +0200
Subject: [PATCH] move all operations, that should lead to a gapless model into
one function
The idea is that this whole procedure can then be applied to only
a part of the model by providing chain_idx and resnum ranges.
---
doc/tests/scripts/modelling_steps.py | 28 ++-----
modelling/doc/pipeline.rst | 13 ++-
modelling/pymod/_closegaps.py | 32 +++++++-
modelling/pymod/_pipeline.py | 118 +++++++++++++++++++++++----
4 files changed, 143 insertions(+), 48 deletions(-)
diff --git a/doc/tests/scripts/modelling_steps.py b/doc/tests/scripts/modelling_steps.py
index 93b4a65d..58fe1a87 100644
--- a/doc/tests/scripts/modelling_steps.py
+++ b/doc/tests/scripts/modelling_steps.py
@@ -13,35 +13,23 @@ aln = io.LoadAlignment('data/1crn.fasta')
aln.AttachView(1, tpl.CreateFullView())
mhandle = modelling.BuildRawModel(aln)
-# perform loop modelling to close all gaps
-modelling.SetupDefaultBackboneScorer(mhandle)
+# we're not modelling termini
modelling.RemoveTerminalGaps(mhandle)
-modelling.CloseSmallDeletions(mhandle)
-for distance in range(merge_distance):
- modelling.MergeGapsByDistance(mhandle, distance)
- modelling.FillLoopsByDatabase(mhandle, fragment_db,
- structure_db, torsion_sampler,
- min_loops_required=-1,
- max_res_extension=6)
-# if above fails, try DB-fill with less restrictions
-modelling.FillLoopsByDatabase(mhandle, fragment_db,
- structure_db, torsion_sampler,
- min_loops_required=-1)
-modelling.FillLoopsByDatabase(mhandle, fragment_db,
- structure_db, torsion_sampler)
-# if above fails on some gaps, use Monte Carlo
-modelling.FillLoopsByMonteCarlo(mhandle, torsion_sampler)
-# as a last resort, try to close large deletions
-modelling.CloseLargeDeletions(mhandle, structure_db)
+
+# perform loop modelling to close all gaps
+modelling.CloseGaps(mhandle)
+
# build sidechains
modelling.BuildSidechains(mhandle, merge_distance,
fragment_db, structure_db,
torsion_sampler)
+
# minimize energy of final model using molecular mechanics
modelling.MinimizeModelEnergy(mhandle)
+
# check final model and report issues
modelling.CheckFinalModel(mhandle)
# extract final model
final_model = mhandle.model
-io.SavePDB(final_model, 'model.pdb')
\ No newline at end of file
+io.SavePDB(final_model, 'model.pdb')
diff --git a/modelling/doc/pipeline.rst b/modelling/doc/pipeline.rst
index 06692d9e..c37d7e75 100644
--- a/modelling/doc/pipeline.rst
+++ b/modelling/doc/pipeline.rst
@@ -9,7 +9,9 @@ A protein homology modelling pipeline has the following main steps:
- Perform loop modelling to close (or remove) all gaps (see functions
:func:`CloseSmallDeletions`, :func:`RemoveTerminalGaps`,
:func:`MergeGapsByDistance`, :func:`FillLoopsByDatabase`,
- :func:`FillLoopsByMonteCarlo`, :func:`CloseLargeDeletions`)
+ :func:`FillLoopsByMonteCarlo`, :func:`CloseLargeDeletions` or
+ :func:`CloseGaps` that calls all these functions using predefined
+ heuristics)
- Build sidechains (see :func:`BuildSidechains` function)
- Minimize energy of final model using molecular mechanics
(see :func:`MinimizeModelEnergy` function)
@@ -23,13 +25,6 @@ purposes:
.. literalinclude:: ../../../tests/doc/scripts/modelling_steps.py
-In the default pipeline above, we call :func:`FillLoopsByDatabase` multiple
-times. First, we try to close "easy" gaps which require few extensions (we wish
-to limit the damage we do on the template) and for which we have plenty of loop
-candidates. If some gaps cannot be closed like this, we try less restrictive
-options. This approach is helpful if neighboring gaps are close together and the
-one closer to the C-terminus is easier to close. Several variants of the
-pipeline were evaluated on 1752 target-template-pairs and this one worked best.
Build Raw Modelling Handle
--------------------------------------------------------------------------------
@@ -207,6 +202,8 @@ Modelling Steps
:return: Number of gaps which were removed.
:rtype: :class:`int`
+.. autofunction:: CloseGaps
+
.. autofunction:: CloseSmallDeletions
.. autofunction:: MergeGapsByDistance
diff --git a/modelling/pymod/_closegaps.py b/modelling/pymod/_closegaps.py
index 7033ffd2..a30dad17 100644
--- a/modelling/pymod/_closegaps.py
+++ b/modelling/pymod/_closegaps.py
@@ -393,7 +393,8 @@ def CloseSmallDeletions(mhandle, max_extension=9, clash_thresh=1.0,
if not success:
current_gap_index += 1
-def MergeGapsByDistance(mhandle, distance):
+def MergeGapsByDistance(mhandle, distance, chain_idx = None,
+ resnum_range = None):
'''Merge 2 neighbouring gaps by deleting residues in-between.
Check if two neighbouring gaps are at max. *distance* residues apart from
@@ -411,6 +412,14 @@ def MergeGapsByDistance(mhandle, distance):
:param distance: The max. no. of residues between two gaps up to which
merge happens.
:type distance: :class:`int`
+
+ :param chain_idx: If not None, only gaps from chain with given index get
+ processed
+ :type chain_idx: :class:`int`
+
+ :param resnum_range: If not None, two gaps only get merged if they're
+ both in this resnum range.
+ :type resnum_range: :class:`tuple` containing two :class:`int`
'''
prof_name = 'closegaps::MergeGapsByDistance'
@@ -438,6 +447,10 @@ def MergeGapsByDistance(mhandle, distance):
# check that we are on the same chain
if current_gap.GetChain() != next_gap.GetChain():
continue
+ # check for range (if given)
+ if not (_InRange(current_gap, chain_idx, resnum_range) and\
+ _InRange(next_gap, chain_idx, resnum_range)):
+ continue
# no merging of gaps at the end AND the start :)
if current_gap.IsNTerminal() and next_gap.IsCTerminal():
continue
@@ -1013,7 +1026,8 @@ def ModelTermini(mhandle, torsion_sampler, fragger_handles=None,
def CloseLargeDeletions(mhandle, structure_db, linker_length=8,
num_fragments=500, use_scoring_extender=True,
- use_full_extender=True):
+ use_full_extender=True, chain_idx=None,
+ resnum_range=None):
'''Try to close large deletions.
@@ -1046,6 +1060,14 @@ def CloseLargeDeletions(mhandle, structure_db, linker_length=8,
:class:`GapExtender`.
See :func:`CloseSmallDeletions`.
:type use_full_extender: :class:`bool`
+
+ :param chain_idx: If not None, only gaps from chain with given index get
+ processed
+ :type chain_idx: :class:`int`
+
+ :param resnum_range: If not None, only gaps within this resnum range get
+ processed
+ :type resnum_range: :class:`tuple` containing two :class:`int`
'''
prof_name = 'closegaps::CloseLargeDeletions'
@@ -1067,6 +1089,12 @@ def CloseLargeDeletions(mhandle, structure_db, linker_length=8,
# keep copy of original gap
gap_orig = mhandle.gaps[gap_idx].Copy()
+
+ # check whether we are in the desired range
+ if not _InRange(gap_orig, chain_idx, resnum_range):
+ gap_idx += 1
+ continue
+
# get chain for gap
actual_chain_idx = gap_orig.GetChainIndex()
actual_chain = mhandle.model.chains[actual_chain_idx]
diff --git a/modelling/pymod/_pipeline.py b/modelling/pymod/_pipeline.py
index b31adee2..f7c3e2ec 100644
--- a/modelling/pymod/_pipeline.py
+++ b/modelling/pymod/_pipeline.py
@@ -162,6 +162,99 @@ def _GetSimEntity(sim):
return ent
###############################################################################
+def CloseGaps(mhandle, merge_distance=4, fragment_db=None, structure_db=None,
+ torsion_sampler=None, chain_idx = None, resnum_range = None):
+ '''
+ Tries to close all gaps in a model, except termini. It will go through
+ following steps:
+
+ - Try to close small deletions by relaxing them
+ (see :func:`CloseSmallDeletions`)
+ - Iteratively merge gaps up to a distance **merge_distance**
+ (see :func:`MergeGapsByDistance`) and try to fill them with a database
+ approach (see :func:`FillLoopsByDatabase`)
+ - Try to fill remaining gaps using a Monte Carlo approach
+ (see :func:`FillLoopsByMonteCarlo`)
+ - Large deletions get closed using a last resort approach
+ (see :func:`CloseLargeDeletions`)
+
+ :param mhandle: Modelling handle on which to apply change.
+ :type mhandle: :class:`ModellingHandle`
+
+ :param merge_distance: Max. merge distance when performing the database
+ approach
+ :type merge_distance: :class:`int`
+ :param fragment_db: Database for searching fragments in database
+ approach, must be consistent with provided
+ **structure_db**. A default is loaded if None.
+ :type fragment_db: :class:`~promod3.loop.FragDB`
+ :param structure_db: Structure db from which the **fragment_db** gets
+ it's structural information. A default is loaded
+ if None.
+ :type structure_db: :class:`~promod3.loop.StructureDB`
+ :param torsion_sampler: Used as parameter for :func:`FillLoopsByDatabase`
+ and :func:`FillLoopsByMonteCarlo` A default one is
+ loaded if None.
+ :param chain_idx: If not None, only gaps from chain with given index get
+ processed
+ :type chain_idx: :class:`int`
+
+ :param resnum_range: If not None, only gaps within this resnum range get
+ processed.
+ :type resnum_range: :class:`tuple` containing two :class:`int`
+
+ '''
+
+
+ # load stuff if needed
+ if not IsBackboneScorerSet(mhandle) or not IsBackboneScorerEnvSet(mhandle):
+ SetupDefaultBackboneScorer(mhandle)
+ if fragment_db is None:
+ fragment_db = loop.LoadFragDB()
+ if structure_db is None:
+ structure_db = loop.LoadStructureDB()
+ if torsion_sampler is None:
+ torsion_sampler = loop.LoadTorsionSamplerCoil()
+
+ #try to close small deletions by relaxing them
+ CloseSmallDeletions(mhandle, chain_idx = chain_idx,
+ resnum_range = resnum_range)
+
+ # iteratively merge gaps of distance i and fill loops by database
+ for distance in range(merge_distance):
+ MergeGapsByDistance(mhandle, distance, chain_idx = chain_idx,
+ resnum_range = resnum_range)
+ FillLoopsByDatabase(mhandle, fragment_db, structure_db,
+ torsion_sampler, min_loops_required=-1,
+ max_res_extension=6, chain_idx = chain_idx,
+ resnum_range = resnum_range)
+
+ # if above fails, try DB-fill with less restrictions
+ FillLoopsByDatabase(mhandle, fragment_db, structure_db,
+ torsion_sampler, min_loops_required=-1,
+ chain_idx = chain_idx, resnum_range = resnum_range)
+ FillLoopsByDatabase(mhandle, fragment_db, structure_db,
+ torsion_sampler, chain_idx = chain_idx,
+ resnum_range = resnum_range)
+
+ # close remaining gaps by Monte Carlo
+ FillLoopsByMonteCarlo(mhandle, torsion_sampler, chain_idx = chain_idx,
+ resnum_range = resnum_range)
+
+ # last resort approach to close large deletions
+ CloseLargeDeletions(mhandle, structure_db, chain_idx = chain_idx,
+ resnum_range = resnum_range)
+
+ #In the function above, we call :func:`FillLoopsByDatabase` multiple
+ #times. First, we try to close "easy" gaps which require few extensions
+ #(we wish to limit the damage we do on the template) and for which we have
+ #plenty of loop candidates. If some gaps cannot be closed like this, we try
+ #less restrictive options. This approach is helpful if neighboring gaps are
+ #close together and the one closer to the C-terminus is easier to close.
+ #Several variants were evaluated on 1752 target-template-pairs and this one
+ #worked best.
+
+
def BuildSidechains(mhandle, merge_distance=4, fragment_db=None,
structure_db=None, torsion_sampler=None):
'''Build sidechains for model.
@@ -467,24 +560,13 @@ def BuildFromRawModel(mhandle, use_amber_ff=False, extra_force_fields=list()):
else:
SetupDefaultBackboneScorer(mhandle)
- # remove terminal gaps and close small deletions
+ # remove terminal gaps
RemoveTerminalGaps(mhandle)
- CloseSmallDeletions(mhandle)
- # iteratively merge gaps of distance i and fill loops by database
- for distance in range(merge_distance):
- MergeGapsByDistance(mhandle, distance)
- FillLoopsByDatabase(mhandle, fragment_db, structure_db,
- torsion_sampler, min_loops_required=-1,
- max_res_extension=6)
- FillLoopsByDatabase(mhandle, fragment_db, structure_db,
- torsion_sampler, min_loops_required=-1)
- FillLoopsByDatabase(mhandle, fragment_db, structure_db,
- torsion_sampler)
-
- # close remaining gaps by Monte Carlo
- FillLoopsByMonteCarlo(mhandle, torsion_sampler)
- CloseLargeDeletions(mhandle, structure_db)
+ #close gaps
+ CloseGaps(mhandle, merge_distance = merge_distance,
+ fragment_db = fragment_db, structure_db = structure_db,
+ torsion_sampler = torsion_sampler)
# build sidechains
BuildSidechains(mhandle, merge_distance, fragment_db,
@@ -501,5 +583,5 @@ def BuildFromRawModel(mhandle, use_amber_ff=False, extra_force_fields=list()):
return mhandle.model
# these methods will be exported into module
-__all__ = ('BuildFromRawModel', 'BuildSidechains', 'MinimizeModelEnergy',
- 'CheckFinalModel')
+__all__ = ('BuildFromRawModel', 'CloseGaps','BuildSidechains',
+ 'MinimizeModelEnergy', 'CheckFinalModel')
--
GitLab