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

Charges dependent on the conformation #40

Open
maciejwisniewski-drugdiscovery opened this issue Jul 1, 2024 · 3 comments
Open

Charges dependent on the conformation #40

maciejwisniewski-drugdiscovery opened this issue Jul 1, 2024 · 3 comments

Comments

@maciejwisniewski-drugdiscovery

Hi,
Are the Espaloma charges dependent on the conformation of the chemical compound?
Because after using this code, I get the same charges for each conformation, and I thought that the AM1 method inherently takes conformation into account, unlike, for example, Gasteiger.

Code:

import rdkit
from rdkit import Chem
from rdkit.Chem import AllChem
from rdkit.Chem import rdDistGeom, rdFMCS
from openff.toolkit import Molecule
from espaloma_charge.openff_wrapper import EspalomaChargeToolkitWrapper
import numpy as np
import time
from espaloma_charge import charge

toolkit_registry = EspalomaChargeToolkitWrapper()
def get_conformations(molecule):
    all_partial_charges = []
    molecule.RemoveAllConformers()
    molecule = Chem.AddHs(molecule)
    params = rdDistGeom.ETKDGv3()
    params.randomSeed = 0xd06f00d
    params.numThreads = 4
    params.maxAttempts = 100
    # print how many atoms there is
    cids = rdDistGeom.EmbedMultipleConfs(molecule, 32, params)
    if len(cids) == 0:
        print("Failed to generate conformations for molecule. Skipping.")
        return None, None
    molecule = Chem.RemoveHs(molecule)
    # Calculating Partial Charges with AM1-BCC
    conformers = Molecule.from_rdkit(molecule,
                                     allow_undefined_stereo=True,
                                     hydrogens_are_explicit=False).conformers
    molecule.RemoveAllConformers()
    for conformer in conformers:
        start = time.time()
        template = Molecule.from_rdkit(molecule,
                                       allow_undefined_stereo=True,
                                       hydrogens_are_explicit=False)
        template.add_conformer(conformer)
        
        template.assign_partial_charges(
            partial_charge_method = 'espaloma-am1bcc',
            toolkit_registry=toolkit_registry,
            use_conformers = template.conformers
        )
        end = time.time()
        print('DURATION:\t',(end - start)/60)
        all_partial_charges.append(template.partial_charges)
        
    return all_partial_charges
rdkit_mol = Chem.MolFromSmiles('COC(=O)C(C1CCCCN1)C2=CC=CC=C2')
charges = get_conformations(rdkit_mol)

@ljmartin
Copy link

I wondered this too - looks like the graph definition does not include conformer information:

def from_rdkit_mol(mol, use_fp=True):

@jgreener64
Copy link

Espaloma charges do not depend on conformation, just topology, though Espaloma was trained on conformation-dependent charges (I am not involved in the project).

@pavankum
Copy link

The paper mentions training to AM1-BCC ELF10 charges, which is a method to diminish conformer dependence, https://pubs.acs.org/doi/full/10.1021/acs.jpca.4c01287.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants