Merge pull request #3 from 9527567/dev

Dev

Author: 9527567

.github/workflows/cmake-single-platform.yml

@@ -0,0 +1,40 @@
+# This starter workflow is for a CMake project running on a single platform. There is a different starter workflow if you need cross-platform coverage.
+# See: https://github.com/actions/starter-workflows/blob/main/ci/cmake-multi-platform.yml
+name: CMake on a single platform
+
+on:
+  push:
+    branches: [ "dev" , "main" ]
+  pull_request:
+    branches: [ "dev" , "main" ]
+
+env:
+  # Customize the CMake build type here (Release, Debug, RelWithDebInfo, etc.)
+  BUILD_TYPE: Release
+
+jobs:
+  build:
+    # The CMake configure and build commands are platform agnostic and should work equally well on Windows or Mac.
+    # You can convert this to a matrix build if you need cross-platform coverage.
+    # See: https://docs.github.com/en/free-pro-team@latest/actions/learn-github-actions/managing-complex-workflows#using-a-build-matrix
+    runs-on: ubuntu-latest
+
+    steps:
+    - uses: actions/checkout@v3
+      with:
+          submodules: true
+    - name: Configure CMake
+      # Configure CMake in a 'build' subdirectory. `CMAKE_BUILD_TYPE` is only required if you are using a single-configuration generator such as make.
+      # See https://cmake.org/cmake/help/latest/variable/CMAKE_BUILD_TYPE.html?highlight=cmake_build_type
+      run: cmake -B ${{github.workspace}}/build -DCMAKE_BUILD_TYPE=${{env.BUILD_TYPE}}
+
+    - name: Build
+      # Build your program with the given configuration
+      run: cmake --build ${{github.workspace}}/build --config ${{env.BUILD_TYPE}}
+
+    - name: Test
+      working-directory: ${{github.workspace}}/build
+      # Execute tests defined by the CMake configuration.
+      # See https://cmake.org/cmake/help/latest/manual/ctest.1.html for more detail
+      run: ctest -C ${{env.BUILD_TYPE}}
+

AmberMDrun/equil.py

@@ -1,17 +1,26 @@
 from . import pyamber
 import os
 import pandas as pd
+from pathlib import Path
+
+
 def density():
-    input = f"for FILE in final_?.out,final_??.out \n\
-readdata \$FILE name MD \n\
-done \n\
-evalplateau *[Density] name EQ out Eval.agr resultsout Eval.results\n\
-go\n\
-quit"
+    path = Path.cwd()
+    final_file = []
+    for out_file in path.glob('final_*.out'):
+        final_file.append(out_file)
+    input = (f"for FILE in {' '.join([item.name for item in final_file])} \n"
+             "readdata \$FILE name MD \n"
+             "done \n"
+             "evalplateau *[Density] name EQ out Eval.agr resultsout Eval.results\n"
+             "go\n"
+             "quit"
+             )
+
     with open("cpptraj.in", "w") as f:
         f.write(input)
     os.system(f'cpptraj -i cpptraj.in')
-    result = pd.read_csv("Eval.results",sep="\s+")
+    result = pd.read_csv("Eval.results", sep="\s+")
     if result["EQ[result]"][0] == "yes":
         return 0
     elif result["EQ[result]"][0] == "no":
@@ -52,10 +61,12 @@ def prep(rst7, s, temp, heavymask, backbonemask, loop=20):
     ref = "step9.rst7"
     for i in range(loop):
         final = pyamber.NPT(f"final_{i}", systemInfo=s, ref=ref, temp=temp,
-                         refc="step5.rst7", irest=True, dt=0.002, nscm=1000, nstlim=500000, ntwx=5000)
+                            refc="step5.rst7", irest=True, dt=0.002, nscm=1000, nstlim=500000, ntwx=5000)
         final.Run()
         result = density()
         if result == 0:
-            return f'final_{i}.rst7'    
+            return f'final_{i}.rst7'
         ref = f'final_{i}.rst7'
-    raise RuntimeError("More than 20 iterations of final density equil required. Bailing out.")
+    raise RuntimeError(
+        "More than 20 iterations of final density equil required. Bailing out.")
+

AmberMDrun/extern/pybind11

@@ -6,6 +6,7 @@
 import logging
 from logging import getLogger
 import subprocess
+from typing import List
 
 
 def runCMD(inCmd, *, raise_on_fail: bool = True, logger: logging.Logger = getLogger("mmpbsa"), **kwargs):
