Skip to content

New dih #607

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.

Sign up for GitHub

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

Merged
carlocamilloni merged 3 commits into from
Jan 29, 2026
Merged

New dih #607

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
950612f
lys
carlocamilloni Dec 11, 2025
11c7d3a
Merge branch 'multi-ego:new_dih' into new_dih
carlocamilloni Dec 11, 2025
ed864d5
preliminary changes for the new model
carlocamilloni Jan 29, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion multi-ego-basic.ff/ffbonded.itp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -536,7 +536,7 @@
#define gd_44K -120.000 1.00 1
; Backbone dihedral angle -C-N-CA-C- LYS
;
#define gd_45K 180.000 2.00 1
#define gd_45K 180.000 1.40 1
; Backbone dihedral angle -N-CA-C-N- LYS
;
;;;;;;;;;;;;;;;;;;;;
Expand Down
35 changes: 31 additions & 4 deletions src/multiego/ensemble.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -349,6 +349,10 @@ def init_meGO_matrices(ensemble, args, custom_dict):
reference_contact_matrices[name] = io.read_molecular_contacts(
path, ensemble["molecules_idx_sbtype_dictionary"], reference["reference"], path.endswith(".h5")
)

# TODO: probabilities of equivalent atoms should be averaged


reference_contact_matrices[name] = reference_contact_matrices[name].add_prefix("rc_")
reference_contact_matrices[name]["c6_i"] = [lj_data_dict[x][0] for x in reference_contact_matrices[name]["rc_ai"]]
reference_contact_matrices[name]["c6_j"] = [lj_data_dict[x][0] for x in reference_contact_matrices[name]["rc_aj"]]
Expand Down Expand Up @@ -1035,14 +1039,30 @@ def set_sig_epsilon(meGO_LJ, parameters):
# by comparing them with RC_probabilities so that the resulting epsilon is between eps_min and eps_0
condition = (
meGO_LJ["probability"] > meGO_LJ["limit_rc_att"] * np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"])
) & (meGO_LJ["probability"] > meGO_LJ["md_threshold"])
) & (meGO_LJ["probability"] > meGO_LJ["md_threshold"]) & (meGO_LJ["epsilon_prior"] > 0.)

