ob3.1 compatibility update, added Dockerfile for running automated un…

…it tests

.github/workflows/main.yml

@@ -0,0 +1,18 @@
+name: CI
+on:  
+  pull_request:    
+    branches:
+    - main  
+  push:
+    branches:
+    - main
+jobs:  
+  build:     
+   runs-on: ubuntu-latest
+
+   steps:    
+   - uses: actions/checkout@v1
+   - name: Build and test with Docker
+     run: |
+       docker build --tag ci -f ./DOCKERFILE .
+       docker run --rm ci

Dockerfile

@@ -0,0 +1,55 @@
+# EVOLVE Dockerfile for CI
+
+# Github Actions does not support GPU acceleration and the parallized amber pmemd engine is not available for free via conda. 
+# The unittests regarding this are not executed
+
+
+FROM continuumio/miniconda3
+
+# Setup a spot for the code
+WORKDIR /evolve
+# Install Python dependencies
+COPY environment.yml environment.yml
+RUN conda env create -f environment.yml
+
+# Activate environment
+SHELL ["conda", "run", "-n", "evolve", "/bin/bash", "-c"]
+
+
+# Copy necessary files
+COPY share share/
+COPY unittests unittests/
+COPY main.py .
+COPY evaluators.py .
+
+# This is necessary due to a bug in Docker
+RUN true
+
+COPY src ./src/
+
+# Set Environment variables
+ENV PYTHONPATH="/evolve/:${PYTHONPATH}"
+
+ENV PATH="/opt/conda/envs/evolve/bin:$PATH"
+
+ENV CONDA_DEFAULT_ENV evolve
+
+ENV EVOLVE_GPU=FALSE
+
+SHELL ["conda", "run", "-n", "evolve", "/bin/bash", "-c"]
+
+# Chceck that openbabel is available
+RUN echo "Make sure openbabel is installed:"
+RUN python -c "from openbabel import openbabel"
+
+# Change to unittest directory
+WORKDIR /evolve/unittests
+
+
+# Execute tests
+RUN ["python", "test_gaapi.py"]
+RUN ["python", "test_openbabel.py"]
+RUN ["python", "test_StartGA_Run.py"]
+
+
+RUN ["/bin/bash"]

environment.yml

@@ -0,0 +1,27 @@
+name: evolve
+channels:
+  - conda-forge
+  - omnia
+  - defaults
+dependencies:
+  - openbabel=3.1.0
+  - numpy=1.18.4=py37h8960a57_0
+  - openmm=7.4.1=py37_cuda101_rc_1
+  - openmmtools=0.19.0=py37_1
+  - pagmo=2.15.0=hec06268_1
+  - pip=20.1.1=pyh9f0ad1d_0
+  - ambertools=20.9
+  - pygmo=2.15.0=py37h2566bc3_0
+  - python=3.7.6=h8356626_5_cpython
+  - six=1.15.0=pyh9f0ad1d_0
+  - sphinx=3.0.3=py_0
+  - sphinxcontrib-applehelp=1.0.2=py_0
+  - sphinxcontrib-bibtex=1.0.0=py_0
+  - sphinxcontrib-devhelp=1.0.2=py_0
+  - sphinxcontrib-htmlhelp=1.0.3=py_0
+  - sphinxcontrib-jsmath=1.0.1=py_0
+  - sphinxcontrib-napoleon=0.7=py_0
+  - sphinxcontrib-qthelp=1.0.3=py_0
+  - sphinxcontrib-serializinghtml=1.1.4=py_0
+  - future
+

evaluators.py

@@ -15,7 +15,7 @@
 from src.gaapi import generators
 from src import constants as cnts
 import numpy as np
-import openbabel
+from openbabel import openbabel
 import src.outprocesses as op
 import os
 import subprocess

main.py

@@ -15,7 +15,7 @@
 from past.utils import old_div
 import numpy as np
 import argparse
-import openbabel
+from openbabel import openbabel
 import time
 import json
 import sys

src/JobConfig.py

@@ -19,7 +19,7 @@
 from src.gaapi import operator_types
 import sys
 import os
-import openbabel
+from openbabel import openbabel
 import collections
 from src import constants as cnts
 

src/MoleculeCreator.py

@@ -11,17 +11,19 @@
 
 from builtins import range
 from past.utils import old_div
-import openbabel
-from openbabel import OBResidue
-from openbabel import OBAtom
-from openbabel import OBMol
-from openbabel import OBAtomTyper
-# from openbabel import openbabel
-# from openbabel.openbabel import OBResidue
-# from openbabel.openbabel import OBAtom
-#from openbabel.openbabel import OBMol
-#from openbabel.openbabel import OBAtomTyper
+# import openbabel
+# from openbabel import OBResidue
+# from openbabel import OBAtom
+# from openbabel import OBMol
+# from openbabel import OBAtomTyper
+from openbabel import openbabel
+from openbabel.openbabel import OBResidue
+from openbabel.openbabel import OBAtom
+from openbabel.openbabel import OBMol
+from openbabel.openbabel import OBAtomTyper
 import numpy as np