@@ -61,36 +62,45 @@ def split_pdb(pdb: str):
     return "pro.pdb", "mol.mol2"
 
 
-def run_tleap(protein: str, mol: str, charge: int, multiplicity: int):
+def run_tleap(protein: str, mol_list: List, charge: List, multiplicity: List, guess_charge: bool):
     cmdline = f'pdb4amber -i {protein} -o _{str(protein)} -y -d -p'
     runCMD(cmdline)
     protein_path = Path(protein).absolute()
-    mol_path = Path(mol).absolute()
-    cmdline = f'acpype -i {str(mol_path)} -n {charge} -m {multiplicity}'
-    runCMD(cmdline,
-           message="Perhaps you should check the charge of the ligand and the correctness of the hydrogen atom.")
-    leapin = f"source leaprc.protein.ff14SB\n\
-            source leaprc.DNA.OL15\n\
-            source leaprc.RNA.OL3\n\
-            source leaprc.water.tip3p\n\
-            source leaprc.gaff2\n\
-            pro = loadpdb _{protein}\n\
-            loadamberparams {mol_path.stem}.acpype/{mol_path.stem}_AC.frcmod\n\
-            mol = loadmol2 {mol_path.stem}.acpype/{mol_path.stem}_bcc_gaff2.mol2\n\
-            com = combine{{pro mol}}\n\
-            solvatebox com TIP3PBOX 10.0\n\
-            addions2 com Na+ 0\n\
-            addions2 com Cl- 0\n\
-            saveamberparm com {protein_path.stem}_{mol_path.stem}.parm7 {protein_path.stem}_{mol_path.stem}.rst7\n\
-            quit"
+    mol_list = [Path(mol) for mol in mol_list]
+    if guess_charge:
+        for mol in mol_list:
+            cmdline = f'acpype -i {str(Path(mol).absolute())}'
+            runCMD(cmdline,
+                   message="Perhaps you should check the charge of the ligand and the correctness of the hydrogen atom.")
+    else:
+        for mol, c, spin in zip(mol_list, charge, multiplicity):
+            cmdline = f'acpype -i {str(Path(mol).absolute())} -n {c} -m {spin}'
+            runCMD(cmdline,
+                   message="Perhaps you should check the charge of the ligand and the correctness of the hydrogen atom.")
+
+    mol_frcmod = f"".join(
+        f'loadamberparams {mol_path.stem}.acpype/{mol_path.stem}_AC.frcmod\n' for mol_path in mol_list)
+    mol_load = f"".join(
+        f'{mol_path.stem} = loadmol2 {mol_path.stem}.acpype/{mol_path.stem}_bcc_gaff2.mol2\n' for mol_path in mol_list)
+    combine = f'com = combine{{pro {" ".join(mol_path.stem for mol_path in mol_list)}}}\n'
+    leapin = (f"""source leaprc.protein.ff14SB
+source leaprc.DNA.OL15
+source leaprc.RNA.OL3
+source leaprc.water.tip3p
+source leaprc.gaff2
+pro = loadpdb _{protein}\n""") + mol_frcmod + mol_load + combine + (f"""solvatebox com TIP3PBOX 10.0
+addionsrand com Na+ 0
+addionsrand com Cl- 0
+saveamberparm com com.parm7 com.rst7
+quit""")
     with open("leap.in", "w") as f:
         for i in leapin:
             f.write(i)
     runCMD('tleap -f leap.in')
-    return f'{protein_path.stem}_{mol_path.stem}.parm7', f'{protein_path.stem}_{mol_path.stem}.rst7'
+    return f'com.parm7', f'com.rst7'
 
 
-def mmpbsa(parm7: str, rst7: str, netcdf: str, system: pyamber.SystemInfo):
+def run_mmpbsa(parm7: str, rst7: str, netcdf: str, system: pyamber.SystemInfo, mol_list: List):
     parm7 = Path(parm7).absolute()
     rst7 = Path(rst7).absolute()
 
@@ -138,25 +148,36 @@ def mmpbsa(parm7: str, rst7: str, netcdf: str, system: pyamber.SystemInfo):
     with open("mmpbsa.in", 'w') as f:
         for i in mmpbsa_in:
             f.write(i)