meGO_LJ.loc[condition, "epsilon"] = np.maximum(0.0, meGO_LJ["epsilon_prior"]) - (
(meGO_LJ["epsilon_0"] - np.maximum(0.0, meGO_LJ["epsilon_prior"])) / np.log(meGO_LJ["rc_threshold"])
) * (np.log(meGO_LJ["probability"] / (np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"]))))
meGO_LJ.loc[condition, "learned"] = 1
meGO_LJ.loc[condition, "sigma"] = meGO_LJ["distance"] / 2.0 ** (1.0 / 6.0)

condition = (
meGO_LJ["probability"] > meGO_LJ["limit_rc_att"] * np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"])
) & (meGO_LJ["probability"] > meGO_LJ["md_threshold"]) & (meGO_LJ["epsilon_prior"] < 0.) & (meGO_LJ["rc_probability"] > meGO_LJ["md_threshold"])
meGO_LJ.loc[condition, "epsilon"] = (-meGO_LJ["rep"] * (meGO_LJ["distance"] / meGO_LJ["rc_distance"]) ** 12).clip(
lower=-20 * meGO_LJ["rep"], upper=-0.05 * meGO_LJ["rep"]
)[condition]
meGO_LJ.loc[condition, "learned"] = 1

#condition = (
# meGO_LJ["probability"] > meGO_LJ["limit_rc_att"] * np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"])
#) & (meGO_LJ["probability"] > meGO_LJ["md_threshold"]) & (meGO_LJ["epsilon_prior"] < 0.) & (meGO_LJ["rc_probability"] <= meGO_LJ["md_threshold"])
#meGO_LJ.loc[condition, "epsilon"] = -meGO_LJ["rep"] * (
# 1.0 + (np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"]) - meGO_LJ["probability"])
#)
#meGO_LJ.loc[condition, "learned"] = 1

# Not-attractive interactions
# this is used only when MD_th < MD_p < limit_rc_att*RC_p
# negative epsilon are used to identify non-attractive interactions
Expand All @@ -1051,14 +1071,21 @@ def set_sig_epsilon(meGO_LJ, parameters):
& (meGO_LJ["probability"] > meGO_LJ["md_threshold"])
& (meGO_LJ["rc_probability"] > meGO_LJ["md_threshold"])
)
# meGO_LJ.loc[condition, "epsilon"] = -meGO_LJ["rep"] * (
# 1.0 + (np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"]) - meGO_LJ["probability"])
# )
meGO_LJ.loc[condition, "epsilon"] = (-meGO_LJ["rep"] * (meGO_LJ["distance"] / meGO_LJ["rc_distance"]) ** 12).clip(
lower=-20 * meGO_LJ["rep"], upper=-0.05 * meGO_LJ["rep"]
)[condition]
meGO_LJ.loc[condition, "learned"] = 1

#condition = (
# (meGO_LJ["probability"] <= meGO_LJ["limit_rc_att"] * np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"]))
# & (meGO_LJ["probability"] > meGO_LJ["md_threshold"])
# & (meGO_LJ["rc_probability"] <= meGO_LJ["md_threshold"])
#)
#meGO_LJ.loc[condition, "epsilon"] = -meGO_LJ["rep"] * (
# 1.0 + (np.maximum(meGO_LJ["rc_probability"], meGO_LJ["rc_threshold"]) - meGO_LJ["probability"])
#)
#meGO_LJ.loc[condition, "learned"] = 1

# 1-4 interactions are special and cannot become attractive because they are part of the bonded-interactions
condition = (meGO_LJ["bond_distance"] < 4) & (meGO_LJ["same_chain"])
meGO_LJ.loc[condition, "epsilon"] = -meGO_LJ["rep"]
Expand Down
66 changes: 49 additions & 17 deletions src/multiego/resources/type_definitions.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,13 +6,15 @@
mg_HH_c12_rep = 1.2e-8
mg_ON_c12_rep = 7.5e-7
mg_NN_c12_rep = 2.5e-5

mg_HO_sigma = 0.169500
mg_eps_ch3 = 0.15
mg_eps_HO = 0.15 # hydrogen bond strength
mg_eps_ch_aromatic = 0.15
mg_eps_ch3 = 0.10
mg_eps_HO = 0.13 # hydrogen bond strength
mg_eps_ch_aromatic = 0.13
mg_eps_ch2 = 0.10
mg_eps_pol = 0.12
mg_eps_ch1 = 0.09
mg_eps_pol_nz = 0.11
mg_eps_pol = 0.09
mg_eps_ch1 = 0.07
mg_eps_bkbn_O_CB = 0.10

# Dataframe with atom types and associated parameters
Expand Down Expand Up @@ -80,7 +82,7 @@
4.0 * 0.35722**12 * mg_eps_pol, # "NT", sig=0.35722
4.0 * 0.31365**12 * mg_eps_pol, # "NL", sig=0.31365
4.0 * 0.33411**12 * mg_eps_pol, # "NR", sig=0.33411
4.0 * 0.31365**12 * mg_eps_pol, # "NZ", sig=0.31365
4.0 * 0.31365**12 * mg_eps_pol_nz, # "NZ", sig=0.31365
4.0 * 0.31365**12 * mg_eps_pol, # "NE", sig=0.31365
4.0 * 0.35812**12 * mg_eps_pol, # "C", sig=0.35812
4.0 * 0.35812**12 * mg_eps_ch_aromatic, # "CH", sig=0.35812
Expand All @@ -106,7 +108,7 @@
4.0 * 0.35722**6 * mg_eps_pol, # "NT",
4.0 * 0.31365**6 * mg_eps_pol, # "NL",
4.0 * 0.33411**6 * mg_eps_pol, # "NR",
4.0 * 0.31365**6 * mg_eps_pol, # "NZ",
4.0 * 0.31365**6 * mg_eps_pol_nz, # "NZ",
4.0 * 0.31365**6 * mg_eps_pol, # "NE",
4.0 * 0.35812**6 * mg_eps_pol, # "C",
4.0 * 0.35812**6 * mg_eps_ch_aromatic, # "CH"
Expand Down Expand Up @@ -198,17 +200,29 @@ def lj14_generator(df):

# Special non-local interactions different from basic mg combination rules
# PROTEIN
polar_sbtype = ["OA", "OM", "N", "NL", "NT", "NR", "NZ", "NE", "C", "S", "P", "OE", "CR1"]
hyd_sbtype = ["CH3", "CH3p", "CH2", "CH2r", "CH1"]
polar_sbtype = ["O", "OA", "OM", "N", "NL", "NT", "NR", "NZ", "NE", "C", "S"]
hyd_sbtype = ["CH3", "CH2", "CH2r", "CH1", "CAH", "CAH2"]
special_non_local = [
#{
# "atomtypes": (["O", "OM"], ["O", "OM"]), # charged oxygen-oxygen repulsion
# "interaction": "rep",
# "sigma": None, # not needed for repulsion
# "epsilon": mg_OO_c12_rep,
#},
{
"atomtypes": (["O", "OM"], ["O", "OM"]), # charged oxygen-oxygen repulsion
"atomtypes": (["O"], ["O"]), # charged oxygen-oxygen repulsion
"interaction": "rep",
"sigma": None, # not needed for repulsion
"epsilon": mg_OO_c12_rep,
},
#{
# "atomtypes": (["NZ", "NL"], ["NZ", "NL"]), # charged nitrogen-nitrogen repulsion
# "interaction": "rep",
# "sigma": None, # not needed for repulsion
# "epsilon": mg_NN_c12_rep,
#},
{
"atomtypes": (["NZ", "NL"], ["NZ", "NL"]), # charged nitrogen-nitrogen repulsion
"atomtypes": (["NZ"], ["NZ"]), # charged nitrogen-nitrogen repulsion
"interaction": "rep",
"sigma": None, # not needed for repulsion
"epsilon": mg_NN_c12_rep,
Expand All @@ -219,18 +233,36 @@ def lj14_generator(df):
"sigma": None, # not needed for repulsion
"epsilon": mg_HH_c12_rep,
},
{
"atomtypes": (polar_sbtype, hyd_sbtype), # polar - hydrophobic repulsion
"interaction": "rep",
"sigma": None, # not needed for repulsion
"epsilon": None,
}, # If None use default rc c12 repulsion
#{
# "atomtypes": (polar_sbtype, hyd_sbtype), # polar - hydrophobic repulsion
# "interaction": "att",
# "sigma": None, # not needed for repulsion
# "epsilon": 0.09,
#}, # If None use default rc c12 repulsion
{
"atomtypes": (["O", "OM", "OA"], ["H"]), # hydrogen bond attraction
"interaction": "att",
"sigma": mg_HO_sigma,
"epsilon": mg_eps_HO,
},
{
"atomtypes": (["NL"], ["N", "NT", "NL", "NR", "NZ", "NE", "C", "CH", "CH1", "CAH", "CH2", "CAH2", "CH3", "CH2r", "S"]), # hydrogen bond attraction
"interaction": "rep",
"sigma": None,
"epsilon": None,
},
{
"atomtypes": (["OM"], ["O", "OM", "C", "CH", "CH1", "CAH", "CH2", "CAH2", "CH3", "CH2r", "S"]), # hydrogen bond attraction
"interaction": "rep",
"sigma": None,
"epsilon": None,
},
#{
# "atomtypes": (["CH1", "CAH"], ["CH1", "CAH"]), # hydrogen bond attraction
# "interaction": "rep",
# "sigma": None,
# "epsilon": None,
#},
]

# List of amino acids and nucleic acids
Expand Down
Loading