+import os
+
 from . import MoleculeInfo as mi
 
 
@@ -116,17 +118,19 @@ def add_fragment(mol, fragment, id_start):
         new_atom.SetId(id)
         new_atom.SetVector(frag_atom.GetVector())
         new_atom.SetAtomicNum(frag_atom.GetAtomicNum())
-
         new_atom.SetType(frag_atom.GetType())
         new_atom.SetSpinMultiplicity(frag_atom.GetSpinMultiplicity())
         new_atom.SetFormalCharge(frag_atom.GetFormalCharge())
-        # new_atom.SetImplicitValence(frag_atom.GetImplicitValence())
+
+        #new_atom.SetImplicitValence(frag_atom.GetImplicitValence())
+        new_atom.SetImplicitHCount(frag_atom.GetImplicitHCount())
         new_atom.SetPartialCharge(frag_atom.GetPartialCharge())
 
         for data in frag_atom.GetData():
             new_atom.SetData(data)
 
         id += 1
+        #print(debugAtom(new_atom))
 
     for obbond in openbabel.OBMolBondIter(fragment):
         begin_atom = obbond.GetBeginAtom()
@@ -135,7 +139,7 @@ def add_fragment(mol, fragment, id_start):
 
 
 def debugAtom(obatom):
-    return "AN:", obatom.GetAtomicNum(), "T:", obatom.GetType(), "idx:", obatom.GetIdx(), "id:", obatom.GetId(), "BO(1):", obatom.CountBondsOfOrder(1)
+    return "AN:", obatom.GetAtomicNum(), "T:", obatom.GetType(), obatom.GetResidue().GetName(), "idx:", obatom.GetIdx(), "id:", obatom.GetId(), "BO(1):", obatom.CountBondsOfOrder(1)
 
 
 def debugBond(obbond):