-    mmpbsa = f"mpirun -np {cpu_count() // 2} gmx_MMPBSA MPI -O -i mmpbsa.in -cs {str(parm7.with_suffix('.pdb'))} -ci index.ndx -cg 1 13 -ct {str(parm7.with_suffix('.xtc'))}  -cp \
-    {str(parm7.with_suffix('.top'))} -nogui"
-    runCMD(mmpbsa)
+    mol_number = 13
+    for mol in range(len(mol_list)):
+        mol_path = Path.cwd().joinpath(f'lig{mol}')
+        mol_path.mkdir(exist_ok=True)
+        os.chdir(str(mol_path))
+        mmpbsa = f"mpirun -np {cpu_count() // 2} gmx_MMPBSA MPI -O -i ../mmpbsa.in -cs {str(parm7.with_suffix('.pdb'))} -ci ../index.ndx -cg 1 {mol_number} -ct {str(parm7.with_suffix('.xtc'))}  -cp \
+        {str(parm7.with_suffix('.top'))} -nogui"
+        runCMD(mmpbsa)
+        os.chdir('..')
+        mol_number += 1
 
 
 def arg_parse():
-    parser = argparse.ArgumentParser(description='Tools for automating the operation of MMPBSA')
+    parser = argparse.ArgumentParser(
+        description='Tools for automating the operation of MMPBSA')
     parser.add_argument('--protein', '-p', type=str,
                         required=True, help="pdb file for protein")
