Ambertools-OpenMM-MD_GBSA.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "2156ddee",
      "metadata": {
        "id": "2156ddee"
      },
      "source": [
        "# Ambertools-OpenMM Protein-ligand complex MD Setup tutorial using Jupyter Notebook\n",
        "This notebook is wrritten by quantaosun@gmail.com in Shanghai, China, 2021. This notebook took PDB ID 7L10, the protease of SARS-CoV-2 with an inhibitor, as an example, relative work was firstly published on ACS Cent. Sci. 2021, 7, 3, 467–475, by William L. Jorgensen, Yale University.\n",
        "\n",
        "Ambertools (https://ambermd.org/AmberTools.php) and Openmm (https://openmm.org/) should be properly refer to if any simulation result would be used in any publications."
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Please modify first 4 cells if you use local Linux instead of Google Colab\n"
      ],
      "metadata": {
        "id": "Jzn6LXx4EmRb"
      },
      "id": "Jzn6LXx4EmRb"
    },
    {
      "cell_type": "code",
      "source": [
        "!pip install condacolab\n",
        "import condacolab\n",
        "condacolab.install_mambaforge()"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "NXpwM8al_V7e",
        "outputId": "c8a26d38-5eb6-433b-cd4f-33ed97fce53e"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Collecting condacolab\n",
            "  Downloading condacolab-0.1.7-py3-none-any.whl (7.2 kB)\n",
            "Installing collected packages: condacolab\n",
            "Successfully installed condacolab-0.1.7\n",
            "⏬ Downloading https://github.com/conda-forge/miniforge/releases/download/23.1.0-1/Mambaforge-23.1.0-1-Linux-x86_64.sh...\n",
            "📦 Installing...\n",
            "📌 Adjusting configuration...\n",
            "🩹 Patching environment...\n",
            "⏲ Done in 0:00:12\n",
            "🔁 Restarting kernel...\n"
          ]
        }
      ],
      "id": "NXpwM8al_V7e"
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "d545a63b"
      },
      "outputs": [],
      "source": [
        "#install openmm\n",
        "!mamba install -c conda-forge openmm -y\n",
        "# install openmm-setup for protein fix and simulation configuration\n",
        "!mamba install -c conda-forge openmm-setup -y\n",
        "#install open babel\n",
        "!mamba install -c openbabel openbabel -y\n",
        "# install Amber tools\n",
        "!mamba install -c conda-forge ambertools -y\n",
        "#install visualization\n",
        "! pip install py3Dmol"
      ],
      "id": "d545a63b"
    },
    {
      "cell_type": "code",
      "source": [
        "amberhome = \"source /usr/local/amber.sh\""
      ],
      "metadata": {
        "id": "swFTi8hLDpAS"
      },
      "execution_count": null,
      "outputs": [],
      "id": "swFTi8hLDpAS"
    },
    {
      "cell_type": "code",
      "source": [
        "\n",
        "!MMPBSA.py -h"
      ],
      "metadata": {
        "id": "26CpONxKDsCf"
      },
      "execution_count": null,
      "outputs": [],
      "id": "26CpONxKDsCf"
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "2d70a9e8",
      "metadata": {
        "id": "2d70a9e8",
        "outputId": "d0af2830-4617-4c4c-fcb6-eff61e91084d"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "# conda environments:\n",
            "#\n",
            "base                     /Users/quantaosun/opt/anaconda3\n",
            "Ambertools-OpenMM-MD  *  /Users/quantaosun/opt/anaconda3/envs/Ambertools-OpenMM-MD\n",
            "\n"
          ]
        }
      ],
      "source": [
        "!conda env list"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "965bcdaf",
      "metadata": {
        "id": "965bcdaf"
      },
      "outputs": [],
      "source": [
        "import sys\n",
        "sys.path.append('/home/sqt/anaconda3/envs/Ambertools-OpenMM-MD/lib/python3.7/site-packages')"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "4aa960a3",
      "metadata": {
        "id": "4aa960a3"
      },
      "source": [
        "# Download your PDB, if the command failed, just mannuly download and skip this."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "41750bad",
      "metadata": {
        "id": "41750bad"
      },
      "outputs": [],
      "source": []
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "cb8bf286",
      "metadata": {
        "id": "cb8bf286",
        "outputId": "e9db277e-a5f1-4a48-dcb5-b0bbca3f7eb2"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "zsh:1: command not found: wget\r\n"
          ]
        }
      ],
      "source": [
        "PDB_ID = \"7L10\" #@param {type:\"string\"}\n",
        "pdb = PDB_ID + \".pdb\"\n",
        "wget https://files.rcsb.org/download/$pdb"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "744c6f2b",
      "metadata": {
        "id": "744c6f2b"
      },
      "source": [
        "# Fix your protein with PDBfixer with OpenMM Setup"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "458651e0",
      "metadata": {
        "id": "458651e0",
        "outputId": "50fe7e45-a6fe-4183-d760-e9e8531f12f2"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "PDBFixer running: http://localhost:8000 \n",
            "Warning: importing 'simtk.openmm' is deprecated.  Import 'openmm' instead.\n",
            "^C\n",
            "Traceback (most recent call last):\n",
            "  File \"/Users/quantaosun/opt/anaconda3/envs/Ambertools-OpenMM-MD/bin/pdbfixer\", line 10, in <module>\n",
            "    sys.exit(main())\n",
            "  File \"/Users/quantaosun/opt/anaconda3/envs/Ambertools-OpenMM-MD/lib/python3.7/site-packages/pdbfixer/pdbfixer.py\", line 1212, in main\n",
            "    ui.launchUI()\n",
            "  File \"/Users/quantaosun/opt/anaconda3/envs/Ambertools-OpenMM-MD/lib/python3.7/site-packages/pdbfixer/ui.py\", line 282, in launchUI\n",
            "    time.sleep(0.5)\n",
            "KeyboardInterrupt\n"
          ]
        }
      ],
      "source": [
        "!pdbfixer\n",
        "#Click the square pause button to force quit after finished."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "41d7fae7",
      "metadata": {
        "id": "41d7fae7"
      },
      "source": [
        "# Rename processed pdb to \"output.pdb\""
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "ed717745",
      "metadata": {
        "id": "ed717745"
      },
      "outputs": [],
      "source": [
        "#!mv 7L10-processed.pdb output.pdb"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6eea0905",
      "metadata": {
        "id": "6eea0905"
      },
      "source": [
        "# Go to H++, register and process the output.pdb structure, download top and crd. If you leave everything default, you don't need to change anything below. http://newbiophysics.cs.vt.edu/H++/uploadpdb.php"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "03f24598",
      "metadata": {
        "id": "03f24598"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2011-02-08.png![image-2.png](attachment:image-2.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ba83d6a1",
      "metadata": {
        "id": "ba83d6a1"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2011-03-21.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "32dd1d44",
      "metadata": {
        "id": "32dd1d44"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2011-44-00.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6f115124",
      "metadata": {
        "id": "6f115124"
      },
      "source": [
        "# It is assumed you have save the processed pdb with a name of \"output.pdb\", otherwise please modifiy commands below."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2cc109fd",
      "metadata": {
        "id": "2cc109fd"
      },
      "source": [
        "# Generate protein pdb"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "7dd3d96a",
      "metadata": {
        "id": "7dd3d96a"
      },
      "outputs": [],
      "source": [
        "!ambpdb -p 0.15_80_10_pH6.5_output.top -c 0.15_80_10_pH6.5_output.crd > 0.pdb"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "f3367d14",
      "metadata": {
        "id": "f3367d14"
      },
      "outputs": [],
      "source": [
        "!awk '$4==\"XEY\"' 7L10.pdb > ligand1.pdb"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "38cd9803",
      "metadata": {
        "id": "38cd9803"
      },
      "source": [
        "# Add hydrogens to you ligand"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "bdeeef5e",
      "metadata": {
        "id": "bdeeef5e"
      },
      "outputs": [],
      "source": [
        "!obabel -ipdb ligand1.pdb -opdb -O ligand2.pdb -h"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "b2b10aa2",
      "metadata": {
        "id": "b2b10aa2"
      },
      "source": [
        "# Beautify pdb format"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "30f01d70",
      "metadata": {
        "id": "30f01d70"
      },
      "outputs": [],
      "source": [
        "!pdb4amber -i ligand2.pdb -o ligand3.pdb"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "255837f3",
      "metadata": {
        "id": "255837f3"
      },
      "source": [
        "# Generate ligand mol2 file"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "b7ec9c9a",
      "metadata": {
        "id": "b7ec9c9a"
      },
      "outputs": [],
      "source": [
        "!antechamber -fi pdb -i ligand3.pdb -fo mol2 -o ligand4.mol2 -c bcc -pf y"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "eb8e40ff",
      "metadata": {
        "id": "eb8e40ff"
      },
      "source": [
        "# Make sure you see \"calculation completed\", otherwise something may have been wrong"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "c2564b16",
      "metadata": {
        "id": "c2564b16"
      },
      "outputs": [],
      "source": [
        "!tail  sqm.out"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "f9e71278",
      "metadata": {
        "id": "f9e71278"
      },
      "source": [
        "# Generate ligand prepi file"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "6d0107dc",
      "metadata": {
        "id": "6d0107dc"
      },
      "outputs": [],
      "source": [
        "!antechamber -i ligand4.mol2 -fi mol2 -o ligand5.prepi -fo prepi -pf y"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "22bd5bf8",
      "metadata": {
        "id": "22bd5bf8"
      },
      "source": [
        "# Generate ligand frcmod file"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "8805b93b",
      "metadata": {
        "id": "8805b93b"
      },
      "outputs": [],
      "source": [
        "!parmchk2 -f prepi -i ligand5.prepi -o ligand6.frcmod"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "20aa5c15",
      "metadata": {
        "id": "20aa5c15"
      },
      "source": [
        "# Combine protein with ligand again"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "31d546e1",
      "metadata": {
        "id": "31d546e1"
      },
      "outputs": [],
      "source": [
        "!cat 0.pdb ligand3.pdb > complex.pdb"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c0beb0f9",
      "metadata": {
        "id": "c0beb0f9"
      },
      "source": [
        "# (Not working properly for now) Beautify PDB format, but we might just not need this anyway, skip for now."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "d8cd3645",
      "metadata": {
        "id": "d8cd3645"
      },
      "outputs": [],
      "source": [
        "#!pdb4amber -i complex.pdb -o complex_clean.pdb\n",
        "# This step is problematic for now, it has been observed it is not necessary or even could bring errors to next step."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0be7464f",
      "metadata": {
        "id": "0be7464f"
      },
      "source": [
        "# Create and run a file called tleap.in with Ambertools"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "04e397f5",
      "metadata": {
        "id": "04e397f5"
      },
      "outputs": [],
      "source": [
        "com_file = open('tleap.in','w')\n",
        "com_file.write('''\n",
        "source leaprc.protein.ff14SB #Source leaprc file for ff14SB protein force field\n",
        "source leaprc.gaff #Source leaprc file for gaff\n",
        "source leaprc.water.tip3p #Source leaprc file for TIP3P water model\n",
        "loadamberprep ligand5.prepi #Load the prepi file for the ligand\n",
        "loadamberparams ligand6.frcmod #Load the additional frcmod file for ligand\n",
        "mol = loadpdb complex.pdb #Load PDB file for protein-ligand complex\n",
        "solvatebox mol TIP3PBOX 8 #Solvate the complex with a cubic water box\n",
        "addions mol Cl- 0 #Add Cl- ions to neutralize the system\n",
        "saveamberparm mol complex.prmtop complex.inpcrd #Save AMBER topology and coordinate files\n",
        "quit #Quit tleap program\n",
        "''')\n",
        "com_file.close()\n",
        "!tleap -s -f tleap.in > tleap.out"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "38a0922f",
      "metadata": {
        "id": "38a0922f"
      },
      "source": [
        "# Upload prmtop and inpcrd to openmm-setup, and start the simulation! If you have a proper GPU, select CUDA, otherwise select CPU."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "02255657",
      "metadata": {
        "id": "02255657"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2011-00-12.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0d5e7887",
      "metadata": {
        "id": "0d5e7887"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2011-01-00.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6db5d460",
      "metadata": {
        "id": "6db5d460"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2012-02-57.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "9b142011",
      "metadata": {
        "id": "9b142011"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2012-03-57.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "8dd40cf9",
      "metadata": {
        "id": "8dd40cf9"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2012-04-37.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a8815182",
      "metadata": {
        "id": "a8815182"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2012-05-08.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a13daa24",
      "metadata": {
        "id": "a13daa24"
      },
      "source": [
        "# You could leave all parameters by default and click \"Run Simulation\""
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0e6745de",
      "metadata": {
        "id": "0e6745de"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2012-05-36.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "3c867c26",
      "metadata": {
        "id": "3c867c26"
      },
      "source": [
        "# Click \"Start\""
      ]
    },
    {
      "cell_type": "markdown",
      "id": "53ecc7b7",
      "metadata": {
        "id": "53ecc7b7"
      },
      "source": [
        "file:///home/sqt/Pictures/Screenshot%20from%202021-12-03%2012-06-44.png![image.png](attachment:image.png)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "4fe6496a",
      "metadata": {
        "id": "4fe6496a"
      },
      "outputs": [],
      "source": [
        "!openmm-setup\n",
        "#Click the square pause button to force quit after finished."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "90eef1c1",
      "metadata": {
        "id": "90eef1c1"
      },
      "source": [
        "# If all goes well, you should see something like below (CPU only)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "21473154",
      "metadata": {
        "id": "21473154"
      },
      "outputs": [],
      "source": [
        "Building system...\n",
        "Performing energy minimization...\n",
        "Equilibrating...\n",
        "Simulating...\n",
        "#\"Progress (%)\"\t\"Step\"\t\"Potential Energy (kJ/mole)\"\t\"Temperature (K)\"\t\"Speed (ns/day)\"\n",
        "0.1%\t1000\t-488259.139341292\t298.57624202824957\t0\n",
        "0.2%\t2000\t-492082.04384769656\t300.09275305741824\t16.6\n",
        "0.3%\t3000\t-494725.1462343714\t302.3090730320394\t16.6\n",
        "0.4%\t4000\t-496075.06792005355\t302.3939891905172\t16.4\n",
        "0.5%\t5000\t-498015.23323377833\t300.135566206444\t16.4\n",
        "0.6%\t6000\t-498488.50188332796\t298.5545999644654\t16.3\n",
        "0.7%\t7000\t-499572.9765659773\t300.7757759109986\t16.3\n",
        "0.8%\t8000\t-499030.79245316726\t297.2492234858693\t16.2\n",
        "0.9%\t9000\t-499823.2656625096\t300.60293949805\t16.2\n",
        "1.0%\t10000\t-500134.98822511156\t297.5748675220697\t16.2\n",
        "1.1%\t11000\t-499849.7045947542\t300.5532217535929\t16.2\n",
        "1.2%\t12000\t-500542.2799043103\t299.2530600741315\t16.2"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "61e7f2e4",
      "metadata": {
        "id": "61e7f2e4"
      },
      "source": [
        "# Alternatively, you could save the simulation input files and run it on a cloud or HPC based OpenMM with GPU acceleration, an example script would be :"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "6322eb4f",
      "metadata": {
        "id": "6322eb4f"
      },
      "outputs": [],
      "source": [
        "# This script was generated by OpenMM-Setup on 2021-12-15.\n",
        "\n",
        "from simtk.openmm import *\n",
        "from simtk.openmm.app import *\n",
        "from simtk.unit import *\n",
        "\n",
        "# Input Files\n",
        "\n",
        "prmtop = AmberPrmtopFile('complex.prmtop')\n",
        "inpcrd = AmberInpcrdFile('complex.inpcrd')\n",
        "\n",
        "# System Configuration\n",
        "\n",
        "nonbondedMethod = PME\n",
        "nonbondedCutoff = 1.0*nanometers\n",
        "ewaldErrorTolerance = 0.0005\n",
        "constraints = HBonds\n",
        "rigidWater = True\n",
        "constraintTolerance = 0.000001\n",
        "hydrogenMass = 1.5*amu\n",
        "\n",
        "# Integration Options\n",
        "\n",
        "dt = 0.004*picoseconds\n",
        "temperature = 300*kelvin\n",
        "friction = 1.0/picosecond\n",
        "pressure = 1.0*atmospheres\n",
        "barostatInterval = 25\n",
        "\n",
        "# Simulation Options\n",
        "\n",
        "steps = 1000000\n",
        "equilibrationSteps = 1000\n",
        "platform = Platform.getPlatformByName('CPU')\n",
        "dcdReporter = DCDReporter('trajectory.dcd', 10000)\n",
        "dataReporter = StateDataReporter('log.txt', 1000, totalSteps=steps,\n",
        "    step=True, speed=True, progress=True, potentialEnergy=True, temperature=True, separator='\\t')\n",
        "checkpointReporter = CheckpointReporter('checkpoint.chk', 10000)\n",
        "\n",
        "# Prepare the Simulation\n",
        "\n",
        "print('Building system...')\n",
        "topology = prmtop.topology\n",
        "positions = inpcrd.positions\n",
        "system = prmtop.createSystem(nonbondedMethod=nonbondedMethod, nonbondedCutoff=nonbondedCutoff,\n",
        "    constraints=constraints, rigidWater=rigidWater, ewaldErrorTolerance=ewaldErrorTolerance, hydrogenMass=hydrogenMass)\n",
        "system.addForce(MonteCarloBarostat(pressure, temperature, barostatInterval))\n",
        "integrator = LangevinMiddleIntegrator(temperature, friction, dt)\n",
        "integrator.setConstraintTolerance(constraintTolerance)\n",
        "simulation = Simulation(topology, system, integrator, platform)\n",
        "simulation.context.setPositions(positions)\n",
        "if inpcrd.boxVectors is not None:\n",
        "    simulation.context.setPeriodicBoxVectors(*inpcrd.boxVectors)\n",
        "\n",
        "# Minimize and Equilibrate\n",
        "\n",
        "print('Performing energy minimization...')\n",
        "simulation.minimizeEnergy()\n",
        "print('Equilibrating...')\n",
        "simulation.context.setVelocitiesToTemperature(temperature)\n",
        "simulation.step(equilibrationSteps)\n",
        "\n",
        "# Simulate\n",
        "\n",
        "print('Simulating...')\n",
        "simulation.reporters.append(dcdReporter)\n",
        "simulation.reporters.append(dataReporter)\n",
        "simulation.reporters.append(checkpointReporter)\n",
        "simulation.currentStep = 0\n",
        "simulation.step(steps)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "398dfaba",
      "metadata": {
        "id": "398dfaba"
      },
      "source": [
        "If you could access AI Studio, you could consider to fork my project there https://aistudio.baidu.com/aistudio/projectdetail/2694076?contributionType=1 , it provide free Tesla V100 GPU at the time when this tutorial is written."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "61dccdda",
      "metadata": {
        "id": "61dccdda"
      },
      "source": [
        "# Below this cell was added in 2022."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "cb860867",
      "metadata": {
        "id": "cb860867"
      },
      "source": [
        "## MMPBSA to estimate the binding energy of protein and small molecule."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "3aaaa963",
      "metadata": {
        "id": "3aaaa963"
      },
      "outputs": [],
      "source": [
        "ante-MMPBSA.py  -p complex.prmtop -c com.prmtop -r rec.prmtop -l ligand.prmtop -s :WAT:Na+:Cl-:Mg+:K+ -n :LIG --radii mbondi2"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "84e121b5",
      "metadata": {
        "id": "84e121b5"
      },
      "outputs": [],
      "source": [
        "f = open(\"mmpbsa.in\", \"w\")\n",
        "f.write(\"\"\"&general \"\"\"  + \"\\n\"\n",
        "\"\"\"  endframe=\"\"\" + \"\"\",  interval=\"\"\" + \"\"\", strip_mask=:WAT:Na+:Cl-:Mg+:K+, \"\"\" + \"\\n\"\n",
        "\"\"\"/ \"\"\" + \"\\n\"\n",
        "\"\"\"&gb \"\"\" + \"\\n\"\n",
        "\"\"\" igb=\"\"\" +  \"\"\", saltcon=\"\"\" +  \"\"\", \"\"\" + \"\\n\"\n",
        "\"\"\"/\"\"\")\n",
        "f.close()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "552413ee",
      "metadata": {
        "id": "552413ee"
      },
      "source": [
        "Please use Pymol to visulize how many dcd frames you got, then change the mmpbsa.in accordingly. mmpbsa.in examples can be found https://ambermd.org/doc12/Amber21.pdf"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "ebda2cee",
      "metadata": {
        "id": "ebda2cee"
      },
      "outputs": [],
      "source": [
        "MMPBSA.py -O -i mmpbsa.in -o FINAL_RESULTS_MMPBSA_decomposition.dat -sp  complex.prmtop -cp com.prmtop -rp rec.prmtop -lp ligand.prmtop -y prot_lig_prod_all.dcd"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7d269ad9",
      "metadata": {
        "id": "7d269ad9"
      },
      "source": [
        "# The above mmpbsa.in only runs GBSA calculation, for the sake of time."
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.8.8"
    },
    "colab": {
      "provenance": []
    }
  },
  "nbformat": 4,
  "nbformat_minor": 5
}