openbabel / openbabel

Open Babel is a chemical toolbox designed to speak the many languages of chemical data.
http://openbabel.org/
GNU General Public License v2.0
1.06k stars 410 forks source link

SetTorsion results in destroyed molecule #2434

Open duerrsimon opened 2 years ago

duerrsimon commented 2 years ago

Environment Information

Open Babel version: Open Babel 2.4.1 3.1.0 and 3.1.1 Operating system and version: Linux

Expected Behavior

I want to change the chi1 torsion angles in a protein sidechain using SetTorsion Only the chi1 dihedral should change.

Actual Behavior

Openbabel also changes some of the backbone dihedrals and thus completely destroys the structure. The only difference between the input pdb files are a few coordinates.

$ diff rotamer1.pdb rotamer2.pdb
< ATOM   1784  HB2 GLU A 117       2.988   8.056  41.088  1.00  0.00           H  
< ATOM   1785  HB3 GLU A 117       3.539   6.939  39.816  1.00  0.00           H  
< ATOM   1786  CG  GLU A 117       1.452   7.430  39.726  1.00  0.00           C  
< ATOM   1787  HG2 GLU A 117       1.208   6.525  39.170  1.00  0.00           H  
< ATOM   1788  HG3 GLU A 117       0.658   7.642  40.442  1.00  0.00           H  
< ATOM   1789  CD  GLU A 117       1.588   8.594  38.764  1.00  0.00           C  
< ATOM   1790  OE1 GLU A 117       2.605   9.314  38.840  1.00  0.00           O1-
< ATOM   1791  OE2 GLU A 117       0.676   8.785  37.931  1.00  0.00           O1-
---
> ATOM   1784  HB2 GLU A 117       3.023   8.049  41.082  1.00  0.00           H  
> ATOM   1785  HB3 GLU A 117       3.515   6.922  39.795  1.00  0.00           H  
> ATOM   1786  CG  GLU A 117       1.439   7.460  39.761  1.00  0.00           C  
> ATOM   1787  HG2 GLU A 117       0.654   7.633  40.497  1.00  0.00           H  
> ATOM   1788  HG3 GLU A 117       1.591   8.365  39.174  1.00  0.00           H  
> ATOM   1789  CD  GLU A 117       1.029   6.326  38.842  1.00  0.00           C  
> ATOM   1790  OE1 GLU A 117       1.917   5.737  38.193  1.00  0.00           O1-
> ATOM   1791  OE2 GLU A 117      -0.182   6.026  38.773  1.00  0.00           O1-

Steps to Reproduce

I made a reproduction here which also allows to view the structures: https://colab.research.google.com/drive/1JU7UA1uRWsS_WYqvMjQpdf4EJ5V6VaKz?usp=sharing

welcome[bot] commented 2 years ago

Thanks for opening your first issue here! Be sure to follow the issue template!

nbehrnd commented 2 years ago

For a moment, I thought the reformulation

Instead of an iteration on all single bonds (section 9.2) how to constrain confab to act on a specific bond?

might describe what you aim for. Within reason, the larger the protein and increasing number of flexible bonds, the more this approach could save computational work ahead.

Does the reformulation match your intent? (Maybe suitable for a post on https://mattermodeling.stackexchange.com/).

duerrsimon commented 2 years ago

@nbehrnd thanks for your comment. I do not really get how confab is related to the problem here?

The SetTorsion function works in the vast majority of cases as intended and just rotates the specifc dihedral (here chi1). In maybe 2 of 100 cases I will produce faulty structures (with systems that only differ in the coordinates and residues/atoms being equal). I double checked that the atoms we are choosing are correct and I also reversed the order of atoms for SetTorsion but nothing gives me a satisfactory result.

I also tested this in openbabel 2.4.1 and the problem exists there too.

I am not an openbabel expert but upon reading ObRotorList Documentation there is an argument for fixed atoms but I am not sure how to use it. Any pointers on how to use this? Perhaps if the bb_nitrogen atom is fixed the rotation would work correctly by rotating the sidechain and not the backbone atoms.

duerrsimon commented 2 years ago

I also checked what happens in other programs. I I do the same change in e.g PyMol the dihedral is correctly changed changed so this seems like a problem specific to OpenBabel.

nbehrnd commented 2 years ago

Any pointers on how to use this?

No, I don't have pointers to offer because so far, my use of OpenBabel for the generation of conformers with confab iterated over all single bonds. Count me out.