Skip to content

Commit

Permalink
Merge pull request #422 from xikunliu/region_function
Browse files Browse the repository at this point in the history 
include region
  • Loading branch information
@pckroon
pckroon committed May 23, 2022
2 parents 4182683 + 61125ee commit 69fe046
Show file tree
Hide file tree
Showing 3 changed files with 229 additions and 4 deletions.
19 changes: 15 additions & 4 deletions bin/martinize2
View file
Original file line number Diff line number Diff line change
Expand Up @@ -333,7 +333,6 @@ def maxwarn(value):
return (splitted[0], count)
raise argparse.ArgumentTypeError(msg)


def entry():
"""
Parses commandline arguments and performs the logic.
Expand Down Expand Up @@ -441,14 +440,16 @@ def entry():
help=('Establish what is the structural unit for the '
'elastic network. Bonds are only created within'
' a unit.'))

rb_group.add_argument('-region', dest='rb_region', default='molecule',
help=('Apply rubber bands only within the specified regions. '
'The format is <start_resid_1>:<end_resid_1>, <start_resid_2>:<end_resid_2>...'
'Start and end resids are included'))
go_group = parser.add_argument_group('Virtual site based GoMartini')
go_group.add_argument('-govs-includes', action='store_true', default=False,
help='Write include statements to use Vitrual Site Go Martini.')
go_group.add_argument('-govs-moltype', default='molecule_0',
help=('Set the name of the molecule when using '
'Virtual Sites GoMartini.'))

prot_group = parser.add_argument_group('Protein description')
prot_group.add_argument('-scfix', dest='scfix',
action='store_true', default=False,
Expand Down Expand Up @@ -726,6 +727,16 @@ def entry():
message = 'Unknown value for -eunit: "{}".'.format(args.rb_unit)
LOGGER.critical(message)
raise ValueError(message)
if args.rb_region == 'molecule':
domain_criterion = vermouth.processors.apply_rubber_band.always_true
else:
regions = [tuple([int(i) for i in apair.split(':')]) for apair in args.rb_region.split(',')]
if any(len(region) != 2 for region in regions):
message = 'Unknown value for -region: "{}".'.format(args.rb_region)
LOGGER.critical(message)
raise ValueError(message)
else:
domain_criterion = vermouth.processors.apply_rubber_band.make_same_region_criterion(regions)
if args.rb_selection is not None:
selector = functools.partial(
selectors.proto_select_attribute_in,
Expand All @@ -743,7 +754,7 @@ def entry():
minimum_force=args.rb_minimum_force,
selector=selector,
domain_criterion=domain_criterion,
res_min_dist=args.res_min_dist
res_min_dist=args.res_min_dist,
)
rubber_band_processor.run_system(system)

Expand Down
40 changes: 40 additions & 0 deletions vermouth/processors/apply_rubber_band.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,11 +18,13 @@

import numpy as np
import networkx as nx
import copy

from .processor import Processor
from .. import selectors
from ..graph_utils import make_residue_graph


# the bond type of the RB
DEFAULT_BOND_TYPE = 6
# the minimum distance between the resids
Expand Down Expand Up @@ -365,7 +367,45 @@ def same_chain(graph, left, right):
node_right = graph.nodes[right]
return node_left.get('chain') == node_right.get('chain')

def make_same_region_criterion(regions):
"""
Returns ``True`` is the nodes are part of the same region.
Nodes are considered part of the same region if their value
under the "resid" attribute are within the same residue range.
Parameters
----------
graph: networkx.Graph
A graph the nodes are part of.
left:
A node key in 'graph'.
right:
A node key in 'graph'.
regions:
[(resid_start_1,resid_end_1),(resid_start_2,resid_end_2),...] resid_start and resid_end are included)
Returns
-------
bool
``True`` if the nodes are part of the same region.
"""

regions = copy.deepcopy(regions)

def same_region(graph, left, right):
node_left = graph.nodes[left]
node_right = graph.nodes[right]
left_resid = node_left.get('resid')
right_resid = node_right.get('resid')
for region in regions:
lower = min(region)
upper = max(region)
if lower <= left_resid <= upper and lower <= right_resid <= upper:
return True
return False
return same_region

class ApplyRubberBand(Processor):
"""
Add an elastic network to a system between particles fulfilling the
Expand Down
174 changes: 174 additions & 0 deletions vermouth/tests/test_apply_rubber_band.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@
from vermouth import selectors
from vermouth.processors import apply_rubber_band
from vermouth.processors.apply_rubber_band import (same_chain,
make_same_region_criterion,
are_connected,
build_connectivity_matrix)

