-
Notifications
You must be signed in to change notification settings - Fork 24
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
Relative Free Energy of solvation #1120
Comments
Hi @mathilfiker, Could you provide more details about how you did your rhfe calculation? I.e. did you calculate both a vacuum and solvent leg? What were your components for the end states of each calculation? |
Hello @IAlibay ,
Maybe there is something I didn't understand from the tutorial, so any input can be useful :). Best, |
So there's two separate issues I can see here:
Issue 1This will need more information:
Issue 2In order to do a relative hydration free energy calculation, you will need a vacuum transformation. This should be something like: ligands = [ligandA,ligandB]
ligand_network = network_planner(
ligands=ligands,
mappers=[mapper],
scorer=scorer
)
# Create your settings
solvent_settings = RelativeHybridTopologyProtocol.default_settings()
vacuum_settings = RelativeHybridTopologyProtocol.default_settings()
vacuum_settings.forcefield_settings.nonbonded_method = "nocutoff"
# Create the Protocol objects
protocols = {
'solvent': RelativeHybridTopologyProtocol(solvent_settings),
'vacuum': RelativeHybridTopologyProtocol(vacuum_settings)
}
transformations = []
for mapping in ligand_network.edges:
for leg in ['solvent', 'vacuum']:
sysA_dict = {'ligand': mapping.componentA,}
sysB_dict = {'ligand': mapping.componentB,}
if leg == 'solvent':
sysA_dict['solvent'] = solvent
sysB_dict['solvent'] = solvent
sysA = openfe.ChemicalSystem(sysA_dict, name=f"{mapping.componentA.name}_{leg}")
sysB = openfe.ChemicalSystem(sysB_dict, name=f"{mapping.componentB.name}_{leg}")
transformation = openfe.Transformation(
stateA=sysA,
stateB=sysB,
mapping={'ligand': mapping},
protocol=protocols[leg]
name=f"{sysA.name}_{sysB.name}_{leg}"
)
transformations.append(transformation)
network = openfe.AlchemicalNetwork(transformations) By calculating both the vacuum and solvent transformation, you can then take their difference to get the ddG of hydration. |
Hello, I tried adding the vacuum leg to the simulation, but both the transformations return a dG of zero. The way I feed partial charges is the following:
This method proved successful when I computed ABHFEs. |
This is rather odd - could you tell us what your python environment looks like @mathilfiker ? Something like |
This is the list of all installed packages:
|
Hi @mathilfiker , thank you for sending this along!
Please let me know if you have any questions! |
Hi @hannahbaumann, thank you for answering!
Is this .db file supposed to look like this? |
Thanks @mathilfiker ! Yes, that's how the db.json would look like, the charges should be in the ffxml, maybe you could just check it for the first ligand atom?
That's very odd indeed! Just to clarify, this is from the |
Where should I find the ffxml? it is not in the working directory. Yes, the DGs = 0 come from simulation_real_time_analysis.yaml |
Thanks @mathilfiker ! This is indeed very weird! We will try to replicate this issue on our end to see whether there could be a bug somewhere. Apologies for the delays on this! |
@mathilfiker - could you clarify if within |
Hi @mathilfiker just a brief update, it looks like this is a legitimate (albeit obscure) bug in the way we parameterize our systems. You can see a very messy bugfix for it here: #1116 FIxing it properly is going to require some rework of some of our Protocol internals and a lot of validation, so we need to discuss things internally (probably next week) before we move forward with this. If you want to try out completly experimental code, then that bugfix may work. |
Dear developers,
Thank you for your work and your dedication.
I'm trying to calculate the relative free energy of solvation of two molecules. They are the same molecules but with different sets of charges. I tried following your RBFE tutorial but neglected the protein. The transformation is correctly created and the calculations run, but the RFE results to be zero. Do you have any suggestions on how to tackle this problem?
Best wishes,
Mathias
The text was updated successfully, but these errors were encountered: