Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Correct the judgment of hydrogen bonds #370

Open
wants to merge 10 commits into
base: master
Choose a base branch
from
60 changes: 30 additions & 30 deletions graphein/protein/edges/distance.py
Original file line number Diff line number Diff line change
Expand Up @@ -334,8 +334,30 @@ def add_disulfide_interactions(


def add_hydrogen_bond_interactions(
G: nx.Graph, rgroup_df: Optional[pd.DataFrame] = None
G: nx.Graph, rgroup_df: Optional[pd.DataFrame] = None
):
rgroup_df = G.graph["raw_pdb_df"]
HBOND_ATOMS = [
"ND1", # histidine
"ND2", # asparagine
"NE", # arginine
"NE1", # tryptophan
"NE2", # glutamate, histidine
"NH1", # arginine
"NH2", # arginine
"NZ", # lysine
"OD1", # asparagine, aspartic
"OD2", # aspartic
"OE1", # glutamate
"OE2", # glutamate
"OG", # serine
"OG1", # threonine
"OH", # tyrosine, threonine, serine
"SD", # methionine
"SG", # cysteine
"N", # backbone
"O", # backbone
]
DONORS_ATOM: List[str] = [
"ND2", # asparagine
"NE", # arginine
Expand All @@ -361,32 +383,11 @@ def add_hydrogen_bond_interactions(
]
"""Add all hydrogen-bond interactions."""
# For these atoms, find those that are within 3.5A of one another.
if rgroup_df is None:
rgroup_df = G.graph["rgroup_df"]
# if rgroup_df is None:
# rgroup_df = G.graph["rgroup_df"]
Copy link
Owner

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is the purpose of not using the supplied rgroup_df?

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, you're right. I commented out the code like that just to facilitate debugging and determine the hydrogen bonds on the main chain.

rgroup_df = filter_dataframe(rgroup_df, "node_id", list(G.nodes()), True)
HBOND_ATOMS = [
"ND1", # histidine
"ND2", # asparagine
"NE", # arginine
"NE1", # tryptophan
"NE2", # glutamate, histidine
"NH1", # arginine
"NH2", # arginine
"NZ", # lysine
"OD1", # asparagine, aspartic
"OD2", # aspartic
"OE1", # glutamate
"OE2", # glutamate
"OG", # serine
"OG1", # threonine
"OH", # tyrosine, threonine, serine
"SD", # methionine
"SG", # cysteine
"N", # backbone
"O",
]

hbond_df = filter_dataframe(rgroup_df, "atom_name", HBOND_ATOMS, True)

if len(hbond_df.index) > 0:
distmat = compute_distmat(hbond_df)
interacting_atoms = get_interacting_atoms(3.5, distmat)
Expand All @@ -404,20 +405,19 @@ def add_hydrogen_bond_interactions(
for r1, r2 in get_edges_by_bond_type(G, "hbond"):
# print(G.nodes.data())
condition1 = (
G.nodes[r1]["atom_type"] in ACCEPTORS_ATOM
and G.nodes[r2]["atom_type"] in DONORS_ATOM
G.nodes[r1]["atom_type"] in ACCEPTORS_ATOM
and G.nodes[r2]["atom_type"] in DONORS_ATOM
)
condition2 = (
G.nodes[r2]["atom_type"] in ACCEPTORS_ATOM
and G.nodes[r1]["atom_type"] in DONORS_ATOM
G.nodes[r2]["atom_type"] in ACCEPTORS_ATOM
and G.nodes[r1]["atom_type"] in DONORS_ATOM
)
is_ionic = condition1 or condition2
if not is_ionic:
G.edges[r1, r2]["kind"].remove("hbond")
if len(G.edges[r1, r2]["kind"]) == 0:
G.remove_edge(r1, r2)


def add_ionic_interactions(
G: nx.Graph, rgroup_df: Optional[pd.DataFrame] = None
):
Expand Down