Expand Down Expand Up @@ -157,8 +158,107 @@ def test_same_chain(nodes, edges, chain, outcome):
graph.add_nodes_from(nodes)
graph.add_edges_from(edges)
nx.set_node_attributes(graph, chain, "chain")
print (same_chain(graph, 1, 2))
assert same_chain(graph, 1, 2) == outcome

@pytest.mark.parametrize('regions, left, right, nodes, chain, resid, edges, outcome', (
([(2, 3)],
2,
3,
[1, 2, 3, 4],
{1: "A", 2: "A", 3: "A", 4: "A"},
{1: 1, 2: 2, 3: 3, 4: 4},
[(1, 2), (2, 3), (3, 4)],
True),
([(2, 3)],
1,
2,
[1, 2, 3, 4],
{1: "A", 2: "A", 3: "A", 4: "A"},
{1: 1, 2: 2, 3: 3, 4: 4},
[(1, 2), (2, 3), (3, 4)],
False),
([(3, 2)],
2,
3,
[1, 2, 3, 4],
{1: "A", 2: "A", 3: "A", 4: "A"},
{1: 1, 2: 2, 3: 3, 4: 4},
[(1, 2), (2, 3), (3, 4)],
True),
([(2, 3)],
2,
3,
[1, 2, 3, 4],
{1: "A", 2: "A", 3: "B", 4: "B"},
{1: 1, 2: 2, 3: 3, 4: 4},
[(1, 2), (2, 3), (3, 4)],
True),
([(2, 3)],
1,
2,
[1, 2, 3, 4],
{1: "A", 2: "A", 3: "B", 4: "B"},
{1: 1, 2: 2, 3: 3, 4: 4},
[(1, 2), (2, 3), (3, 4)],
False),
([(2, 3)],
2,
5,
[1, 2, 3, 4, 5, 6],
{1: "A", 2: "A", 3: "A", 4: "B", 5: "B", 6: "B"},
{1: 1, 2: 2, 3: 3, 4: 1, 5: 2, 6: 3},
[(1, 2), (2, 3), (3, 4)],
True),
([(2, 3)],
5,
6,
[1, 2, 3, 4, 5, 6],
{1: "A", 2: "A", 3: "A", 4: "B", 5: "B", 6: "B"},
{1: 1, 2: 2, 3: 3, 4: 1, 5: 2, 6: 3},
[(1, 2), (2, 3), (3, 4)],
True),
([(2, 3)],
2,
6,
[1, 2, 3, 4, 5, 6],
{1: "A", 2: "A", 3: "A", 4: "B", 5: "B", 6: "B"},
{1: 1, 2: 2, 3: 3, 4: 1, 5: 2, 6: 3},
[(1, 2), (2, 3), (3, 4)],
True),
([(2, 3)],
2,
4,
[1, 2, 3, 4, 5, 6],
{1: "A", 2: "A", 3: "A", 4: "B", 5: "B", 6: "B"},
{1: 1, 2: 2, 3: 3, 4: 1, 5: 2, 6: 3},
[(1, 2), (2, 3), (3, 4)],
False),
([(3, 3)],
3,
3,
[1, 2, 3, 4, 5, 6],
{1: "A", 2: "A", 3: "A", 4: "B", 5: "B", 6: "B"},
{1: 1, 2: 2, 3: 3, 4: 1, 5: 2, 6: 3},
[(1, 2), (2, 3), (3, 4)],
True),
([(3, 3)],
3,
6,
[1, 2, 3, 4, 5, 6],
{1: "A", 2: "A", 3: "A", 4: "B", 5: "B", 6: "B"},
{1: 1, 2: 2, 3: 3, 4: 1, 5: 2, 6: 3},
[(1, 2), (2, 3), (3, 4)],
True),
))
def test_make_same_region_criterion(regions, left, right, nodes, edges, chain, resid, outcome):
graph = nx.Graph()
graph.add_nodes_from(nodes)
graph.add_edges_from(edges)
nx.set_node_attributes(graph, chain, "chain")
nx.set_node_attributes(graph, resid, "resid")
same_region = make_same_region_criterion(regions)
assert same_region(graph=graph, left=left, right=right) == outcome