-    parser.add_argument('--mol2', '-m', type=str,
+    parser.add_argument('--mol2', '-m', type=str, nargs='+',
                         required=False, help='mol2 file for mol')
     parser.add_argument('--temp', '-t', type=float,
                         required=False, help='Temperature', default=303.15)
     parser.add_argument("--ns", '-n', type=int,
                         help="time for MD(ns)", default=100)
-    parser.add_argument("--charge", type=int,
-                        default=0, help="charge of mol")
-    parser.add_argument("--multiplicity", type=int,
-                        default=1, help="multiplicity of mol")
+    parser.add_argument('-g', '--guess_charge',
+                        action='store_true', help='guess charge')
+    parser.add_argument('-c', "--charge", type=int, nargs='+',
+                        default=[0], help="charge of mol")
+
+    parser.add_argument("--multiplicity", type=int, nargs='+',
+                        default=[1], help="multiplicity of mol")
     parser.add_argument("--MIN", type=str,
                         default="pmemd.cuda_DPFP", help="Engine for MIN")
     parser.add_argument("--MD", type=str,
@@ -168,11 +189,18 @@ def arg_parse():
 def mmpbsa():
     args = arg_parse()
     protein = args.protein
-    mol = args.mol2
+    mol_list = args.mol2
     temp = args.temp
-    if mol is None:
+    if not args.guess_charge:
+        if len(mol_list) != len(args.charge) and len(mol_list) != len(args.multiplicity):
+            raise ValueError(
+                "If the charge is not guessed, it is necessary to specify the charge and spin multiplicity for each ligand.")
+
+    if mol_list is None:
         protein, mol = split_pdb(protein)
-    parm7, rst7 = run_tleap(protein, mol, args.charge, args.multiplicity)
+        mol_list = [mol]
+    parm7, rst7 = run_tleap(protein, mol_list, args.charge,
+                            args.multiplicity, args.guess_charge)
     s = pyamber.SystemInfo(parm7, rst7, runMin=args.MIN, runMd=args.MD)
     heavymask = "\"" + s.getHeavyMask() + "\""
     backbonemask = "\"" + s.getBackBoneMask() + "\""
@@ -181,7 +209,8 @@ def mmpbsa():
     md = pyamber.NPT("md", s, rst7, rst7, ntwx=50000,
                      irest=True, nscm=1000, nstlim=args.ns * 500000)
     md.Run()
-    mmpbsa(parm7, rst7, "md.nc", s)
+    rst7 = 'final_0.rst7'
+    run_mmpbsa(parm7, rst7, "md.nc", s, mol_list)
 
 
 if __name__ == '__main__':

AmberMDrun/mmpbsa.py

@@ -1 +1 @@
-__version__ = '0.0.4'
+__version__ = '0.0.5'

AmberMDrun/version.py

@@ -3,6 +3,8 @@
 
 # AmberMDrun 
 Easy to use, easy to expand, high-performance Amber simulation package
+## Update 
+v0.0.5 Added support for multiple ligands.
 ## Install
 This software only supports **Linux** because some Linux system functions are called.**Mac OS X** and **Windows** are not supported.
 ### Necessary
@@ -77,19 +79,22 @@ options:
 ## How to calculate MM-PB (GB) SA between small molecules and proteins of a single drug
 
 ~~~bash
-usage: mmpbsa [-h] --protein PROTEIN [--mol2 MOL2] [--temp TEMP] [--ns NS] [--charge CHARGE] [--multiplicity MULTIPLICITY] [--MIN MIN] [--MD MD]
+usage: mmpbsa [-h] --protein PROTEIN [--mol2 MOL2 [MOL2 ...]] [--temp TEMP] [--ns NS] [-g] [-c CHARGE [CHARGE ...]] [--multiplicity MULTIPLICITY [MULTIPLICITY ...]] [--MIN MIN] [--MD MD]
 
 Tools for automating the operation of MMPBSA
 
 options:
   -h, --help            show this help message and exit
   --protein PROTEIN, -p PROTEIN
                         pdb file for protein
-  --mol2 MOL2, -m MOL2  mol2 file for mol
+  --mol2 MOL2 [MOL2 ...], -m MOL2 [MOL2 ...]
+                        mol2 file for mol
   --temp TEMP, -t TEMP  Temperature
   --ns NS, -n NS        time for MD(ns)
-  --charge CHARGE       charge of mol
-  --multiplicity MULTIPLICITY
+  -g, --guess_charge    guess charge
+  -c CHARGE [CHARGE ...], --charge CHARGE [CHARGE ...]
+                        charge of mol
+  --multiplicity MULTIPLICITY [MULTIPLICITY ...]
                         multiplicity of mol
   --MIN MIN             Engine for MIN
   --MD MD               Engine for MD
@@ -98,6 +103,11 @@ Typically, the complex structure after molecular docking is used to perform MMPB
 ~~~bash
 mmpbsa -p complex.pdb
 ~~~
+## V0.0.5 added support for multiple ligands
+Just follow the files of multiple ligands after -m, and add an option `-g` to guess the static charge of small molecules, or manually specify the static charge, for example:
+~~~bash
+mmpbsa -p pro.pdb -m lig1.mol2 lig2.mol2 -g -n 100
+~~~
 ## How to extend code through inheritance classes
 Will be described in the near future
 
@@ -116,4 +126,20 @@ URL = {https://www.mdpi.com/2218-273X/13/4/635},
 ISSN = {2218-273X},
 DOI = {10.3390/biom13040635}
 }
+~~~
+## If you are interested, you can also cite this article
+~~~tex
+@article{CUI2023134812,
+title = {A TastePeptides-Meta system including an umami/bitter classification model Umami_YYDS, a TastePeptidesDB database and an open-source package Auto_Taste_ML},
+journal = {Food Chemistry},
+volume = {405},
+pages = {134812},
+year = {2023},
+issn = {0308-8146},
+doi = {https://doi.org/10.1016/j.foodchem.2022.134812},
+url = {https://www.sciencedirect.com/science/article/pii/S0308814622027741},
+author = {Zhiyong Cui and Zhiwei Zhang and Tianxing Zhou and Xueke Zhou and Yin Zhang and Hengli Meng and Wenli Wang and Yuan Liu},
+keywords = {Peptides, Umami prediction, TastePeptidesDB, Machine learning},
+abstract = {Taste peptides with umami/bitterness play a role in food attributes. However, the taste mechanisms of peptides are not fully understood, and the identification of these peptides is time-consuming. Here, we created a taste peptide database by collecting the reported taste peptide information. Eight key molecular descriptors from di/tri-peptides were selected and obtained by modeling screening. A gradient boosting decision tree model named Umami_YYDS (89.6\% accuracy) was established by data enhancement, comparison algorithm and model optimization. Our model showed a great prediction performance compared to other models, and its outstanding ability was verified by sensory experiments. To provide a convenient approach, we deployed a prediction website based on Umami_YYDS and uploaded the Auto_Taste_ML machine learning package. In summary, we established the system TastePeptides-Meta, containing a taste peptide database TastePeptidesDB an umami/bitter taste prediction model Umami_YYDS and an open-source machine learning package Auto_Taste_ML, which were helpful for rapid screening of umami peptides.}
+}
 ~~~

README.md

@@ -1,5 +1,7 @@
 # AmberMDrun 
 易于使用、易于扩展、高性能的Amber模拟软件包。
+## 版本更新
+v0.0.5 添加了多配体的支持。
 ## 安装
 此软件仅支持**Linux**，因为某些Linux系统功能被调用。**Mac OS X**和**Windows**不受支持。
 ### 必要的依赖
@@ -76,19 +78,22 @@ options:
 ~~~
 ## 如何计算单个小分子和蛋白质之间的MM-PB（GB）SA
 ~~~bash
-usage: mmpbsa [-h] --protein PROTEIN [--mol2 MOL2] [--temp TEMP] [--ns NS] [--charge CHARGE] [--multiplicity MULTIPLICITY] [--MIN MIN] [--MD MD]
+usage: mmpbsa [-h] --protein PROTEIN [--mol2 MOL2 [MOL2 ...]] [--temp TEMP] [--ns NS] [-g] [-c CHARGE [CHARGE ...]] [--multiplicity MULTIPLICITY [MULTIPLICITY ...]] [--MIN MIN] [--MD MD]
 
 Tools for automating the operation of MMPBSA
 
 options:
   -h, --help            show this help message and exit
   --protein PROTEIN, -p PROTEIN
                         pdb file for protein
-  --mol2 MOL2, -m MOL2  mol2 file for mol
+  --mol2 MOL2 [MOL2 ...], -m MOL2 [MOL2 ...]
+                        mol2 file for mol
   --temp TEMP, -t TEMP  Temperature
   --ns NS, -n NS        time for MD(ns)
-  --charge CHARGE       charge of mol
-  --multiplicity MULTIPLICITY
+  -g, --guess_charge    guess charge
+  -c CHARGE [CHARGE ...], --charge CHARGE [CHARGE ...]
+                        charge of mol
+  --multiplicity MULTIPLICITY [MULTIPLICITY ...]
                         multiplicity of mol
   --MIN MIN             Engine for MIN
   --MD MD               Engine for MD
@@ -97,6 +102,12 @@ options:
 ~~~bash
 mmpbsa -p complex.pdb
 ~~~
+
+## V0.0.5 添加了多配体的支持
+只需要在-m 后跟多个配体的文件即可，添加了一个选项`-g`用于猜测小分子的静电荷，或者手动指定静电荷，例如：
+~~~bash
+mmpbsa -p pro.pdb -m lig1.mol2 lig2.mol2 -g -n 100
+~~~
 ## 如何通过继承扩展代码
 我们将在不久的将来进行描述。
 
@@ -115,4 +126,21 @@ URL = {https://www.mdpi.com/2218-273X/13/4/635},
 ISSN = {2218-273X},
 DOI = {10.3390/biom13040635}
 }
+~~~
+## 如果您感兴趣的话，也可以引用这篇文章
+bibtex:
+~~~tex
+@article{CUI2023134812,
+title = {A TastePeptides-Meta system including an umami/bitter classification model Umami_YYDS, a TastePeptidesDB database and an open-source package Auto_Taste_ML},
+journal = {Food Chemistry},
+volume = {405},
+pages = {134812},
+year = {2023},
+issn = {0308-8146},
+doi = {https://doi.org/10.1016/j.foodchem.2022.134812},
+url = {https://www.sciencedirect.com/science/article/pii/S0308814622027741},
+author = {Zhiyong Cui and Zhiwei Zhang and Tianxing Zhou and Xueke Zhou and Yin Zhang and Hengli Meng and Wenli Wang and Yuan Liu},
+keywords = {Peptides, Umami prediction, TastePeptidesDB, Machine learning},
+abstract = {Taste peptides with umami/bitterness play a role in food attributes. However, the taste mechanisms of peptides are not fully understood, and the identification of these peptides is time-consuming. Here, we created a taste peptide database by collecting the reported taste peptide information. Eight key molecular descriptors from di/tri-peptides were selected and obtained by modeling screening. A gradient boosting decision tree model named Umami_YYDS (89.6\% accuracy) was established by data enhancement, comparison algorithm and model optimization. Our model showed a great prediction performance compared to other models, and its outstanding ability was verified by sensory experiments. To provide a convenient approach, we deployed a prediction website based on Umami_YYDS and uploaded the Auto_Taste_ML machine learning package. In summary, we established the system TastePeptides-Meta, containing a taste peptide database TastePeptidesDB an umami/bitter taste prediction model Umami_YYDS and an open-source machine learning package Auto_Taste_ML, which were helpful for rapid screening of umami peptides.}
+}
 ~~~

README.zh.md

@@ -125,7 +125,7 @@ def build_extension(self, ext: CMakeExtension) -> None:
 # logic and declaration, and simpler if you include description/version in a file.
 setup(
     name="AmberMDrun",
-    version="0.0.4",
+    version="0.0.5",
     author="ZhiWei Zhang",
     author_email="[email protected]",
     description="A scripting tool for running Amber MD in an easy way",