|
|
|
|
@@ -1,23 +1,33 @@
|
|
|
|
|
"""Functions to manipulate Determinant objects.
|
|
|
|
|
|
|
|
|
|
from __future__ import division
|
|
|
|
|
from __future__ import print_function
|
|
|
|
|
|
|
|
|
|
import math, time
|
|
|
|
|
|
|
|
|
|
import propka.iterative, propka.lib, propka.vector_algebra
|
|
|
|
|
import propka.calculations
|
|
|
|
|
TODO - it is confusing to have both `determinant.py` and `determinants.py`.
|
|
|
|
|
Should these be merged?
|
|
|
|
|
"""
|
|
|
|
|
import math
|
|
|
|
|
import propka.iterative
|
|
|
|
|
import propka.lib
|
|
|
|
|
import propka.vector_algebra
|
|
|
|
|
from propka.calculations import squared_distance, get_smallest_distance
|
|
|
|
|
from propka.calculations import angle_distance_factors, hydrogen_bond_energy
|
|
|
|
|
from propka.determinant import Determinant
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
# TODO - it is confusing to have both `determinant.py` and `determinants.py`.
|
|
|
|
|
# Should these be merged?
|
|
|
|
|
# Cutoff for angle factor
|
|
|
|
|
# TODO - this constant appears elsewhere in the package.
|
|
|
|
|
# It should be moved to a configuration file.
|
|
|
|
|
FANGLE_MIN = 0.001
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def setDeterminants(propka_groups, version=None, options=None):
|
|
|
|
|
def set_determinants(propka_groups, version=None, options=None):
|
|
|
|
|
"""Add side-chain and coulomb determinants/perturbations to all residues.
|
|
|
|
|
|
|
|
|
|
NOTE - backbone determinants are set separately
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
propka_groups: groups to adjust
|
|
|
|
|
version: version object
|
|
|
|
|
options: options object
|
|
|
|
|
"""
|
|
|
|
|
adding side-chain and coulomb determinants/perturbations to all residues - note, backbone determinants are set separately
|
|
|
|
|
"""
|
|
|
|
|
|
|
|
|
|
iterative_interactions = []
|
|
|
|
|
# --- NonIterative section ---#
|
|
|
|
|
for group1 in propka_groups:
|
|
|
|
|
@@ -27,172 +37,198 @@ def setDeterminants(propka_groups, version=None, options=None):
|
|
|
|
|
# do not calculate interactions for coupled groups
|
|
|
|
|
if group2 in group1.covalently_coupled_groups:
|
|
|
|
|
break
|
|
|
|
|
|
|
|
|
|
distance = propka.calculations.distance(group1, group2)
|
|
|
|
|
|
|
|
|
|
if distance < version.parameters.coulomb_cutoff2:
|
|
|
|
|
interaction_type = version.parameters.interaction_matrix.get_value(group1.type,group2.type)
|
|
|
|
|
interaction_type = version.parameters.interaction_matrix.get_value(group1.type,
|
|
|
|
|
group2.type)
|
|
|
|
|
if interaction_type == 'I':
|
|
|
|
|
propka.iterative.addtoDeterminantList(group1, group2, distance, iterative_interactions, version=version)
|
|
|
|
|
propka.iterative.addtoDeterminantList(group1, group2,
|
|
|
|
|
distance,
|
|
|
|
|
iterative_interactions,
|
|
|
|
|
version=version)
|
|
|
|
|
elif interaction_type == 'N':
|
|
|
|
|
addDeterminants(group1, group2, distance, version)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
add_determinants(group1, group2, distance, version)
|
|
|
|
|
# --- Iterative section ---#
|
|
|
|
|
propka.iterative.addDeterminants(iterative_interactions, version, options=options)
|
|
|
|
|
propka.iterative.add_determinants(iterative_interactions, version,
|
|
|
|
|
options=options)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def addDeterminants(group1, group2, distance, version):
|
|
|
|
|
def add_determinants(group1, group2, distance, version):
|
|
|
|
|
"""Add determinants and perturbations for distance(R1, R2) < coulomb_cutoff.
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
group1: first group to add
|
|
|
|
|
group2: second group to add
|
|
|
|
|
distance: distance between groups
|
|
|
|
|
version: version object
|
|
|
|
|
"""
|
|
|
|
|
adding determinants/perturbations, distance(R1, R2) < coulomb_cutoff always
|
|
|
|
|
"""
|
|
|
|
|
|
|
|
|
|
# side-chain determinant
|
|
|
|
|
addSidechainDeterminants(group1, group2, version)
|
|
|
|
|
|
|
|
|
|
add_sidechain_determinants(group1, group2, version)
|
|
|
|
|
# Coulomb determinant
|
|
|
|
|
addCoulombDeterminants(group1, group2, distance, version)
|
|
|
|
|
add_coulomb_determinants(group1, group2, distance, version)
|
|
|
|
|
|
|
|
|
|
return
|
|
|
|
|
|
|
|
|
|
def addSidechainDeterminants(group1, group2, version=None):
|
|
|
|
|
def add_sidechain_determinants(group1, group2, version=None):
|
|
|
|
|
"""Add side-chain determinants and perturbations.
|
|
|
|
|
|
|
|
|
|
NOTE - res_num1 > res_num2
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
group1: first group to add
|
|
|
|
|
group2: second group to add
|
|
|
|
|
version: version object
|
|
|
|
|
"""
|
|
|
|
|
adding side-chain determinants/perturbations
|
|
|
|
|
Note, res_num1 > res_num2
|
|
|
|
|
"""
|
|
|
|
|
|
|
|
|
|
hbond_interaction = version.hydrogen_bond_interaction(group1, group2)
|
|
|
|
|
|
|
|
|
|
if hbond_interaction:
|
|
|
|
|
|
|
|
|
|
if group1.charge == group2.charge:
|
|
|
|
|
# acid pair or base pair
|
|
|
|
|
if group1.model_pka < group2.model_pka:
|
|
|
|
|
newDeterminant1 = Determinant(group2, -hbond_interaction)
|
|
|
|
|
newDeterminant2 = Determinant(group1, hbond_interaction)
|
|
|
|
|
new_determinant1 = Determinant(group2, -hbond_interaction)
|
|
|
|
|
new_determinant2 = Determinant(group1, hbond_interaction)
|
|
|
|
|
else:
|
|
|
|
|
newDeterminant1 = Determinant(group2, hbond_interaction)
|
|
|
|
|
newDeterminant2 = Determinant(group1, -hbond_interaction)
|
|
|
|
|
new_determinant1 = Determinant(group2, hbond_interaction)
|
|
|
|
|
new_determinant2 = Determinant(group1, -hbond_interaction)
|
|
|
|
|
else:
|
|
|
|
|
newDeterminant1 = Determinant(group2, hbond_interaction*group1.charge)
|
|
|
|
|
newDeterminant2 = Determinant(group1, hbond_interaction*group2.charge)
|
|
|
|
|
new_determinant1 = Determinant(group2, hbond_interaction*group1.charge)
|
|
|
|
|
new_determinant2 = Determinant(group1, hbond_interaction*group2.charge)
|
|
|
|
|
group1.determinants['sidechain'].append(new_determinant1)
|
|
|
|
|
group2.determinants['sidechain'].append(new_determinant2)
|
|
|
|
|
|
|
|
|
|
group1.determinants['sidechain'].append(newDeterminant1)
|
|
|
|
|
group2.determinants['sidechain'].append(newDeterminant2)
|
|
|
|
|
|
|
|
|
|
return
|
|
|
|
|
def add_coulomb_determinants(group1, group2, distance, version):
|
|
|
|
|
"""Add non-iterative Coulomb determinants and perturbations.
|
|
|
|
|
|
|
|
|
|
def addCoulombDeterminants(group1, group2, distance, version):
|
|
|
|
|
Args:
|
|
|
|
|
group1: first group to add
|
|
|
|
|
group2: second group to add
|
|
|
|
|
distance: distance between groups
|
|
|
|
|
version: version object
|
|
|
|
|
"""
|
|
|
|
|
adding NonIterative Coulomb determinants/perturbations
|
|
|
|
|
"""
|
|
|
|
|
|
|
|
|
|
coulomb_interaction = version.electrostatic_interaction(group1, group2, distance)
|
|
|
|
|
|
|
|
|
|
coulomb_interaction = version.electrostatic_interaction(group1, group2,
|
|
|
|
|
distance)
|
|
|
|
|
if coulomb_interaction:
|
|
|
|
|
Q1 = group1.charge
|
|
|
|
|
Q2 = group2.charge
|
|
|
|
|
|
|
|
|
|
q1 = group1.charge
|
|
|
|
|
q2 = group2.charge
|
|
|
|
|
# assigning the Coulombic interaction
|
|
|
|
|
if Q1 < 0.0 and Q2 < 0.0:
|
|
|
|
|
""" both are acids """
|
|
|
|
|
addCoulombAcidPair(group1, group2, coulomb_interaction)
|
|
|
|
|
elif Q1 > 0.0 and Q2 > 0.0:
|
|
|
|
|
""" both are bases """
|
|
|
|
|
addCoulombBasePair(group1, group2, coulomb_interaction)
|
|
|
|
|
if q1 < 0.0 and q2 < 0.0:
|
|
|
|
|
# both are acids
|
|
|
|
|
add_coulomb_acid_pair(group1, group2, coulomb_interaction)
|
|
|
|
|
elif q1 > 0.0 and q2 > 0.0:
|
|
|
|
|
# both are bases
|
|
|
|
|
add_coulomb_base_pair(group1, group2, coulomb_interaction)
|
|
|
|
|
else:
|
|
|
|
|
""" one of each """
|
|
|
|
|
addCoulombIonPair(group1, group2, coulomb_interaction)
|
|
|
|
|
|
|
|
|
|
return
|
|
|
|
|
# one of each
|
|
|
|
|
add_coulomb_ion_pair(group1, group2, coulomb_interaction)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def addCoulombAcidPair(object1, object2, value):
|
|
|
|
|
def add_coulomb_acid_pair(object1, object2, value):
|
|
|
|
|
"""Add the Coulomb interaction (an acid pair).
|
|
|
|
|
|
|
|
|
|
The higher pKa is raised.
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
object1: first part of pair
|
|
|
|
|
object2: second part of pair
|
|
|
|
|
value: determinant value
|
|
|
|
|
"""
|
|
|
|
|
Adding the Coulomb interaction (an acid pair):
|
|
|
|
|
the higher pKa is raised
|
|
|
|
|
"""
|
|
|
|
|
|
|
|
|
|
if object1.model_pka > object2.model_pka:
|
|
|
|
|
newDeterminant = Determinant(object2, value)
|
|
|
|
|
object1.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
new_determinant = Determinant(object2, value)
|
|
|
|
|
object1.determinants['coulomb'].append(new_determinant)
|
|
|
|
|
else:
|
|
|
|
|
newDeterminant = Determinant(object1, value)
|
|
|
|
|
object2.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
new_determinant = Determinant(object1, value)
|
|
|
|
|
object2.determinants['coulomb'].append(new_determinant)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def addCoulombBasePair(object1, object2, value):
|
|
|
|
|
"""
|
|
|
|
|
Adding the Coulomb interaction (a base pair):
|
|
|
|
|
the lower pKa is lowered
|
|
|
|
|
def add_coulomb_base_pair(object1, object2, value):
|
|
|
|
|
"""Add the Coulomb interaction (a base pair).
|
|
|
|
|
|
|
|
|
|
The lower pKa is lowered.
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
object1: first part of pair
|
|
|
|
|
object2: second part of pair
|
|
|
|
|
value: determinant value
|
|
|
|
|
"""
|
|
|
|
|
if object1.model_pka < object2.model_pka:
|
|
|
|
|
newDeterminant = Determinant(object2, -value)
|
|
|
|
|
object1.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
new_determinant = Determinant(object2, -value)
|
|
|
|
|
object1.determinants['coulomb'].append(new_determinant)
|
|
|
|
|
else:
|
|
|
|
|
newDeterminant = Determinant(object1, -value)
|
|
|
|
|
object2.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
new_determinant = Determinant(object1, -value)
|
|
|
|
|
object2.determinants['coulomb'].append(new_determinant)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def addCoulombIonPair(object1, object2, value):
|
|
|
|
|
def add_coulomb_ion_pair(object1, object2, value):
|
|
|
|
|
"""Add the Coulomb interaction (an acid-base pair).
|
|
|
|
|
|
|
|
|
|
The pKa of the acid is lowered & the pKa of the base is raised.
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
object1: first part of pair
|
|
|
|
|
object2: second part of pair
|
|
|
|
|
value: determinant value
|
|
|
|
|
"""
|
|
|
|
|
Adding the Coulomb interaction (an acid-base pair):
|
|
|
|
|
the pKa of the acid is lowered & the pKa of the base is raised
|
|
|
|
|
"""
|
|
|
|
|
|
|
|
|
|
# residue1
|
|
|
|
|
Q1 = object1.charge
|
|
|
|
|
newDeterminant = Determinant(object2, Q1*value)
|
|
|
|
|
object1.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
|
|
|
|
|
q1 = object1.charge
|
|
|
|
|
new_determinant = Determinant(object2, q1*value)
|
|
|
|
|
object1.determinants['coulomb'].append(new_determinant)
|
|
|
|
|
# residue2
|
|
|
|
|
Q2 = object2.charge
|
|
|
|
|
newDeterminant = Determinant(object1, Q2*value)
|
|
|
|
|
object2.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
q2 = object2.charge
|
|
|
|
|
new_determinant = Determinant(object1, q2*value)
|
|
|
|
|
object2.determinants['coulomb'].append(new_determinant)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def set_ion_determinants(conformation_container, version):
|
|
|
|
|
"""Add ion determinants and perturbations.
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
def setIonDeterminants(conformation_container, version):
|
|
|
|
|
"""
|
|
|
|
|
adding ion determinants/perturbations
|
|
|
|
|
Args:
|
|
|
|
|
conformation_container: conformation to set
|
|
|
|
|
version: version object
|
|
|
|
|
"""
|
|
|
|
|
for titratable_group in conformation_container.get_titratable_groups():
|
|
|
|
|
for ion_group in conformation_container.get_ions():
|
|
|
|
|
squared_distance = propka.calculations.squared_distance(titratable_group, ion_group)
|
|
|
|
|
if squared_distance < version.parameters.coulomb_cutoff2_squared:
|
|
|
|
|
weight = version.calculate_pair_weight(titratable_group.Nmass, ion_group.Nmass)
|
|
|
|
|
# the pKa of both acids and bases are shifted up by negative ions (and vice versa)
|
|
|
|
|
value = (-ion_group.charge) * version.calculate_coulomb_energy(math.sqrt(squared_distance), weight)
|
|
|
|
|
newDeterminant = Determinant(ion_group, value)
|
|
|
|
|
titratable_group.determinants['coulomb'].append(newDeterminant)
|
|
|
|
|
dist_sq = squared_distance(titratable_group, ion_group)
|
|
|
|
|
if dist_sq < version.parameters.coulomb_cutoff2_squared:
|
|
|
|
|
weight = version.calculate_pair_weight(titratable_group.Nmass,
|
|
|
|
|
ion_group.Nmass)
|
|
|
|
|
# the pKa of both acids and bases are shifted up by negative
|
|
|
|
|
# ions (and vice versa)
|
|
|
|
|
value = (-ion_group.charge) \
|
|
|
|
|
* version.calculate_coulomb_energy(math.sqrt(dist_sq),
|
|
|
|
|
weight)
|
|
|
|
|
new_det = Determinant(ion_group, value)
|
|
|
|
|
titratable_group.determinants['coulomb'].append(new_det)
|
|
|
|
|
|
|
|
|
|
return
|
|
|
|
|
|
|
|
|
|
def setBackBoneDeterminants(titratable_groups, backbone_groups, version):
|
|
|
|
|
def set_backbone_determinants(titratable_groups, backbone_groups, version):
|
|
|
|
|
"""Set determinants between titrable and backbone groups.
|
|
|
|
|
|
|
|
|
|
Args:
|
|
|
|
|
titratable_groups: list of titratable groups
|
|
|
|
|
backbone_groups: list of backbone groups
|
|
|
|
|
version: version object
|
|
|
|
|
"""
|
|
|
|
|
for titratable_group in titratable_groups:
|
|
|
|
|
titratable_group_interaction_atoms = titratable_group.interaction_atoms_for_acids
|
|
|
|
|
titratable_group_interaction_atoms \
|
|
|
|
|
= titratable_group.interaction_atoms_for_acids
|
|
|
|
|
if not titratable_group_interaction_atoms:
|
|
|
|
|
continue
|
|
|
|
|
|
|
|
|
|
# find out which backbone groups this titratable is interacting with
|
|
|
|
|
for backbone_group in backbone_groups:
|
|
|
|
|
# find the interacting atoms
|
|
|
|
|
backbone_interaction_atoms = backbone_group.get_interaction_atoms(titratable_group)
|
|
|
|
|
backbone_interaction_atoms \
|
|
|
|
|
= backbone_group.get_interaction_atoms(titratable_group)
|
|
|
|
|
if not backbone_interaction_atoms:
|
|
|
|
|
continue
|
|
|
|
|
|
|
|
|
|
# find the smallest distance
|
|
|
|
|
[backbone_atom, distance, titratable_atom] = propka.calculations.get_smallest_distance(backbone_interaction_atoms,
|
|
|
|
|
titratable_group_interaction_atoms)
|
|
|
|
|
[backbone_atom, distance, titratable_atom] \
|
|
|
|
|
= get_smallest_distance(backbone_interaction_atoms, \
|
|
|
|
|
titratable_group_interaction_atoms)
|
|
|
|
|
# get the parameters
|
|
|
|
|
parameters = version.get_backbone_hydrogen_bond_parameters(backbone_atom, titratable_atom)
|
|
|
|
|
parameters = version.get_backbone_hydrogen_bond_parameters(backbone_atom,
|
|
|
|
|
titratable_atom)
|
|
|
|
|
if not parameters:
|
|
|
|
|
continue
|
|
|
|
|
[dpKa_max, [cutoff1, cutoff2]] = parameters
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
[dpka_max, [cutoff1, cutoff2]] = parameters
|
|
|
|
|
if distance < cutoff2:
|
|
|
|
|
# calculate angle factor
|
|
|
|
|
f_angle = 1.0
|
|
|
|
|
@@ -206,19 +242,20 @@ def setBackBoneDeterminants(titratable_groups, backbone_groups, version):
|
|
|
|
|
# ||
|
|
|
|
|
# C
|
|
|
|
|
if backbone_group.type == 'BBC':
|
|
|
|
|
if titratable_group.type in version.parameters.angular_dependent_sidechain_interactions:
|
|
|
|
|
if titratable_group.type \
|
|
|
|
|
in version.parameters.angular_dependent_sidechain_interactions:
|
|
|
|
|
if titratable_atom.element == 'H':
|
|
|
|
|
heavy_atom = titratable_atom.bonded_atoms[0]
|
|
|
|
|
heavy_atom = titratable_atom.bonded_atoms[0]
|
|
|
|
|
hydrogen_atom = titratable_atom
|
|
|
|
|
[d1, f_angle, d2] = propka.calculations.angle_distance_factors(atom1=heavy_atom,
|
|
|
|
|
atom2=hydrogen_atom,
|
|
|
|
|
atom3=backbone_atom)
|
|
|
|
|
[_, f_angle, _] = angle_distance_factors(atom1=heavy_atom,
|
|
|
|
|
atom2=hydrogen_atom,
|
|
|
|
|
atom3=backbone_atom)
|
|
|
|
|
else:
|
|
|
|
|
# Either the structure is corrupt (no hydrogen), or the heavy atom is closer to
|
|
|
|
|
# the titratable atom than the hydrogen. In either case we set the angle factor
|
|
|
|
|
# to 0
|
|
|
|
|
# Either the structure is corrupt (no hydrogen),
|
|
|
|
|
# or the heavy atom is closer to the titratable
|
|
|
|
|
# atom than the hydrogen. In either case we set the
|
|
|
|
|
# angle factor to 0
|
|
|
|
|
f_angle = 0.0
|
|
|
|
|
|
|
|
|
|
# for BBN groups, the hydrogen is on the backbone group
|
|
|
|
|
#
|
|
|
|
|
# Titra.
|
|
|
|
|
@@ -229,23 +266,21 @@ def setBackBoneDeterminants(titratable_groups, backbone_groups, version):
|
|
|
|
|
# / \
|
|
|
|
|
if backbone_group.type == 'BBN':
|
|
|
|
|
if backbone_atom.element == 'H':
|
|
|
|
|
backbone_N = backbone_atom.bonded_atoms[0]
|
|
|
|
|
backbone_H = backbone_atom
|
|
|
|
|
[d1, f_angle, d2] = propka.calculations.angle_distance_factors(atom1=titratable_atom,
|
|
|
|
|
atom2=backbone_H,
|
|
|
|
|
atom3=backbone_N)
|
|
|
|
|
backbone_n = backbone_atom.bonded_atoms[0]
|
|
|
|
|
backbone_h = backbone_atom
|
|
|
|
|
[_, f_angle, _] = angle_distance_factors(atom1=titratable_atom,
|
|
|
|
|
atom2=backbone_h,
|
|
|
|
|
atom3=backbone_n)
|
|
|
|
|
else:
|
|
|
|
|
# Either the structure is corrupt (no hydrogen), or the heavy atom is closer to
|
|
|
|
|
# the titratable atom than the hydrogen. In either case we set the angle factor
|
|
|
|
|
# to 0
|
|
|
|
|
# Either the structure is corrupt (no hydrogen), or the
|
|
|
|
|
# heavy atom is closer to the titratable atom than the
|
|
|
|
|
# hydrogen. In either case we set the angle factor to 0
|
|
|
|
|
f_angle = 0.0
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
if f_angle > 0.001:
|
|
|
|
|
value = titratable_group.charge * propka.calculations.hydrogen_bond_energy(distance, dpKa_max, [cutoff1,cutoff2], f_angle)
|
|
|
|
|
|
|
|
|
|
newDeterminant = Determinant(backbone_group, value)
|
|
|
|
|
titratable_group.determinants['backbone'].append(newDeterminant)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
return
|
|
|
|
|
if f_angle > FANGLE_MIN:
|
|
|
|
|
value = titratable_group.charge * hydrogen_bond_energy(distance,
|
|
|
|
|
dpka_max,
|
|
|
|
|
[cutoff1, cutoff2],
|
|
|
|
|
f_angle)
|
|
|
|
|
new_determinant = Determinant(backbone_group, value)
|
|
|
|
|
titratable_group.determinants['backbone'].append(new_determinant)
|
|
|
|
|
|