@pytest.fixture
def test_molecule():
Expand Down Expand Up @@ -304,3 +404,77 @@ def test_apply_rubber_bands(test_molecule, chain_attribute, atom_names, res_min_
res_min_dist=res_min_dist)
process.run_molecule(test_molecule)
assert test_molecule.interactions['bonds'] == outcome


@pytest.mark.parametrize('regions, chain_attribute, atom_names, res_min_dist, outcome',
(([(1, 4)],
{0: 'A', 1: 'A', 2: 'A',
3: 'A', 4: 'A', 5: 'A',
6: 'A', 7: 'A', 8: 'A'},
['BB'],
2,
[vermouth.molecule.Interaction(
atoms=(0, 6),
meta={'group': 'Rubber band'},
parameters=[6, 1.5, 1000])]
),
# different min_res and different regions
([(1, 3)],
{0: 'A', 1: 'A', 2: 'A',
3: 'A', 4: 'A', 5: 'A',
6: 'A', 7: 'A', 8: 'A'},
['BB'],
1,
[vermouth.molecule.Interaction(
atoms=(0, 5),
meta={'group': 'Rubber band'},
parameters=[6, 1.0, 1000])]
),
# select more than only BB atoms
([(1, 3)],
{0: 'A', 1: 'A', 2: 'A',
3: 'A', 4: 'A', 5: 'A',
6: 'A', 7: 'A', 8: 'A'},
['BB', 'SC1'],
2,
[]
),
# change chain identifier
([(1, 3)],
{0: 'A', 1: 'A', 2: 'A',
3: 'B', 4: 'B', 5: 'B',
6: 'B', 7: 'B', 8: 'B'},
['BB'],
1,
[vermouth.molecule.Interaction(
atoms=(0, 5),
meta={'group': 'Rubber band'},
parameters=[6, 1.0, 1000])]
)))
def test_apply_rubber_bands_same_regions(test_molecule, regions, chain_attribute, atom_names, res_min_dist, outcome):
"""
Takes molecule and sets the chain attributes. Based on chain, minimum distance
between residues, and atom names elagible it is tested if rubber bands are applied
for the correct atoms in molecule.
"""
selector = functools.partial(
selectors.proto_select_attribute_in,
attribute='atomname',
values=atom_names)

domain_criterion = vermouth.processors.apply_rubber_band.make_same_region_criterion(regions)
nx.set_node_attributes(test_molecule, chain_attribute, 'chain')

process = vermouth.processors.apply_rubber_band.ApplyRubberBand(
selector=selector,
lower_bound=0.0,
upper_bound=10.,
decay_factor=0,
decay_power=0.,
base_constant=1000,
minimum_force=1,
bond_type=6,
domain_criterion=domain_criterion,
res_min_dist=res_min_dist)
process.run_molecule(test_molecule)
assert test_molecule.interactions['bonds'] == outcome

0 comments on commit 69fe046

Please sign in to comment.