@@ -168,12 +172,11 @@ def swapsidechain(settings, mol, res_index, aa_mol):
     """
 
     mol.BeginModify()
-
+    aa_mol.BeginModify()
     curr = mol.GetResidue(res_index)
 
     new = aa_mol.GetResidue(0)
-    print(mol, mol.NumResidues())
-    print(mol.GetResidue(4))
+
     mol_CA = mi.getAlphaCarbon(curr)
     mol_CB = mi.getBetaAtom(curr)
     mol_N = mi.getBBNitrogen(curr)
@@ -183,11 +186,12 @@ def swapsidechain(settings, mol, res_index, aa_mol):
     aa_bb_nitrogen = mi.getBBNitrogen(new)
     aa_gamma_atom = mi.getChi1DihedralAtom(new)
 
-    frag_res = aa_CB.GetResidue()
+    frag_res = new
     frag_name = frag_res.GetName()
 
-    res_name =mi.getResType(curr)
 
+    res_name =mi.getResType(curr)
+
     aa_tor = 0
     if (aa_gamma_atom is not None):
         aa_tor = aa_mol.GetTorsion(aa_gamma_atom, aa_CA, aa_CB, aa_bb_nitrogen)
@@ -232,6 +236,7 @@ def swapsidechain(settings, mol, res_index, aa_mol):
     frag_atom_ids_del = [aa_mol.GetAtom(i).GetId() for i in frag_atoms_del]
     frag_atom_ids_del.append(aa_CA.GetId())
 
+
     # delete unnecessary atoms
     for i in range (0, len(frag_atom_ids_del)):
 
@@ -243,7 +248,7 @@ def swapsidechain(settings, mol, res_index, aa_mol):
             res.RemoveAtom(atom)
 
         aa_mol.DeleteAtom(atom)
-        
+
     for i in range (0, len(mol_atom_ids_del)):
         atom = getAtomByID(mol, mol_atom_ids_del[i])
 
@@ -258,8 +263,8 @@ def swapsidechain(settings, mol, res_index, aa_mol):
 
     # add fragment to the protein
     add_fragment(mol, aa_mol, orig_num_atoms)
-    
-    mol.EndModify()
+
+    #mol.EndModify()
 
     renum_atoms = openbabel.vectorInt(prev_atoms + aa_mol.NumAtoms())
 
@@ -282,21 +287,30 @@ def swapsidechain(settings, mol, res_index, aa_mol):
         if (atom.GetId() == aa_CB.GetId()):
             corr_frag_cb = atom
 
+
     mol.AddBond(mol_CA.GetIdx(), corr_frag_cb.GetIdx(), 1)
 
     frag_res_type = mi.getResType(frag_res)
-
+    
     if settings.use_res_type == True:
         curr.SetName(frag_res_type)
         pass
     else:
         curr.SetName(frag_name)
-        
+
     for obatom in openbabel.OBResidueAtomIter(frag_res):
         frag_atom = getAtomByID(mol, obatom.GetId())
         curr.AddAtom(frag_atom)
         curr.SetAtomID(frag_atom, frag_res.GetAtomID(obatom))
+
+    # need to renumber IDs now to be consecutive
+    mol.RenumberAtoms(renum_atoms)
 
+    for i in range (1, mol.NumAtoms() + 1) :
+        mol.GetAtom(i).SetId(i - 1) 
+
+    aa_mol.EndModify()
+
 
     # For several residues in order for the amber settings to work, one needs to rename a few hydrogens which is done here.
     # if old resname = glycine we need to rename the remaining H to HA, the other has been replaced by the sidechain
@@ -335,24 +349,34 @@ def swapsidechain(settings, mol, res_index, aa_mol):
         for obatom in openbabel.OBResidueAtomIter(curr):
             if curr.GetAtomID(obatom) in ["HE3", "HE4", "HE5", "HE6", "HE7", "HE8"]:
                 curr.SetAtomID(obatom, "HE1")
-               
+
 
+    # Fix to get around openbabel screwing up the structure when calling EndModify()
+    # This is necessary because we need to set the Chi1 dihedral correctly
+    obConversion = openbabel.OBConversion()
+    obConversion.SetInAndOutFormats("pdb", "pdb")            
+    obConversion.WriteFile(mol, "temp.pdb")
 
+    obConversion = openbabel.OBConversion() 
+    obConversion.SetInAndOutFormats("pdb", "pdb") 
+    obConversion.ReadFile(mol, "temp.pdb")    
+    os.remove('temp.pdb')
 
 
+
+    # preparation for setting the chi1 dihedral
+    curr = mol.GetResidue(res_index)
+
     alpha_carbon = mi.getAlphaCarbon(curr)
     bb_nitrogen = mi.getBBNitrogen(curr)
     beta_atom = mi.getBetaAtom(curr)
     chi_atom = mi.getChi1DihedralAtom(curr)
-
+
+    ## Here the correct dihedral needs to be set! 
     if (chi_atom is not None and aa_gamma_atom is not None):
         mol.SetTorsion(bb_nitrogen, alpha_carbon, beta_atom, chi_atom, aa_tor * (old_div(np.pi, 180.0)))
-
-    # need to renumber IDs now to be consecutive
-    mol.RenumberAtoms(renum_atoms)
 
-    for i in range (1, mol.NumAtoms() + 1) :
-        mol.GetAtom(i).SetId(i - 1)
+
 
 
 

src/MoleculeInfo.py

@@ -12,12 +12,12 @@
 from builtins import str
 from builtins import range
 from past.utils import old_div
-import openbabel
-from openbabel import OBResidue
-from openbabel import OBAtom
-# from openbabel import openbabel
-# from openbabel.openbabel import OBResidue
-# from openbabel.openbabel import OBAtom
+# import openbabel
+# from openbabel import OBResidue
+# from openbabel import OBAtom
+from openbabel import openbabel
+from openbabel.openbabel import OBResidue
+from openbabel.openbabel import OBAtom
 import numpy as np
 from . import constants
 
@@ -86,7 +86,7 @@ def getResType(obres):
         Query residue
     """
     res = obres.GetName()
-    
+
     # check richardson rotamers and return rotamer_type
     for i in range (0, len(constants.rotamers)):
         if (constants.rotamers[i] == res):
@@ -185,7 +185,7 @@ def getBetaAtom(obres):
     obres : OBResidue
         Query residue
     """
-    
+
     alpha_carbon = getAlphaCarbon(obres)
 
     bbcarboxyl = getBBCarboxyl(obres)
@@ -206,7 +206,7 @@ def getBetaAtom(obres):
     carboxyl_vec = np.asarray([bbcarboxyl.GetX(), bbcarboxyl.GetY(), bbcarboxyl.GetZ()])
 
     res = getResType(obres)
-
+    
     for obatom in openbabel.OBAtomAtomIter(alpha_carbon):
         if (res == "GLY"):
             if (obatom.GetType() == 'H'):
@@ -604,6 +604,7 @@ def getChi1DihedralAtom(obres):
     beta_atom = getBetaAtom(obres)
     gamma_atom = None
 
+
     if beta_atom is None:
         return None
 

src/gaapi/Individual.py

@@ -13,7 +13,7 @@
 from src import MoleculeInfo as mi
 from src import MoleculeCreator as mcr
 import copy
-import openbabel
+from openbabel import openbabel
 import numpy as np
 import pygmo as pg
 

src/gaapi/generators.py

@@ -14,7 +14,7 @@
 import numpy as np
 
 from src import MoleculeInfo as mi
-import openbabel
+from openbabel import openbabel
 
 
 def initialisePopulation(settings):