Alternative to CGenFF for generating large ligand topology

Question

I am now preparing a ligand topology following GROMACS protein-ligand tutorial. However, I have a large ligand (446 atoms).

$ wc -l my_ligand.pdb
446 my_ligand.pdb

When I tried to upload my_ligand_fix.mol2 to CGenFF server, I noticed it the maximum number of atoms allowable is 384.

Note that my_ligand_fix.mol2 has gone through processing right up to sort_mol2_bonds.pl output.

Is there a way to get around this? Or what other alternative to CGenFF that I can use instead?

-Sca

Answer 1

I'd use OpenMM's openmmforcefields package, which is available on GitHub.

I haven't used it extensively, but I recommended to a student in a colleague's lab and they've found it useful.

Small molecule force fields:

• GAFF 1.x and 2.x parameters
• OpenFF parameters

Biomolecule force fields:

• Amber-type (albeit not Amber ff19SB yet)
• CHARMM-type

And it will, if needed, generate AM1-BCC partial charges using OpenEye or antechamber tools.

The README gives some example Python code:

# Create an OpenFF Molecule object for benzene from SMILES
from openff.toolkit.topology import Molecule
molecule = Molecule.from_smiles('c1ccccc1')
# Create the GAFF template generator
from openmmforcefields.generators import GAFFTemplateGenerator
gaff = GAFFTemplateGenerator(molecules=molecule)
# Create an OpenMM ForceField object with AMBER ff14SB and TIP3P with compatible ions
from simtk.openmm.app import ForceField
forcefield = ForceField('amber/protein.ff14SB.xml', 'amber/tip3p_standard.xml', 'amber/tip3p_HFE_multivalent.xml')
# Register the GAFF template generator
forcefield.registerTemplateGenerator(gaff.generator)
# You can now parameterize an OpenMM Topology object that contains the specified molecule.
# forcefield will load the appropriate GAFF parameters when needed, and antechamber
# will be used to generate small molecule parameters on the fly.
from simtk.openmm.app import PDBFile
pdbfile = PDBFile('t4-lysozyme-L99A-with-benzene.pdb')
system = forcefield.createSystem(pdbfile.topology)

For example you could load your molecules like this from an .sdf file (again from the README):

molecules = Molecule.from_file('molecules.sdf')
gaff = GAFFTemplateGenerator(molecules=molecules, forcefield='gaff-2.11')