Conformers and Docking with Vina

[dayang-us]

Hi all!

This issue can be viewed in isolation / together with the other Issues I've created so far: Issue #30 : Characterisation of the 2D protein-ligand interaction map Issue #31 : Making a molecule that meets the eNTRy criteria

For my project I will only be working on a small group of selected molecules and not from a library. Over the last week I followed the AutoDock Vina video tutorial step-by-step and managed to dock each molecule separately.

However, I encountered another set of roadblocks:

Ligand Minimisation

I imported the structures from CHEMBL and copied them into ChemSketch 3D to generate mol files, which was converted into pdb via PyMOL. I used AutoDockTools to prepare the ligands and enzyme and exported them as PDBQT file for docking with Vina.

1. I could not generate low-energy conformers of the molecules, so the protein-ligand complex models I have generated from Vina may not be representative of the actual binding modes.

   • I have tried using OpenBabel to combine the molecules (either using SMILES from DataWarrior, or mol files from ChemSketch) into a single SDF file and I've followed the Jupyter notebook for conformation and docking but either got stuck at checking for the molecules' names (Can't read the molecule) or at the generating conformers step (broken process pool)
   • I have also tried using DataWarrior to Generate Conformers but it yielded an empty file. I want to try troubleshooting this step later this week or
   • I might also try generating conformers using UCSF Chimera

There's plenty of points in the steps where the error might occur. I suspect that the coordinates from the SMILES/mol file was not read correctly. I would really appreciate it if anyone can point me to websites/documents/published protocols for converting sketched molecules to the correct mol/sdf/pdb/pdbqt files and troubleshooting the errors in the Jupyter notebook / DataWarrior, thank you! (I was signposted to Dr Chris Swain @drc007, hope you can help? I am using Windows though...)

Macromolecule Preparation

2. AutoDockVina does not automatically assign the charges for the manganese ions in the crystal structure so the metal ion interactions was not correctly computed and identified during docking. For the next iterations I can add the charges manually into the pdbqt file. But...

   • I can see from previous literature that for docking studies the solvent will be removed, because the docking programs focuses on the interactions between the ligand and the protein side chains only
   • So in this case, should I remove/ignored the charges of the manganese ions from the crystal structure too? They were involved in the interaction with the ATP analogue ANP in the crystal structure, but so is water?

3. For docking, I did not select any side chains as flexible. There is an option to do this in AutoDock Vina but is it compulsory / offer better results than modelling a flexible ligand to a rigid receptor?

4. I set the search space centered around the ATP binding site. Will setting a 'blind search' (opening up the search space to the whole molecule) help with making better correlations/conclusions about the best predicted binding modes and their actual binding affinity and/or potency?

If anyone is interested to help out, I am happy to discuss further.

Thank you!

[jarvist]

For using openbabel to generate conformers, I hit some strange bug with automatically (from SMILES) generated SDF files, where it would hang on the presence of the title. obabel x.sdf -O x_conformers.sdf --title "" --conformer --nconf 100 --score rmsd --writeconformers Is taken from a working script I was using. Overwriting the title with nothing fixes this problem. (Here these are selected for RMSD diversity, rather than energy, to try and cover a maximum phase space.)

Here is a working script, using obabel and an old C SMILES->SDF converter program that I find most reliable (https://depth-first.com/articles/2007/12/12/simple-3d-conformer-generation-with-smi23d/), to do high throughput docking with vina: https://github.com/QuantumCorona/AutodockVina_3CL-pro/blob/d9677ed849ebfdd78b589325c1a0401c7449ed6f/0014-MAK-UNK/dock_from_SMILES_local.sh

However I don't really understand why you are trying to generate different conformers. Autodock Vina has it's own heuristics for identifying & sampling rotational degrees of freedom of the ligand. Unless you are sampling some degree of freedom that it doesn't touch, I don't think pre-generating conformers will make any difference to the docking. My understanding of docking in general is that if you increase the exhaustiveness sufficiently, you are likely to identify the correct pose, but that the empirical pseudo-energy based ranking will not suggest it as the best.

I am not sure if Vina is that suitable for metalloenyzmes. http://autodock.scripps.edu/resources/autodockzn-forcefield - might be a better starting point.

[bendndi]

@dayang-us Not suire if this is your question, but If all you have is a datawarrior file with the list of structures, you can easily export this into an sd file:

File -> Save special -> sd file

If however what you're looking to do is take a list of individual smiles strings and create a single sd file then there is a quick / expensive option or a slighlty longer but free option.

Expensive option requires having chemdraw for excel or biovia insight for excel. I'm going to assume you don't have either.

The free option is to use Knime. Install Knime and select the opton to install all available add in / comunity nodes as default (a "full install". Create an excel file with Smiles column and an identifier column, then run this through a knime protocol where you can use e.g. openbabel node to transform entire excel smiles liet (or csv) into an sdf. It's a very straighforward 3 node protocol:

CSV Reader ----Open Babel (select smiles as the input and sdf as the output) --- sdf writer

[holeung]

Yes, Autodock Vina generates its own conformers. You just need to feed it your molecules in a reasonable 3D conformation. This can be done in openbabel with the --gen3d flag.

If your inhibitors are suspected to bind to magnesium, you will have to find special software that can handle that. I believe Autodock (not Autodock Vina) can.

[MelchorSanchez]

1 - To generate conformers there are different alternatives. One of them is to use rdkit. You can find different scripts out there like https://github.com/dkoes/rdkit-scripts/blob/master/rdconf.py or https://gist.github.com/tdudgeon/b061dc67f9d879905b50118408c30aac. Some more info https://pubs-acs-org/doi/10.1021/acs.jcim.5b00654 and https://pubs.acs.org/doi/10.1021/acs.jcim.0c00025

Autodock vina generate its own conformers but are dependent on the initial conformation. You need to have a reasonably good initial 3D conformation. You can generate conformers and select the best one based on internal energy and/or experimental torsional angles matching. Another possibilitty is to generate a 3D structure with openbael or Rdkit and then optimize it.

2 - For prepare the protein you can use pdb2pqr, that includes the usage of propka (but metal ions can be a problem and as you said maybe you have to add the charges manually to the pqr file in the end). Other option is to use some tool like PDBFixer. There exists AutoDockZN. Maybe it can help you, but I believe it should be used together with AutoDock4, not with Vina.

3 - With AutoDock you can perform flexible docking, but I believe not with vina. Another option that is nowadays rather popular is to run small MD simulations with the best docking poses to observe the induced fit events (adaptation form ligand to target and viceversa)

4 - You need to have a big enough docking box. 4-6 Armstrongs from the geometric center of the binding pocket would be enought. However this is the case if you are sure of the binding site you want to explore. If not, you can perform a blind docking. Or just a cavity detection analysis and then a usual docking over the "best" cavities.

[dayang-us]

@jarvist @holeung Thank you so much for the insights! At the beginning I had the same understanding too - that Vina should generate their on conformations and rotate the molecule around using the torsion data provided when preparing the ligand in ADT. I became unsure when I saw that docking without prior ligand energy minimisation might generate false negatives (can't seem to find the source now, sorry!)

I suppose as @MelchorSanchez has said I can generate low energy conformers and use that as a starting point for docking with Vina.

Thanks for helping me out with the rest of the questions too - it will be really helpful for planning the next iteration of my project, esp. using AutoDockZn for metalloenzymes as all of you have mentioned.

And thank you @bendndi for the guide on generating sdf files using DataWarrior and Knime!

[drc007]

Those interested in Open Source tools might be interested in a series of free workshops organised by RSC CICAG. 9 -13 Nov 2020. https://www.rsc.org/events/detail/43180/workshop-on-open-source-tools-for-chemistry.

We would also be interested in hearing about other software tools we could include in future workshops.

[dayang-us]

Those interested in Open Source tools might be interested in a series of free workshops organised by RSC CICAG. 9 -13 Nov 2020. https://www.rsc.org/events/detail/43180/workshop-on-open-source-tools-for-chemistry.

We would also be interested in hearing about other software tools we could include in future workshops.

Thank you @drc007! Looking forward to the workshops.

[jarvist]

At the beginning I had the same understanding too - that Vina should generate their on conformations and rotate the molecule around using the torsion data provided when preparing the ligand in ADT. I became unsure when I saw that docking without prior ligand energy minimisation might generate false negatives (can't seem to find the source now, sorry!)

This is definitely true if you feed in a distorted ligand. But as part of the conversion from SMILES -> SDF its standard to use an empirical forcefield (usually UFF or MM94) to optimise the geometry of the ligand. You can use Babel's 'obminimize' to do this for an individual structure, to see if anything notably changes. (This convenience program is a bit rough around the edges... I remember it is hardwired to spit out a single structure to STDOUT.) Avagadro is a GUI program which has force-fields and an optimiser built in, if you just want to look at a couple of individual examples.

[MelchorSanchez]

1 - To generate conformers there are different alternatives. One of them is to use rdkit. You can find different scripts out there like https://github.com/dkoes/rdkit-scripts/blob/master/rdconf.py or https://gist.github.com/tdudgeon/b061dc67f9d879905b50118408c30aac. Some more info https://pubs-acs-org/doi/10.1021/acs.jcim.5b00654 and https://pubs.acs.org/doi/10.1021/acs.jcim.0c00025

Autodock vina generate its own conformers but are dependent on the initial conformation. You need to have a reasonably good initial 3D conformation. You can generate conformers and select the best one based on internal energy and/or experimental torsional angles matching. Another possibilitty is to generate a 3D structure with openbael or Rdkit and then optimize it.

2 - For prepare the protein you can use pdb2pqr, that includes the usage of propka (but metal ions can be a problem and as you said maybe you have to add the charges manually to the pqr file in the end). Other option is to use some tool like PDBFixer. There exists AutoDockZN. Maybe it can help you, but I believe it should be used together with AutoDock4, not with Vina.

3 - With AutoDock you can perform flexible docking, but I believe not with vina. Another option that is nowadays rather popular is to run small MD simulations with the best docking poses to observe the induced fit events (adaptation form ligand to target and viceversa)

4 - You need to have a big enough docking box. 4-6 Armstrongs from the geometric center of the binding pocket would be enought. However this is the case if you are sure of the binding site you want to explore. If not, you can perform a blind docking. Or just a cavity detection analysis and then a usual docking over the "best" cavities.

I have forget to mention before that there is this software. It seems that can work properly with metalloproteins and it is open source. I have never test it but maybe you can take a look. https://onlinelibrary.wiley.com/doi/full/10.1002/jcc.24847 , https://github.com/insilichem/gaudi

[dayang-us]

At the beginning I had the same understanding too - that Vina should generate their on conformations and rotate the molecule around using the torsion data provided when preparing the ligand in ADT. I became unsure when I saw that docking without prior ligand energy minimisation might generate false negatives (can't seem to find the source now, sorry!)

This is definitely true if you feed in a distorted ligand. But as part of the conversion from SMILES -> SDF its standard to use an empirical forcefield (usually UFF or MM94) to optimise the geometry of the ligand. You can use Babel's 'obminimize' to do this for an individual structure, to see if anything notably changes. (This convenience program is a bit rough around the edges... I remember it is hardwired to spit out a single structure to STDOUT.) Avagadro is a GUI program which has force-fields and an optimiser built in, if you just want to look at a couple of individual examples.

Awesome, it's all clear now, will try Babel and/or Avogadro for setting up the ligand. Thanks @jarvist!

[drc007]

There are graphical user interfaces for Openbabel available if you would prefer not to use the command-line.

Windows http://openbabel.org/wiki/Windows_GUI

MacOSX https://www.macinchem.org/ibabel/version4/iBabel4_0.php

[dayang-us]

@MelchorSanchez Great, thanks for the software recommendation!

@drc007 Thank you! I'm using OpenBabel with GUI at the moment