{ "cells": [ { "cell_type": "markdown", "id": "0e9a6ccd", "metadata": {}, "source": [ "# Introducing FRET restraints" ] }, { "cell_type": "code", "execution_count": 1, "id": "ea263007", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from matplotlib import pyplot as plt\n", "import pandas as pd\n", "import mdtraj as md\n", "import fretraj as ft" ] }, { "cell_type": "code", "execution_count": 2, "id": "b6bd6397", "metadata": {}, "outputs": [], "source": [ "import importlib" ] }, { "cell_type": "markdown", "id": "bedd9c07", "metadata": {}, "source": [ "## Accessible contact volume\n", "\n", "Calculate ACVs for donor and acceptor dye on an **solvated and energy minimized** biomolecule.\n", "> Note: the generated .gro file must first be converted to a .pdb file using `gmx trjconv`" ] }, { "cell_type": "code", "execution_count": 3, "id": "0cdbec73", "metadata": {}, "outputs": [], "source": [ "labels = ft.cloud.labeling_params('gmx/mn_riboswitch_labels.json', verbose=False)" ] }, { "cell_type": "code", "execution_count": 4, "id": "f1af7a93", "metadata": {}, "outputs": [], "source": [ "struct = md.load_pdb('gmx/mn_riboswitch.pdb')\n", "struct_nosolv = struct.remove_solvent()" ] }, { "cell_type": "code", "execution_count": 5, "id": "acfec770", "metadata": {}, "outputs": [], "source": [ "don = ft.cloud.Volume(struct_nosolv, '1-B-U1-O5\\'', labels)\n", "acc = ft.cloud.Volume(struct_nosolv, '1-A-A1-O5\\'', labels)" ] }, { "cell_type": "markdown", "id": "bba17084", "metadata": {}, "source": [ "Compute a FRET trajectory from the ACVs" ] }, { "cell_type": "code", "execution_count": 6, "id": "2d5d1227", "metadata": {}, "outputs": [], "source": [ "fret_traj = ft.cloud.FRET(don, acc, 'A1-U1', labels)" ] }, { "cell_type": "code", "execution_count": 7, "id": "45cfeff6", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
R0 (A)<R_DA> (A)<E_DA><R_DA_E> (A)R_attach (A)R_mp (A)
value54.096.10.0491.986.694.4
stdNaN10.70.0313.0NaNNaN
\n", "
" ], "text/plain": [ " R0 (A) (A) (A) R_attach (A) R_mp (A)\n", "value 54.0 96.1 0.04 91.9 86.6 94.4\n", "std NaN 10.7 0.03 13.0 NaN NaN" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fret_traj.values" ] }, { "cell_type": "markdown", "id": "ea4c67c5", "metadata": {}, "source": [ "Create a Plumed object from the loaded structure and the ACVs. Look for neighboring phosphorus atoms within a range of 15 A." ] }, { "cell_type": "code", "execution_count": 8, "id": "6cd4f8ed", "metadata": {}, "outputs": [], "source": [ "plumed = ft.restraints.Plumed(struct_nosolv, [don, acc], selection='name P', cutoff=15)" ] }, { "cell_type": "markdown", "id": "bbc48092", "metadata": {}, "source": [ "Set your desired FRET value that you want the simulation converge to. The corresponding $\\langle R_{DA,E}\\rangle$ will be computed from the $\\langle E\\rangle$. The $\\langle R_{mp}\\rangle$ will be calculated from $\\langle E\\rangle$ using a third order polynomial (Kalinin et al., *Nat. Methods*, 2012)" ] }, { "cell_type": "code", "execution_count": 9, "id": "877919bd", "metadata": {}, "outputs": [], "source": [ "targetFRET = 0.6\n", "mean_R_DA_E = ft.fret.mean_dist_DA_fromFRET(None,None,targetFRET,54)\n", "targetRmp = ft.fret.R_DAE_to_Rmp(mean_R_DA_E)" ] }, { "cell_type": "markdown", "id": "c35bee9b", "metadata": {}, "source": [ "Compile a plumed input file for the MD simulation" ] }, { "cell_type": "code", "execution_count": 10, "id": "7d881844", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "True" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plumed.write_plumed('gmx/out/plumed_A1-U1.dat', targetRmp, 100, 100)" ] }, { "cell_type": "markdown", "id": "9365e9aa", "metadata": {}, "source": [ "To illustrate the restaints save a visualization script for VMD or PyMOL which can be run with either `vmd_vis -c filename.gro -x trajectory.xtc -v vis.vmd` or `pymol_vis -c filename.gro -x trajectory.xtc -v vis.py`" ] }, { "cell_type": "code", "execution_count": 11, "id": "9f178643", "metadata": {}, "outputs": [], "source": [ "plumed.write_vmd('gmx/out/vis.vmd')" ] }, { "cell_type": "code", "execution_count": 12, "id": "25f75832", "metadata": {}, "outputs": [], "source": [ "plumed.write_pymol('gmx/out/vis.py')" ] }, { "cell_type": "markdown", "id": "9b26ca64", "metadata": {}, "source": [ "Save force field and topology files files of the mean position pseudoatoms" ] }, { "cell_type": "code", "execution_count": 13, "id": "3e6ef4f8", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(1) Copy the \"MP.itp\" file into your force field directory\n", "\n", "(2) Add the following lines to your \"topology.top\" file:\n", "; Include topology for dye mean position\n", "#include \"/MP.itp\"\n", "\n", "#ifdef POSRES\n", "[ position_restraints ]\n", "1 1 1000 1000 1000\n", "#endif\n", "\n", "(3) Add the following line to your \"atomtypes.atp\" file in the force field directory:\n", "MP 0.00000 ;\n", "\n", "(4) Add the following lines to your \"ffnonbonded.itp\" file in the force field directory:\n", "; Dummy mass for dye mean position\n", "MP 35 79.90 0.0000 A 4.64693e-01 2.45414e-01\n" ] } ], "source": [ "plumed.write_pseudo('gmx/out/MP.pdb', 'gmx/out/MP.itp')" ] }, { "cell_type": "code", "execution_count": 14, "id": "549f64cb", "metadata": {}, "outputs": [], "source": [ "don.save_mp('gmx/out/mn_riboswitch_D.dat', units='nm', format='plain')\n", "acc.save_mp('gmx/out/mn_riboswitch_A.dat', units='nm', format='plain')" ] }, { "cell_type": "markdown", "id": "1175676b", "metadata": {}, "source": [ "First, backup the original .gro file. Then, insert the mean position coordinates given in **MP.pdb** at positions specified by **riboswitch_D.dat** and **riboswitch_A.dat** \n", "\n", "```sh\n", "mv em/mn_riboswitch.gro em/mn_riboswitch_original.gro\n", "gmx insert-molecules -ci MP.pdb -f mn_riboswitch_original.gro -o mn_riboswitch.gro -ip out/mn_riboswitch_D.dat -replace water |& tee tmp_insert.dat\n", "gmx insert-molecules -ci MP.pdb -f mn_riboswitch.gro -o mn_riboswitch.gro -ip out/mn_riboswitch_A.dat -replace water |& tee -a tmp_insert.dat\n", "```\n", "\n", "Run the following command will tell you how many solvent molecules have been replaced by the MP pseudoatoms and how to update the topology. \n", "\n", "```sh\n", "awk '$1 == \"Replaced\" {sum += $2}; END {print \"\\nIn your topology.top file under the section [ molecules ] do the following:\\n(1) decrease the number of solvent molecules by \" sum \"\\n(2) add the following line:\\nMP 2\"}' tmp_insert.dat && rm tmp_insert.dat\n", "```" ] } ], "metadata": { "jupytext": { "formats": "ipynb,md:myst", "text_representation": { "extension": ".md", "format_name": "myst" } }, "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.11.11" }, "source_map": [ 11, 15, 23, 25, 32, 36, 41, 44, 48, 52, 54, 58, 60, 64, 68, 72, 75, 79, 83, 85, 89, 93, 96 ] }, "nbformat": 4, "nbformat_minor": 5 }