Exception: Particle coordinate is nan in addMembrane

[nitroamos]

Hi,

I'm getting this exception almost every time I try to add a membrane (I think it only worked once out of the dozens of attempts on a couple different proteins like 1pga (this is the one it worked on), 2x2v, 4yay):

  File "/hpc/grid/adw_hpcc_workspace/users/andera52/PDBCompiler/pdbCompilerAnacondaEnv/lib/python3.6/site-packages/simtk/openmm/app/modeller.py", line 1366, in addMembrane
    integrator.step(20)
  File "/hpc/grid/adw_hpcc_workspace/users/andera52/PDBCompiler/pdbCompilerAnacondaEnv/lib/python3.6/site-packages/simtk/openmm/openmm.py", line 4436, in step
    return _openmm.LangevinIntegrator_step(self, steps)
Exception: Particle coordinate is nan

I'm not sure what I'm doing wrong. I'm doing default options and POPE. Some possibilities:

1. I've been trying the tool memembed and I'm seeing that it doesn't produce exactly the same alignment as what I'm seeing from a file I downloaded from OPM (e.g., 4yay). The DUM residues are coplanar and share the same z-axis centering, but they aren't centered on top of each other in the xy plane (or necessarily at the origin). Does xy plane placement matter? I've tried both, and both give the same error.
2. Does PDBFixer care about the DUM residues or should I remove them before doing addMembrane?
3. One time I got a different exception which seems to imply it's running on a GPU. Should I be doing this on a CPU for robustness? If so, is there a way to force a CPU from PDBFixer? I'm writing my own python script and not using the UI so I could reach in and tweak the fixer object. I do see this error on a CPU, too, so maybe the type of processor doesn't matter.
4. Do you have a specific pdb I could try by way of a positive control?

I have time to help debug this, but I've not used PDBFixer or OpenMM before so I'm not sure where to start digging.

Thanks!

OpenCL internal error: CL_MEM_OBJECT_ALLOCATION_FAILURE error executing CL_COMMAND_NDRANGE_KERNEL on Tesla M60 (Device 0).

  File "/hpc/grid/adw_hpcc_workspace/users/andera52/PDBCompiler/pdbCompilerAnacondaEnv/lib/python3.6/site-packages/simtk/openmm/app/modeller.py", line 1358, in addMembrane
    LocalEnergyMinimizer.minimize(context, 10.0, 30)
  File "/hpc/grid/adw_hpcc_workspace/users/andera52/PDBCompiler/pdbCompilerAnacondaEnv/lib/python3.6/site-packages/simtk/openmm/openmm.py", line 2631, in minimize
    return _openmm.LocalEnergyMinimizer_minimize(context, tolerance, maxIterations)
Exception: Error invoking kernel findBlocksWithInteractions: clEnqueueNDRangeKernel (-4)

Error I got when running on my laptop using the pdbfixer webapp:

  File "/Users/andera52/anaconda/lib/python3.6/site-packages/pdbfixer-1.5-py3.6.egg/pdbfixer/ui.py", line 137, in addHydrogensPageCallback
    fixer.addMembrane(lipidType, 0*unit.nanometer, padding, positiveIon, negativeIon, ionicStrength)
  File "/Users/andera52/anaconda/lib/python3.6/site-packages/pdbfixer-1.5-py3.6.egg/pdbfixer/pdbfixer.py", line 1089, in addMembrane
    modeller.addMembrane(forcefield, lipidType=lipidType, minimumPadding=minimumPadding, positiveIon=positiveIon, negativeIon=negativeIon, ionicStrength=ionicStrength)
  File "/Users/andera52/anaconda/lib/python3.6/site-packages/simtk/openmm/app/modeller.py", line 1366, in addMembrane
    integrator.step(20)
  File "/Users/andera52/anaconda/lib/python3.6/site-packages/simtk/openmm/openmm.py", line 14924, in step
    return _openmm.LangevinIntegrator_step(self, steps)

[nitroamos]

I did just get an OPM aligned 4yay and a memembed aligned 3o7qA to work with POPC. Maybe POPC is more likely to work than POPE?

[peastman]

Could you test it out with the latest OpenMM development code? You can install it with conda install -c omnia/label/dev openmm. We recently fixed a bug in building membranes. This might not be related, but it would be good to rule it out.

If that doesn't fix it, could you post the script you're using? I'll see if I can reproduce it.

Does PDBFixer care about the DUM residues or should I remove them before doing addMembrane?

You should remove them. PDBFixer doesn't use them. All that matters is that the protein is properly oriented and positioned along the z axis, which will be the case for an OPM structure.

[nitroamos]

I've switched to the dev version of openmm (7.3.2) and I'm seeing the same issue. Here's a script:

fixer = PDBFixer(filename="4yay.pdb")
fixer.findMissingResidues()
fixer.findNonstandardResidues()
fixer.replaceNonstandardResidues()
fixer.findMissingAtoms()
fixer.addMissingAtoms()
fixer.addMissingHydrogens(7.0)

