I will be docking a small library of molecules that I have downloaded from PubChem using AutoDock. I thought why not minimize their energy before docking using ChemDraw. I am thinking; what is the best method of optimization/ energy minimization of small molecules. MM2, MMF94 or other? And why?
1 Answer
The difference between different force-fields is not going to be major, it is the side steps which are.
When
If you are starting from a SMILES string, optimisation is a must obviously.
If you are using a 3D conformer from PubChem or even an actual sub-1 Å crystal structure from CSD, optimisation is nice for consistency.
Which
MMF94 is a solid choice. RDkit offers MMF94 (Merck force field 1994) and UFF (universal force field). The latter has some extra features, such as not kicking the bucket at atoms that aren't your usual suspects. ChemDraw offers MM2 and MMF94 IIRC, the former is very old. g16 in Gaussian (paid software) uses Gaussian Electrostatic Model, which is polarisable and will solve partial charges in a nice way (which I'll get to in a minute). OpenBabel uses MMF94 under the hood. You can even optimise with macromolecular forcefields, such as AMBER and CHARMM.
Given that the AutoDock uses a flexible ligand it should not matter too much anyway.
What else?
Autodock does not require conformers by virtue of the flexible ligand business (the sampler does rototorsions). Other programs might if they use rigid ligands (e.g. Fred) or they use a mixed approach and use a conformer library for speed (e.g. Rosetta ligand_dock). Conformer diversity is actually a bigger issue as conformers that are only slightly contorted, often even less than kB (0.6 kcal/mol), may not be sampled.
The thing that messes things is incorrect protonation. OpenBabel can protonate at a given pH, RDKit cannot. Not an issue for ligands with drug screens (cf. Lipinski rule), but is for natural substrates. However, for the protein it is essential to make sure acidic and basic residue are correctly protonated and that histidine have the correct tautomer —there are tools, c.f. MOE, that can do it but by far the best approach is reading up on the mechanism.
The next thing is dodgy Gasteiger Marsili partial charges —for Autodock 4 (and most other full/hybrid physics-based docking programs), partial charges are important: after all PDBQT has a Q in it for charge! For prepare_ligand4.py
you give it a mol2 file. A major consideration is what format you are using for your molecules. This sounds a very silly thing to point out, but you'd be surprised how many trip up. An mol file (or an sdf file, which is the same but multi-entry) does not have an official partial charge column (although they are often added as comments after M END
). So whatever force-field you chose make sure that the program you use gave you partial charges are not formal charges. This does not quite apply to AutoDock Vina, which is charge independent and more reliant on machine learning derived factors.
Verdict
Open babel is nice in that it can start with a SMILES, handles multiple formats including mol2 and PDBQT, can protonate at pH 7, uses MMF94 and does Gasteiger Marsili charges (and its confab algorithm is very good, not need here though).