tdudgeon / simple-simulate-complex

Simple protein-ligand complex simulation with OpenMM
73 stars 16 forks source link

Is the prepareComplex.py code in the "Protein and ligand preparation" section correct? #7

Open Shredderroy opened 2 months ago

Shredderroy commented 2 months ago

Firstly, thank you very much for attempting to simplify the workflow for those who are not computational biologists.

This is the specific excerpt that the title is referring to:

python prepareComplex.py data/Mpro-x0387_0_apo-desolv.pdb Mpro-x0387_0

But the source file prepareComplex.py indicates that one must supply a receptor name, a ligand name, and a complex name. That has not been done here.

Even so, which receptor should be given as the argument here? Should it be one of the outputs of the previous prepareProtein.py step? And if so, should it be the minimised file or the fixed file?

It would be incredibly helpful if you could confirm whether the prepareReceptor.py command should be something like:

$ python prepareComplex.py Protein_minimised.pdb ligand.mol PLComplex

Thanks in advance for your clarification.

Shredderroy commented 2 months ago

I want to add that I have managed to make some progress and actually get a full simulation to run. Here are the steps that I took, just in case they are helpful to someone who is stuck as I was:

  1. I ran prepareProtein.py as shown in the README file.
  2. I used the output PDB file of the docking pose of a ligand, which I converted to MOL using Wolfram Mathematica (but one may also use OpenBabel).
  3. I ran simulateComplex.py with the protein_minimised and ligand.mol files as input.

With the above set-up, I was able to consistently run the simulations for different proteins and ligands.

But when I run simulateComplex.py with the --solvate option, I find that my complex is always off to one corner of the box, with a portion of the ribbon actually sticking out of the cube, like in the image below:

OpenMM-test

Am I doing something incorrectly? Is there a way to force the complex to be in the middle of the solvent cube?

Thanks again for your help/

Shredderroy commented 2 months ago

Just to add a quick note: I tried altering the call to modeller.addSolvent in line 94 of simulateComplex.py by replacing the padding option with boxSize as follows:

modeller.addSolvent(
    system_generator.forcefield,
    model = args.water_model,
    boxSize = Vec3(5.0, 5.0, 5.0) * unit.nanometer,
    positiveIon = args.positive_ion,
    negativeIon = args.negative_ion,
    ionicStrength = args.ionic_strength * unit.molar,
    neutralize = not args.no_neutralize
)

The change was taken from the official OpenMM documentation.

In boxSize I tried varying only one of the axes while keeping the other two at 5.0. But none of the axes seemed to cover that little bit of ribbon sticking out of the box. I could see the visual difference resulting from an increase in the x and z sizes. But an increase in the y size only lengthens the box away from the ribbon.

So, I am not really sure how I can get the solvent box to cover that protruding ribbon.

tdudgeon commented 1 month ago

This looks like a periodic box issue, which is quite common and a bit confusing. Are you using the analyse.py script to re-image the trajectory? This sometimes centres things.

Shredderroy commented 1 month ago

No, I have not been using analyse.py. I could not get that file to run without errors.

The image above was taken from PyMol, which opened both the output_minimised and output_traj files and played a wiggling movie.

In spite of the skewed visualisation, can I at least assume that the simulation was done correctly, and the distribution of water molecules during the calculation was more or less isotropic?

Shredderroy commented 1 month ago

I just want to add that I have managed to get analyse.py and that script appears to have corrected the off-centre image:

OpenMM-Complex-Solvated