print('\nAdding membrane:',args.membrane_lipid_type)
PDBFile.writeFile(fixer.topology, fixer.positions, open("temp.pdb",'w'))
fixer.addMembrane(lipidType=args.membrane_lipid_type,
                  membraneCenterZ=0*unit.nanometer,
                  minimumPadding=1*unit.nanometer,
                  positiveIon="Na+",
                  negativeIon="Cl-",
                  ionicStrength=0.0*unit.molar)

And I got that exception all but one time I tried (after 2-9 minutes). I've attached my input file, which I downloaded the file directly from OPM (chosen for its optimistic name!) and all I did was remove the DUM residues.

I also attached the file I saved right before adding the membrane. The TYR at position 283 (312 before renumbering) looks pretty clashy, don't know if that's related.

It did work one time, and I've attached the output. It has an issue, described in #164 relating to overflowing the resnum field width.

temp.pdb.gz
4yay.pdb.gz
output.pdb.gz

[peastman]

The TYR at position 283 (312 before renumbering) looks pretty clashy, don't know if that's related.

I think you're correct. The whole sidechain is missing in the original file, so it has to get added. But it's jammed right up against a couple of other residues, so it has to get built in just the right conformation or it will clash. PDBFixer is having trouble finding that conformation. But it sounds like it does sometimes manage to find it, just not always. So you can just run it a few times until it succeeds.

[nitroamos]

I've been working away at this. I've integrated Rosetta relax into my project and can use it to refine the structure before trying to add the membrane. I've attached the refined starting point I'm using for adding a membrane.

I've been studying the code in modeller.addMembrane. I think it makes sense to me... so I wrote some debugging code to print out all pairwise distances so that I could see if anything is clashing. Here are some of the distances that I'm seeing right before the LocalEnergyMinimizer.minimize line. There are 123938 atoms and there are 2534 pairs of atoms within 0.09 nm (yes, this uses a lot of memory...). The closest pairs are:

Printing close distances: 0.09 nm 123938
Num close: 2534
 34900 HOH    O :  34978 HOH   H2 : 0.005442828637598851 (-1.7426923461914061, -2.721847720336916, 2.2090001106262207) (-1.746558541870117, -2.7256786300659197, 2.208967781066894)
 33854 HOH   H1 :  35204 HOH   H1 : 0.020743199081561612 (-2.835491532897949, -0.9788476898193377, -2.555000686645508) (-2.824309701538086, -0.9621478988647478, -2.549866104125977)
 3783  HOH   H1 :  34993 HOH    O : 0.04253896908145396 (-1.433692330932617, -2.407848544311525, 2.069000005722046) (-1.4659341522216796, -2.434800334167482, 2.0623947143554684)
 31873 HOH    O :  33561 HOH    O : 0.04502261585893706 (-2.4866927810668944, -2.3648483230590838, 3.500999927520752) (-2.4783385940551756, -2.3206541015625017, 3.498971557617187)
 33854 HOH   H2 :  35201 HOH   H1 : 0.04508281432372642 (-2.742492074584961, -0.9748479797363299, -2.5750001907348636) (-2.7385457702636717, -1.0176488830566424, -2.588600540161133)
 3593  HOH   H1 :  34275 HOH   H2 : 0.04743480970935302 (-3.496492738342285, -2.1958486511230486, -2.731000328063965) (-3.5289453216552733, -2.175949282836916, -2.7027000427246097)
 34685 HOH   H2 :  34253 HOH   H1 : 0.04750599168750323 (-2.7866929718017577, -2.2438480331420916, -2.885000610351563) (-2.813459748840332, -2.2808820678710955, -2.897994422912598)
 3593  HOH   H2 :  34279 HOH    O : 0.047686540224351044 (-3.564492578125, -2.2488481475830095, -2.7730001449584964) (-3.540090913391113, -2.285087771606447, -2.753888511657715)

That is, those distances are, I think, directly from the initial placement. I don't know if distances this close would cause atoms to shoot off to infinity, but before I try to dig much deeper this seems like a good time to check in with you with some questions.

1. Based on what I'm seeing, I'm only seeing waters clash with each other and since starting structure has no water, these must all be waters added by addMembrane. It doesn't look like newly placed waters are being added to your distance overlap check for subsequent waters. Is that intentional? You're not doing that for the POP residues either, but maybe those are already not mutually clashing. When I'm testing it, it's often failing the first integration step where, as I understand, the protein is scaled the smallest, meaning clashes are most likely not associated with the protein.

2. Is it using a Lennard-Jones potential such that we're seeing an inner wall clash? If so, do you have a softer potential we could try here? e.g., something with a capped plateau? This would be a better solution if you knew you wanted to keep all the waters.

3. Another strategy that seems reasonable is to spend more time in the LocalEnergyMinimizer before attempting MD and e.g., setting maxIterations to 0. However, its documentation describes "distance constraints" but doesn't say exactly what is constrained to where which makes me wonder if it's being given sufficient freedom to resolve the clashes. I tried minimizing for 1000 steps (instead of 30) and it didn't clear the clashes.

4. The other solution I've tried is to repeatedly retry adding the membrane. However, it fails too frequently, at least in the 4yay system I'm trying, to be a viable solution. It did work once, so I attached the pdb.

Any suggestions?

Thanks!

4yay_start.pdb.gz
output.cif.gz

[nataliaszulc]

It is solved already ? I've the same problem :-(