{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Contact analysis: number of contacts within a cutoff\n", "\n", "We calculate the number of salt bridges in an enzyme as it transitions from a closed to an open conformation.\n", "\n", "**Last executed:** Feb 06, 2020 with MDAnalysis 0.20.2-dev0\n", "\n", "**Last updated:** January 2020\n", "\n", "**Minimum version of MDAnalysis:** 0.17.0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n", " \n", "**Optional packages for visualisation:**\n", "\n", "* [matplotlib](https://matplotlib.org)\n", "* [pandas](https://pandas.pydata.org)\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import PSF, DCD\n", "from MDAnalysis.analysis import contacts\n", "\n", "import numpy as np\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files\n", "\n", "The test files we will be working with here feature adenylate kinase (AdK), a phosophotransferase enzyme. (Beckstein *et al.*, 2009) The trajectory ``DCD`` samples a transition from a closed to an open conformation." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "u = mda.Universe(PSF, DCD)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Defining the groups for contact analysis" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We define salt bridges as contacts between NH/NZ in ARG/LYS and OE\\*/OD\\* in ASP/GLU. You may not want to use this definition for real work." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "sel_basic = \"(resname ARG LYS) and (name NH* NZ)\"\n", "sel_acidic = \"(resname ASP GLU) and (name OE* OD*)\"\n", "acidic = u.select_atoms(sel_acidic)\n", "basic = u.select_atoms(sel_basic)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Calculating number of contacts within a cutoff" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Below, we define a function that calculates the number of contacts between `group_a` and `group_b` within the `radius` cutoff, for each frame in a trajectory." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "def contacts_within_cutoff(u, group_a, group_b, radius=4.5):\n", " timeseries = []\n", " for ts in u.trajectory:\n", " # calculate distances between group_a and group_b\n", " dist = contacts.distance_array(group_a.positions, group_b.positions)\n", " # determine which distances <= radius\n", " n_contacts = contacts.contact_matrix(dist, radius).sum()\n", " timeseries.append([ts.frame, n_contacts])\n", " return np.array(timeseries)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The results are returned as a numpy array. The first column is the frame, and the second is the number of contacts present in that frame." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(98, 2)" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ca = contacts_within_cutoff(u, acidic, basic, radius=4.5)\n", "ca.shape" ] }, { "cell_type": "code", "execution_count": 6, "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", " \n", " \n", " \n", "
Frame# Contacts
0069
1173
2277
3377
4485
\n", "
" ], "text/plain": [ " Frame # Contacts\n", "0 0 69\n", "1 1 73\n", "2 2 77\n", "3 3 77\n", "4 4 85" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ca_df = pd.DataFrame(ca, columns=['Frame', \n", " '# Contacts'])\n", "ca_df.head()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Plotting" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, '# salt bridges')" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOy9eXxcd3nv//7OrpE0ska2ZFuy7DjO5lhKnD2EkoUAAcIatgRo6O+ylG53aVnCbSlw6Y/ScoHLcltS9hazhBC2BgKklIaU2LEdx3LsrI4tS7YlWSPNjGY0+/f+cc4Zzb5IM1qs5/16+WXPmXNmvj46+pznPKvSWiMIgiCsHmxLvQBBEARhcRHhFwRBWGWI8AuCIKwyRPgFQRBWGSL8giAIqwzHUi+gFtauXau3bNmy1MsQBEFYUezbt++M1npd4fYVIfxbtmxh7969S70MQRCEFYVS6nip7eLqEQRBWGWI8AuCIKwyRPgFQRBWGSvCxy8IwsojmUwyMjJCLBZb6qWc9Xg8Hvr6+nA6nTXtL8IvCEJTGBkZob29nS1btqCUWurlnLVorZmcnGRkZIRzzjmnpmPE1SMIQlOIxWJ0dXWJ6DcZpRRdXV11PVmJ8AuC0DRE9BeHes+zCL+w5DzwxGlOBWeXehmCsGoQ4ReWlGA0yXv+eR/femR4qZcinMXcdddd/PrXv+aHP/whn/jEJ8ru981vfpMdO3YwMDDAzp07+dSnPjWv7ztw4AD333//fJfLZz/7WaLR6LyPr4YIv7CkDI0GAZiKJpZ4JcLZzO7du7nmmmv4zW9+w4te9KKS+/zsZz/js5/9LL/4xS8YGhrikUceoaOjY17fJ8IvCBWwhD84m1zilQhnI+973/sYHBzk0Ucf5dprr+XLX/4y733ve/nYxz5WtO8nPvEJPvWpT7Fx40YA3G4373rXuwBDyK+55hoGBwd53etex9TUFAA33HADH/jAB7jqqqs4//zzeeihh0gkEnz4wx/mu9/9Lpdeeinf/e532bNnD9deey07d+7kBS94AU899RQA6XSav/iLv2DHjh0MDg7y+c9/ns997nOcPHmSG2+8kRtvvJF0Os073vGO7JPIZz7zmQWfF0nnFJaUodFpQIT/bOejP3mCwydDDf3M7Rt9/PWrLq64z9///d/zpje9iW9+85t8+tOf5oYbbuDhhx8uue+hQ4e4/PLLS773+7//+3z+85/n+uuv58Mf/jAf/ehH+exnPwtAKpViz5493H///Xz0ox/lV7/6FR/72MfYu3cvX/jCFwAIhUI89NBDOBwOfvWrX/GhD32Ie++9l7vvvptjx45x4MABHA4HgUAAv9/Ppz/9aX7961+zdu1a9u3bx+joKIcOHQJgenp6vqcsiwi/sKSIxS80m/3793PJJZfw5JNPctFFF9V9fDAYZHp6muuvvx6AO++8kze+8Y3Z91//+tcDcPnll3Ps2LGyn3HnnXfyzDPPoJQimTSu91/96lf84R/+IQ6HIcV+v7/o2K1bt3L06FH+9E//lFe+8pW89KUvrfv/UIgIv7BkTEUSnAgY2Twi/Gc31SzzZnDgwAHe8Y53MDIywtq1a4lGo2itufTSS/nd735HS0tL3v4XX3wx+/bt46abbqrre9xuNwB2u51UKlVyn7/6q7/ixhtv5L777uPYsWPccMMNNX9+Z2cnjz/+OA888AD/+I//yPe+9z2++tWv1rXGQsTHLywZh04a1v7Wda0i/ELDufTSSzlw4ADnn38+hw8f5qabbuKBBx7gwIEDRaIPRubP+973Pk6fPg1AIpHgy1/+Mh0dHXR2dvLQQw8B8M///M9Z678c7e3thMPh7OtgMEhvby8AX//617PbX/KSl/ClL30pe8MIBAJFx585c4ZMJsNtt93Gxz/+cfbv3z/PMzKHCL+wZFhunhduW0toNkkmo5d4RcLZxsTEBJ2dndhsNp588km2b99edt9XvOIV/Mmf/Ak333wzF198MZdddhmhkBGX+MY3vpENFB84cIAPf/jDFb/3xhtv5PDhw9ng7vvf/37uuusudu7cmfdU8M53vpP+/n4GBwe55JJL2LVrFwDvfve7ueWWW7jxxhsZHR3lhhtu4NJLL+Vtb3tbxXTUWlFaL/9ftiuuuEIv9SCW4ckof/Lt/XzlzitZ1+5e0rU0mngqzZ1f3cOfvfg8XnDu2pqPG52e5d3f3MsX77iMLWtb6/7e9/7LPp44GeLt12zmb+4/wsGPvBSfp7YmU8Ly58iRI/PyqQvzo9T5Vkrt01pfUbivWPw18sjRSQ6OBHlseGqpl9JwjpwK88jRAL97brKu4x4/Mc0TJ0N8/T+Pzet7h0aDDPR10OE1xD4YFXePICwGIvw1MhyI5v19NjE0YqSHjYfidR03HjKaQt27f4TZRLquY6ciCUamZhno7aCjxRR+8fMLwqIgwl8jluCfOBuF3/S1j4Xr65s+FjZuFOFYip8ePDmv7xwU4T+rWQmu5LOBes+zCH+NHDcF//hZKfxGAGusTot/LBRjQ4eHc9e18u099fXasYT/YhH+sxaPx8Pk5KSIf5Ox+vF7PJ6aj5E8/ho5cZa6emLJNE+PGWljE3Va/BPhOD0+D7cObuDj/3qEJ0+HuHC9r6Zjh0aCbO7y0tHiJBI3shxE+M8u+vr6GBkZYWJiYqmXctZjTeCqlaYKv1LqvwLvAhTwT1rrzyql/MB3gS3AMeBNWutlHTENx5IEIglcDhsjgVkyGY3Ndnb0GT9yKkQ6o7log48jp0Ik0xmc9toeBMdCMc5Z28ptl/Xxdw88xa7dw3zsNTtqOnZoNMjO/jUAYvGfpTidzponQgmLS9NcPUqpHRiifxVwCXCrUmob8EHgQa31ecCD5utljWXlX7XFTyKd4XTo7Jkharlcbr6oGzCs+FoZCxkWf2eri1fsWM99+0eJJkpXLuYSiCQYnTYCuwBelx2HTYnwC8Ii0Uwf/0XAbq11VGudAn4DvB54DfANc59vAK9t4hoaguXmuW6bkeO+nNw9mYxmaCQ47+OHRoJ0tbq4pM+wvsdK3NQOjQaLiqtiyTTB2STdZk3DHVdvJhxP8dPHT5X8nqfHwvz80Gl+fug033rkOAADfYbwK6XoaHGK8JchEk/x7PjMUi9DOItopvAfAn5PKdWllPICrwA2AT1aa0sdTgM9pQ5WSr1bKbVXKbV3qX2EltC/cBkK/8+fOM2rvvBbnjodrr5zCYZGg+zo7aDHZwSGxgss/qMTM9z6+d/y8ydO5223ngy6zeOu3NLJed1tfKtEkDcST3HbP/wnf/gv+/jDf9nH//7l07gdNnb0zvU67/CK8Jfjq799ntd+8WEJkgoNo2k+fq31EaXUJ4FfABHgAJAu2EcrpUpezVrru4G7wajcbdY6a2E4EGWN18mFG9qx29SySuncf9wIjzw7PsMF69vrOjaWTPPM+Awv2d5Dj8+w3McLLH7L0nxmbAYG5rZbTwbWDUMpxe1X9fOxnx7miZNBLt44J+o/efwk4ViKz92+k23r2gDwt7ryqnQ7WpyERPhLMjo9y0w8RTSRptUt+RjCwmlqOqfW+ita68u11i8CpoCngTGl1AYA8+/xZq6hERyfjNLv9+K029i4xsPxyeUj/AdNH/18nkIOm4HdHb0ddLW5sanilM5yhWvWft057Stuu6wPt8NWlNr57T3DnN/TxqsGN7B9o4/tG32s78hPPetocTItlbslmYwY08nCserxE0GohaYKv1Kq2/y7H8O/vwv4MXCnucudwI+auYZGcCIQZZPfC0C/37tsXD2ZjM4Ot5jPmg6ZN42B3g7sNsW6djfjBSmdJ8oUrln7WRY/GO6aVw5u4IePncymaB4aDfL4SJDbr+pHqfKZUOLjL89UVvjl/AiNodkFXPcqpQ4DPwH+WGs9Dfwt8BKl1DPAzebrZUs6oxmZmqU/K/yty8bV8/xkhBlTYOezpoMjQda2udhgWt/d7Z66LH6nXdHpzW+qdsdV/czE5yp5d+0Zxu2w8fqdlXOMRfjLE7CEPy4Wv9AYmuow1Fr/Xoltk8CLm/m9jeTk9CypjGZzjsU/GUkwE0/RtsT+Viub5/yeNo4HInUff8gM7FqWeI/PzcjUbN4+VqXy6VCMWDKNx2kHjFhAd7unyIq/fHMn5/e0sWv3MK8c3MiPHhvl1sGN2UZs5ehocRKKJc+qGolGEYiKq0doLNKyoQqWJd2fI/xgtGleaoZGg7gdNm66sIeT0zGS6UzNx84mjIrdwZzMmm6fJy+rJ5PRjARmWW+6c0am5v7PY+EY3b7i9tRKKe64qp/HR4J84v4jRBJp7ri6v+p6OlqcaC1WbSGpdCYb+xBXj9AoRPirYLk4LB//5i5v3valZGg0yPaNPrauayWd0Zyarr2w7PCpEBlNXkpld7ubQCRBImXcQMbCMRLpDC88rziNdTwUp6e9dG+Q15lB3m/tHubC9e1cZlboVsJnVu9KZk8+0znnY0YsfqFBSG5YFYYDURw2lfWDWzeA4Xm4ViqRTGd4//cPcmZmzuK+/vx1vPP3tpbcP5PRPDEa5LbL++aeQgJR+s0bE8BTp8N88udPlnwSsPLwrSIqmAvUTszE6V3Tkn2qeeG2tXx/30jeU85YKMa153aVXFtHi5NbBzdy7/4R7ri6clDXYk1O24ZNVfeujy8/dJTfPD1XC7K2zc0nbxvE5Vj+do/l3wdx9QiNY/lf+UvM8UCUvs4WHGb/mo4WJx0tzoZb/MOBKPc9NsqJQJSZeIojp0J87eFjZfc/eiZCJJFmR29HVvgL/fw/PDDKb56eYCaeKvrT4rLzhsv7sm4cIJvLb+XoW/79nf1r8Lrs2dexZJpQLJWX0VPIe284l1cObOB1O3tr+v83s1/P1//zGIdPhpiJpxgLxbjvsVGGRqcb/j3NIF/45WlIaAxi8VchN5XTwkjpnC1zxPyIxo3atv/5yu28ZHsPf/uzJ/nKb4+itS5pMVupmIN9RtWty24ruhkdGg1y4fp27vuj62paQ7fpurEGspwIRLHbFBvXtNDv92bjHeMlcvgL2dbdxhffellN3wtkg7/NyOUPzia57bI+PvLqizkdjHHNJx5kaCTI5Zv9Df+uRpMn/BL/EBqEWPxVGA5Esxa1RX+Xt+EpnVZaZqvbyJrp8blJpjVTZYTw4EgQj9PGtnVt2G2KPn9L3pq01hwcCTKY48qphhWstXL0hwNRNq7x4LTb8uoXrIEt3RUs/npplsWfzmjCsVT283t8bta1u7OFb8sdS/hddpu4eoSGIcJfgeBskuloslj4/V5GpqKkM43rJGEVPFkpopb1XappGhjW/PYNvqwLqt/vzasoHpmaJTibzAveVqOr1Y3dpuZcPZPRvGym4UAUrXVOu4bGDZ1vlvBb7hHr85VSDPR2ZJ+YljuW8Pd1toirR2gYIvwVsCzozV3Fwp9Ma04FG+fuiSQsi98Q/kJ/ey7pjObQyWC2rbG1puHJaLaR18GRuarcWrHbFOva3Hmunqzwd3mJJTNMhOPZ98tl9cyHFqcdp73xrZmtz7OEH4xMpmfHZ2pqIb3UBCIJ2t0O/K2u7FOhICwUEf4KFKZyWmz2Nz6lc6bA4i/XLRPg+TMzRM3ArkW/30s4nsoK3dBoEKdd1d24rdvnZiwcZyaeYjKSoN/fmv18MP7PY+EYLruNNVWKsuqhWa2ZrZhBrvAP9naQ0WTbXSxnApEE/jYXbR6HuHqEhiHCX4Fywr+pCUVckXi+xb+uvXS3TMgZVN43lx/fX3AzOjQa5IL17bgd9rrW0d3uYTwUK1+4FogyHoqzrt1dU5pmPfia0KEza/Hn3KSsFNahFeDuCUQS+FtdtHucIvxCwxDhr8DxySidXmde+2CADR0eHDbVUIs/Ymb1eM2WCB6nnY4WZ8kB6FZg99x1rdltVv7+cdPdMzQaZKC3euFUIT0+N2OhWDZeYAl+b2cLShmfPxaKNdS/b7FmARZ/Mp0hVaJeoZSrp8fnYV27e0EDbBaLQCSB3+uifYksfq018VS6+o7CikKEvwIjU8WpnAAOu43ezhZOTDXQxx9P4XXZ8/rU9PiKu2UCPHEyxEU5gV2ATZ1zFvmJgBHYrce/b9Hd7mEqmuS5CaMPvyX8boedDT4PJwJRxs0h641mIa6e//adA/zZdx4r2l5K+MFw96woi9/tWJLg7o8OnOSa//9BYkkR/7MJEf4KTEUTdLW6Sr63rs3N5Ezt82mrEUmkioZslOqWCXB8MsK55kATi1a3g7Vtbk4Eohw0i5PmI/yWJb/v+JRRrJbjIunvMjJ7xkKxijn886Wjxcn0bKL6jgVorfnts2eMYTEFlBP+Hb0dPDsxk3WxLUe01gSilqvHQTyVybbTWCyeHgszFU3mVZQLKx8R/goEZ5NFgmHR2erKK65ZKDPxdFG3z26fu8jHH0umGQvFi1JMAfr9LQwHogyNBnHZbZy/vq1on2pYlvzeY4GSaaxPj4UJx1INzeG36GhxEpxHAZf1hFPq5xGaTeJ22LJdRS0G+zrQ2uhZtFyJJNIkUpmsjx8Wv3rXOqeNvNaFpUeEvwLT0fLC39Xqyk5GagSReCpbvGXRY3bLzB10Xi7FFOZy+YdG5hfYhbkirlAsVVL4Q6afuVmunnA8VTTYvRqWy2Yqmig6ttzP0HoaWs5+fmsAi7/VlTUKFjulU4T/7ESEvwyFFZ+FdLa6mIokGjYAeyaeotWVb/H3tLtJZTRT0blfunKZRmAI86ngLEMjwboKt3LpzsnNL5fNZOzXeFePz2rNXGcQ03JtZXRxAVi5p7Zun4fudvey9vNP5gh/u8e4NhY7wCvCf3Yiwl8G65HaV8HiT2V01gJeKJESg10sd0qun3+4IM0yl/6uVjJmT/t6WjXk0tXqwm4GmAufKjZ3zWURNcvih/qrd3OrcAPRfIGq5K4b7FveAd5AxPi557p6Qovt6omK8J+NiPCXwRKfNd7Swd1Oc/tUg34hool0UXA3W72bk9kzHIjS6rKXDDrn3gzmE9gFsNlU1pov5eopXFsjsc51PcKvtWZoJMgmfwtQLFCVhH9HbwfPTcws24rYQMQ4D7kW/2L35BeL/+xEhL8M5bJBLPxthkg1ys8/U8LHb7ldJnIt/kkjxbRU8ZQlzC67jfN76qvYzf/e0sLf6XXS5nbgctjKnpeFMB+LfzgQJRRLcf3564D6hH+g1wzwLtMK3nyLf/FdPal0JvuzEOE/uxDhL0M14bcs7lIW/48OjPKJnx3J/tm1e7jq90VK+Pit6t3cfj2luoVadLe7cTlsXLihfUFDRrp9nrzhMxZKKTb5vXQ3oWoX5s61ldKZyWj++ZHjFVMJrZ5E15/fDRQLVGg2WXber/VUdHCkeb35Hz0W4LHhqXkdG4gkcdoVbW5HTVk91vmq5Qnm4WfPsL/KuqZnk1ghrKUQ/qlIgm/tPt6wOJowh/TjL0M14bdcPYW/EFprPnDvQZJpjcOmSGc0qYzm5u3deYHTXDIZXdLV43HaWeN1Zl09WmuGA9GsdVuIzaZ46faeebt5LF64bS12pfIKxCxuvqg7L9jcSAot/oefO8Nf/fAQz4yF+dhrdpQ85pCZunrNVqO3fu7PI5XOEI6XD9B3+zz4W13ZYrVm8NGfPIHX6eB7f3ht3ccGInH8rS6UUjVl9ewfnuKvfngIl13x5isrzzn+qx8eIpXR/Ptf3FB2uH3uuVwK4f/cvz3D1x4+xnXnrmXL2tbqBwg1IxZ/Gapa/GVcPaHZFLFkhrtefiFPffzl7HrXNQAV2wBbnTkLg7tgdMC0umGOh+PEU5mSqZwWX7jjMt5z/bll36+FO1+whX98++Ul3/vzl17Ax187sKDPL0eh8FtPSvftHy3bSXNoNMiFG9pp9zjxuux5AmUF3iu5pXp8c+e3GZwOxjldprV2NQKRJP5W46nP5bDhdlTuyW89/VRrJZLOaE5MRRkORHn4uTMVvn8uq6gwaN5sYsk0P9g/CjDv8yeUR4S/DKW6OubS4rTjdtiKrN/xgiElF2/0odTcL2UprD49hRa/8TlGt0yonMp5NuBx2nDZbQRnk4yHY/zy8BhXbukkHE/x08dPFe0/15PIeMLxFxTVVbt5g9UWoznCn0pnmIzEGQ/H5uWuMCz+ubW3e5wVs8isDKVq0+FOh2Ik08Z6KrkhrXO5bV3bolv89w+dyv78mvXzWc2I8JchNJvEZbfhcZY+RUopo4hrJv8XYizbq96w1FrdDs5d11bR4i+cvpWL1S0T5rqBlvPxr3SUUtkOnffsHSGV0XzytkG2dbfxrT3FAnV8Mko4lsoKf9d8hL/dU3bYzUI5M5NAa4glM/NK+52Kzln8AO0eR0VXT1b4JyNl9wGj5QcYRskvD4+V7AcFc8J/bncbwdlkySZ4zWLX7uFsjKlUh1phYYjwlyE4m8TX4qwYxOxsdRVZ/JaI5LY0GOjtqGjxRyu5enxuJszq3eFAFKWgr/PsFH6AjhYHU5Ek33l0mGu3drF1XRt3XNXP4yemeeJk/jm0xidaxWqFbTRqtfjPzMSbImq5N5T5iNfkTBy/N9fiL9+obSae4rmJGZSq7uqxqr/ff8uFpDKae/aOlNwvK/zrWtHaCPYuBk+Phdl7fIo/uG4LboetaTfm1YwIfxmCs8mqg0b8Jdo2WI+luZWtO3o7GA/Hy/7yW1ac11VK+D2kMkazruFAlI0dLQvK2FnurPG6eOiZCU4EZrn9aiNAedtlfbgdtiK3xKHRIC7HXOrqfFw963weMrpxabm55Loo6nVXJNPGU0KhxV/Ox3/4ZAit4cotfqaiyYqFXsOBKA6b4rpzu7h2axffeXS4ZJuMQCRBu8eRLdZbLHfPrt3DuOw23nD5pmzbEqGxnL0KskAq5X9b+M22DbmMhWK0ux15/vrBKoM/LB9/KYu/OyelczgQzRYqna10tDiJJNL4W1287OIeY5vXySsHN/CjAyfzumkeHJnmovVzqat+bxnhr3AD78kOvGm8uORaqvVardaTpFUvAsb1Ua6Ay7q2XjmwAZiz6ksxHJilt7MFh93GHVf3cyIwy0PPFgd5AxGjO62VurwYwm8EdUe4Zcd6/K2u7HwIobGI8JehVuEv/GUYD8dYV1DVun1D5QBvpJKP3xrBGIpXzOE/W7DO+Rsv78trMvfWq/uZiaf4yeMnASMF9onRUF5PIn+bi9lkmtmEcSMN1eTqqTzUfiGMh2JYnsJS7bUrkc2oyakcN6Zwlbbkh0amWe/zcPnmTqDydLjhyUj2OnrZxevpanWxa/fxkmvobHXRuYjC/9ODpwjFUtxhPu11NznrarWy6oQ/OJvkVZ//LbuPTlbdr5rwd5kDsHMnFI2H4kVDyFvdDrZVCPAWztvNxWqNcHwywkQ4ntcv52zEOudvuSo/D/2y/k7O72njL394iIG/foDBj/6CcDyVV7NgiaSVejgdTeBx2ip2Ke0u0RajXr7462e56wdDRdvHQnHWtrlpczvKBlDLkZtKaVHJ1TM0ajTmsyaxVfLzG0+OZpW3w8YbLu/jV0fGi4S9URb/mZk4L/3MbyomOFh8b+8Jtq5r5epzjLqM7nax+JvBqhP+H+wfYWg0yP7hytWawQotmS06s9W7c1bYWLj0WMKBChOfCuft5mJV7z563KiyPFtTOS3efu1mPvPmSzinoGBHKcUnXj/A71+7hTdesYk3XbGJ995wLi83XRswJ5IBM9Oqlpv32jY3Si3M1fPDx0b52aHidNNx81ow5irM0+LPFX63g5lEcdvqmXiKo2ciDPZ14PM4WeN1lhX+UCzJVDSZ9+T4wvPWks5onh4LF62h0+vK9lCar/AfHJnm6bEZ/umhoxX301pz5GSIF25bm02q6PF5iCTSy7af0kplVVXuaq2zAUKrD0op0hlNOJ4q25nTwrKEJiNx1nd40FozFio9lnBHbwc/eGzUnFeb/34knkIp8LqKLVO3w06n18neYwHg7E3ltDh3XVvRdDGLyzf7uXyzv+yxWeGP1i78TruNrtbSIy5rIWJm01gtoXO/bywUZ0OHh2giXb+Pv6TFb7StjiRS2RYOAE+MBtF6rgVFv99bVviz8xxyrqPNfuMmOxyIcs3WLiBn+lebC5fDRrvHMW/ht9xOPxs6zUdelcgaTIVMR5OE46mSzQDHQzHaylwXQv2sKot/3/Epnhk3yvMDkfJZD7X4hiG3Q2fSPC5FIpXJWum5ZAO8Jfz8kUSaVpejbOpoj29uBOPms1z4F0JW+M2bei3CD9aA+flZ/IdPhbAM8MKA6ng4ZvT9n0eRmJVltKYgnROKG7UNFaS1VhJ+S4Rznxw3rPFgt6m89Wenf5nXeKl4Vq0MB2ax2xSJdIZ795dOHTX2K65Tsdym8/35CKVZVcK/a/cw7W4H565rrWjx15IGCLltG4zPsvzEpSz+7Rt92FTpzB5r0Ho5rBtJu9tRNcV0NdNlpj5aN/XgbPk+PbksxI+ceyPPFc5kOsNkJEF3u9u8cddXvTsVSdDR4sSZ0y+prYzwHxoNsqHDk71O+v1eRqdmS9YmZMU1p+2H025j4xoPx3MCwoVPHAsT/gjndbdxWf8adu0ZLnseSq3NisHM94lMKM2qEf7paIKfDp3itTt76e30Eqgw2zXbi79mi9/4hbDEo5Twe10OtnW3lRT+mRJDWHKxPq9cO2bBoN3jwG5T2Zt6aDZJR0tpt0IuC8kVPzQaxGcK8vEc4T8zE0dr47O7293EUxlCs7X7qSfNwGoulntnJp5/7R4czZ+41u/3kspoTgWLxXI4EKXT68Tnyb+2C58SJguF37sQ4Tey0e64ejNHJyI8cjRQdj9rLRa5WW1C42iq8Cul/rtS6gml1CGl1LeVUh6l1NeVUs8rpQ6Yfy5t5hos7t0/SiKV4far+s3S/hos/irW9RqvC6Xmgl7Zdg1lhpTsMAO8hRaPMW+3kvCX7o8v5GOzKTpzBGo6mqjN4vd55l29e3A0yJVb/PhbXXnCmXstZMWrDqt1KlrsC7dcPbntH8KxJM+fieRlN1nXSalc/nIpwf3+1rz9G2XxWx1l+/1ebh3cgM/jYFeJ9htguKHWtrnzChnb3YeCPlcAACAASURBVA5anHbJ7GkwTRN+pVQv8GfAFVrrHYAdeIv59vu01peafw40aw0WRlD3OJduWsP2jT5DHGbKX8S1unrsNsWaFmfWOso2aCvTfnmwt4OJcLzIXxmJp0vm8FtYn1epK6dg4G91EogkSKYzRBLpml09Whu9derBCuwO9HWwye/NE85s6452T7ZIrB4/9eRMIi+wC4YIQr6r5wmzYncgZ9Sm5So5Xkb4y81rnowkstkz1jVtuc/8bUaHznqbzU3MxIklM/R3efE47bz+sj5+fugUkyVmLBg3iPwCRaVUXqNCoTE029XjAFqUUg7AC5xs8vdleWx4ih8/fpIfP36Su//jKM9NRLJFIV1tLiKJNLFkuuSx0zUKP5jVu2YWyXgoTrvHQUsZf/1AmQre6q4e45fvbE/lbASWZToXoK+euNYzD4scckS3t4N+vzfPR265jnp87pqKxGYTae4fOpW9ZsdCsbziLchx9eQIv5Ubn2vxb+howWFTRQHeVDrD6NRsGYvfzP83/w/WE3Gn2R3U73WRSBk303oobCz41qv7SaY1399XHOQt9zTSk9OosByHRoPZwj2hOk1L59RajyqlPgUMA7PAL7TWv1BK3QH8jVLqw8CDwAe11kW3c6XUu4F3A/T3Vx4qUUgmo3nz3Y+QSM09uvtbXbxqcCOQ45uPJtjQUdwCodasHjAsIqtDZ6lUzVwu2uAD4MlTIV6yvSe7PZKo7OrZ1t2O3aYWPGBlNeBvdfHU6XDN7jrImW1cpx95KEd09w9Pcf/QKVLpDA67jfFQDJuCrjY3bWaBX6Uisf/778/y+X97Nm/b5rX5IjiX1TPn4x8yA7tr2+bci3aboq+zpUj4TwVjpDK65JNjVvgDUbZv9BGIGN1pLYPE3zoXz6pkpBRS6Lc/r6ednf1r+NehU3lzIxKpDKeCs/R39RV9RrfPzRMVxmNGEyle938f5s9uOo8/ffF5Na9tNdM04VdKdQKvAc4BpoF7lFJvA+4CTgMu4G7gA8DHCo/XWt9tvs8VV1xR1/NlNGmkor3n+q288fJNgJFzb1ni/pxKxFLCH5xN4nLY8DjLu18sOludPH/GaHM7Ho7nNWcrxOty0JkzUcvCcPVUEv42Dnz4JXm520JpLIu/VncdzL9tw6HRYNaHv9nfStoMqG7yexkLxVjX7sZuU3hdDtrdjrIBymQ6w3cePcHvnbeWv37VxQDYFGwpqNL2uuzYVL6rZ2gkWNIg6O9qLfLxW69Lunq68uMCgUiczta57rT+bM1Koq4nT6ujbG/n3O/ZVVv8fO3hYyRSmWyfpdHpWTK6dByrx+fh354cL/sdp4PGfIHHTjRvhObZRjNdPTcDz2utJ7TWSeAHwAu01qe0QRz4GnBVo7/YqoTt93vZ1t3Gtu62vECZv0oJei1Vu3Of5c4J7la2+MHqr1/o40/RWiGdExDRrxF/q5vp2WTVQTq5dLW6zOrd+oT/4Mh0VnQ35VjMYBkBc9eCkctf+vMfPDLGRDjOndduyV6vW9e1FY1EtEYwWn74cCzJ0YLArkW/vyXP9QRzPv9S4trR4qSjxcnxgGHE5E7/guIaiVoZDkTZ4PPktc0Y6Osgkc7kVQqXyuix6G53E61QvWu51UolTgilaabwDwPXKKW8yjAbXgwcUUptADC3vRY41OgvrtT7BmoQ/tlk1VTOuc9yMhVNksloxkPxbN5xOQoDVemMZjZZ2eIXasfvNapbj5nDRmpJ53TYbaxtq6/IymqTMNC7BsgJqJpia1Rwz10L3e2esq6kb5lDR264oPQs5VyMKVzGTe3QqOH+yA3sWvT7vQRnkwRz0paHA1GcdlXyKdc6xpreVTj9a+53pr6e/MOTxcFk60aVG+uyhseUckNVeyKztpdKnBBK0zTh11rvBr4P7AeGzO+6G/iWUmrI3LYW+Hijvzvb+6ZEf3ugatOpWis+wbAw0+aQlEQ6U9SgrRBjxuvcBVxp3q5QP37T122532r9Odbb/jfbJqHPiNus93lw2ucCquOhWN4wnp4yFv+JQJSHnjnDm6/cVHK4fSG5jdpKBXYt+s02DCem5qz+4UCUvk4v9jLD1ftzMpMKp38txOIvFPN+vxefx5HXrXY4EMXtsLGurdhwyjbSK/PzyX2CLtcPS8inqVk9Wuu/1lpfqLXeobV+u9Y6rrW+SWs9YG57m9Z6ptHfO1Oh6RkYYmBTFPXSt6hP+I39njxtWF/VLP7ciVpQuUGbUD9WJky9wl/JIi9FYZsEu02xqdMQzkRqrmrXwmq7UeiK+PaeYWwK3nzlppq+t90z15P/4GiQ3jUtdJUQy/4C1xMYN5lK/vn+Li8jU1HSGV00/avN7cBpV3VZ/LOJNOPheJH7RinFQF9HXrdOK8200L0Fc+nM5WIkY6EYLoetbGW8UExV4VdKXaeUajX//Tal1KeVUpubv7T5U2mwCcwV+pSbulSvxQ9w5JThr6zFx29N1DLWKsLfSCzL9PkzEVqc9pqnlZWzyMtxaDTIep8nz4+/yax+PTNjpXLm+vg9JFKZbNAZjKDu9/aOcNOFPWXdL4W0e5yE45arJ8iOXl/J/QpdT9a/C/Pk847xe0mmNScC0aLpX0opM3Be+83RetoodbPZ0dvBk6dD2Zbmw4HZsn2oeqq0bRgPG83wtnW3MTQiAd5aqOW34h+AqFLqEuDPgeeAbzZ1VQuk0mATi8L5rLlY83ZrwbIwj5wyLP7qrp78x9aZ7E2qegaRUB1L+EenZ2u+eYNxQ540C79qobBNAmDm8kdyWnfk+vgt8ZoTzl8eHuPMTJw7rq7N2gfDmAnHUoTMit3BvjVl9+vKqSYORpMEZ5PZTpylsCzzx03xzJ3+BfmJDLVQmMOfy2DvGpJpzdOnZ4zq3slI2aeRNrcDr8te9olsLBSjp93DQO8ahkZDEuCtgVqEP6WNM/ka4Ata6y8C7c1d1sKo5uqB8iXoqXSGmXhtzb1g7pfjydOGxV89uJv/2BqtEo8Q6sMqONK6djcPGNa5Ub1b3aINZ0W3WPhDsRTPjBney9yngVIBym/vGaZ3TQvXn99d8zotH/+hAldTKTb5vfz80Cne/KXf8ftf3Z3dVg5LoA+YaZGFBWRWVXQpYsk0f/69x7MGEMy5mUoND8oN8AYiCSKJdNmWJEqpbKO7UoyHjaSKgV4fZ2binJb2DlWpRfjDSqm7gLcD/6qUsgHLOrcwmqhB+Ms0nbL6oNQs/OYvx3Agis/jqJr7P2f5WRa/uHoaidthz7Y2qM/ir72I6/6hU2gN123ryttuuVf2HjeakOUaAYVFYscnI9mgbrlgaynaPU5mcoS/UlHfW6/uzw6i9zjt3HRhN1edU36ewYYODw6bmhP+1tot/p8fOs29+0f4Qk4R2nAgSpvbqF0pZJO/hY4WJ0Oj0xVTOS3WtZfPuhoPxehu9zBgPv2Uan0u5FOL2rwZuAP4/7TWp5VS/cDfN3dZC8Nyn3griLC/zcXe48UXcT2FPwAtLjstTjuzyXRV/z7MtVi2BCBSw01KqI/OVhfheKqmql2LbNuGGqzFXbuHOb+njcv6O/O2W8L16LEp7DaV7XMDOQFK84b/7T0nsNtUzUFdi3aPg0Q6w95jU/SuaSkS51zeeMUm3nhF7Z/vsNvo7WzJVskWCb+3vMVvDTh64InTTITjrGt3ZwO2pTrKKqWyU+ms4S+VelH1+Dwl/fcz8RSRRJoen5vtG+Zan7/04vW1/adXKVUtfq31aeBewLqKzwD3NXNRC8UqiCqVIWDh97qy+fe5ZFsy1yEa1i9ILcLvdtjxt7qKfPyV4hFCfVg/j/pcPdbs3coW/6HRII+PBLnjqv4iQbPcKM+fibCuzZ1nybe47LR7jOrdRCrD9/ed4MUXdtd0zeRitW145OhkU1p49Pu92VYnpSz+UCxVFAd5ZizMnmMB3nRFH6mM5p59J4DSTddy2dHbwVOnwzxnDkfq66wg/O3ukllRua3QW1x2zutul8yeGqglq+ddGPn4XzI39QI/bOaiFkq1NsdgXNTpjM4Ww1jUa/FbnwVUbNeQS3f73MSnSJViM6F+5iP8XW1ubDVU7+7aM4zbYeN1O4t7yrS5Haw1Yz6lYj2Wn/oXh09zZibB7VfX14MK8lszlyrcWii57pZC48dKXbaaElrs2jOM0654/y0XcvU5fr6z5wTpjJEdVMq/bzHY10EyrfnF4TG6291lmxuCcT5nk2nCBdW7VqzM+t2z0kQlwFuZWnz8fwxcB4QAtNbPALVHo5aAat0uIXd6Vv5FvCDhr9F66/Z5mDAf+SPxFDYFLTX0BRJqYz7Cb7cp1rZVLuKaiaf40WOj3Dq4sawbybL6S7Xm7jb91Lt2G0HdF51XvVK3kDb33Pc2y+IHiqZ/wVzqcq67J5ZMc+++EV528XrWtrm54+p+hgNR7ntslHgqUzGYbK3/ydPhqi3He3ylc/mzrdDN9wd6Ozgzkyg5hEaYoxbhj2utsz9ps8Xysr6dRuIpvFVcJ4XTsyyCpjUzH+EvN4ClkJ4ci38mnqo4b1eon/kIP1SfxPXjAyeJJNLZ9t6lsISz1LXQ4/Pw5KkQ//ncJLdfVV9Q18Ky+KG5wl8qdmBlTOXOsvjXg6cIxVLZc3LLjvV0ep185pdP531eKfo6W7JPFdUav80VceULemHqbLnW50I+tQj/b5RSH8Loq/8S4B7gJ81d1sKIxNNV0yNzuw3mYln8tebx535Wrf7aHp+HiZk46YwmWqUzp1A/8xf+ykPXd+05zoXr27msv3TuPOQKfwmL3+cmkkjjsCneVEfQNRdL+Ps6W4omdDWCTRWEPzvTOMfVs2vPMFvXtnKtGaB1O+y84fI+RqeNnj+VhN8K8FbbD3JjMIXCH6fFac8+4W/f4MNuU5LZU4VahP+DwARGb533APcDf9nMRS2USKK6q6dco7bgbBJ3jS2ZCz+rZh+/z+jvE4gkmElUfzoR6mO+wr+u3cPp4CzPjs8U/fnV4TEOjYa44+rioG4u/VlXT4meM6bVevNFPTW7BQtpN109zZrNYKWklhJ+a9vTY8Y5+Y+nJ9h3fIrbCwLdt19lWP82Bb1rKlck76hR+MvN3h0PG83wrO/3OO2cV2a2NUA8lc77uZ40b1BLTTxVfjBUM6hqamqtM8A/mX9WBLUGd6G08NcrGBvXeIp6jlfCEoCxUIxIDfEIoT42mu0PqhXTFdLX2cJUNMnNn/5NyfdbnHZeu7O34mds624DSrsuNpnXx1uvqT+oa7Gm1YnDpthZ4aljIfg8Tnp8bjZ0FN+YOr1OXHYbn3vwGT734DMAuBw2brs8P9C9dV0bLzi3i5PTs1VbZlgpsdZ5K0eb20Gb21Ek1GMFzfDAuCk++OQ4Wuuim/Rd9w7xg8dG87bd90cvYGdBau5i8xf3HGQqkuBf3nn1onxfVcUxO2kW+vSDwF7g41rryWYsbCHM1OA+8TjteF32hgj/rYMb2bq2reZ+K7m9RyKmj19oHNdt6+Le976AizfWZxX//rWb2dLVSrpMRsg5Xa34qsxF2Nnfyb3vfUFJd9CLL+rh3vdey+WbyxdRVcPncXLfH13HeT2VhXIh7HrXNUVVu2Dk+X/3PddwYmpOfDf7vSWfDv7PW3bmTQorx80XdXPve68t23oilwvWt3P4VP4krvFQrKh6eaCvg3v2jXAyGCt64tj9fICrzvHztms2MxGO879+epjhQHRJhV9rze+eO5NNlS0MqjeDWhTnZ0Aa2GW+fgvG/NzTwNeBVzVlZQvAsKKru0/8ra7i4O5ssq4cfgCn3cYlm2q3wLqz5ftxZuJpetc03le7mlFKcfnm+n+R2z1OXjm4YcHfX+677Ta1ING3aEYaZy7nrit/U9nZ31mTSK5rd2eLFSth/KxqOycDvR1891EjVdRuU2itTVdPscUPRgVvrvAHIglGp2e58wWbefUlGwlEEvyvnx6uq/9QMzgdinHGDJg/PRau22CZD7XcWm7WWt+ltR4y//xP4Hqt9SeBLc1dXv3UM9jE31rcoTM4W3ufnvli9RwfD8VrvkkJwmpnoLeD2WSaoxNGwddMPEXUrNrN5SIrwDuaX+lb2Erbas++1MKfO5dgsYLStQi/XSmVHY+olLoSsJSq9Cy0JSTbAqEG94m/1VVUjBKqozPnfHE5bHS1uhizXD3i4xeEqhSmao5li7fyLf65AG++W6iwsZ3dplhTpmfXYnJoNIjdpmh12RctDbUW4X8n8BWl1PNKqWPAV4B3mT36P9HMxc2Hevrb+1tdTM7k/9Cno4mmW/xguHvGQ7GaMpAEQTBcUC1Oe9ZCtnL6SwXxB/s6GBqZzqvgPTgyzTlr8+M05br0LiZDo0HO625jsG/N8hF+rfWjWusB4FLgEq31oNZ6j9Y6orX+XvOXWB+ROnrfGP165n7oyXSGSCK9OMLf7ubkdIxYMoNXgruCUBW7TXHxRl/WcreK7UrVTAz0djAVTWbrCcCYUVwYCC7XpXex0FozNBJkoLeDgb4OnjwVzvZKaiZlFUcp9T/KbAdAa/3pJq1pQdTT+8bf5iKaMPJnPU47oXm0a5gvPT43jx4z2vdKgzZBqI0dOQFeq2q3VM2E1aL50GiQvk4vkzPxbGA3F3+ri+cmGj79tWZOBWNMRhIM9HXQ6XWRSGd4eixccc5CI6ikjtawlQuAK4Efm69fBexp5qIWQl2uHu9c9W7vmpZ59emZLz0+D9FE5RGRgiDkM9jXwdf/8xjPTcwwForjddlL/v5cuL4dh01xcCTILTs2ZF0oA7352XedrS4Cx5bO4rfcVgO9Hdk2MkMF091+fug0u5+fy5rv9Lr44xu3zavlh0VZxdFafxRAKfUfwGVa67D5+iPAv877G5vMTD0Wf+tcv57eNS3sOz4FlJ4Y1GhyrRQJ7gpCbeSmao6HY/T4PCUrqT1OO+f3zLVotrJlLi6YUdxlJnhkMrpiG/dmYQV2L9rgw+2w0e5xMDQa5Hbz/Ug8xZ9/7wCpjMblsBltXhJprtu2dl4pyxa1BHd7gNxbYsLctiyxsnq8FVq8WhT269m1Z5hz17VW7MXSKHKrDcXiF4Ta2LquDa+Z/TIeilesFRjonWvRPDQaZOva4gI8f6uLjJ7r0bXYWIFdj9M+N5wmJ6Xzx48bjQF3vesahj7yMh788+uBuQyl+VKL8H8T2KOU+ohp7e/GKNxalswNL6/P4j9yKsRjw9NFfUeaRW5ASix+QagNK8A7NBpkzLT4yzHQZwR4R6Zmi9wnFtnWLdHFd/dYN6Tc2c0D5nAaK8C7a/dwXmPA9T4Pa9vcC87+qSWr52+APwCmzD9/oLVedmmcFvWmc4Jh8e/aPWz0HbmseMBGM8gtOpHgriDUzo7eDp44GWQsFKOnisUP8O9PT3AqGCvZ2K5cz67F4GQwRiCSyFvXQF9HNsA7NBI03D45xqjxVOBbcKFXpawen9Y6pJTyA8fMP9Z7fq11YEHf3CSi8RRK1ebq8Xmc2G2K0alZfvjYKK/Ysb4prW5LsbbNjVKgtbh6BKEeBvs6+NrDhkVcqRHfBWaA15oHXKrVxVIKvyXeAzl9iqybwEFT9D1OW1FjwIHeDn7z9ASziXTFqWWVqKQ4u4BbgX3kN2lT5uut8/rGJjNj9uKvxV1jsyk6vU7u3T9COJ7ijqs3Vz2mUTjtRvXumZmE5PELQh3kWsiVXD0ep50L1rdnh8dfvNFXtM+SCv/oNA6b4sL17dlt/X4vPo+DR45O8uCRMV41uLEoy3Cgbw0ZDYdPBefd+6msq0drfasy1PN6rfXWnD/naK2XpeiD1ZK5vl76wdkk27rbuHLL4nbos0rNxeIXhNo5Z21b9om+1IjLXKybxNa1rbSX6Ky6tMIf4rye9rzZH0opBvo6+MnB8tPecjOb5ktFH7826p2XbepmKWYS9fW+sX7wdyxSUDeXbp8x4NvjbH4bVkE4W7DbFDvMDpbVZi5Y7p1yHU3LtWdvNkbF7jSDJeIOA71r0NpoNndpia6/PT6j82lhL6J6qEUh9yulrtRaPzrvb1lE6u1v39XmxuWw8frLKg/YaAbrfR7aPU6ZtysIdTLQ18GeY4Gq404HzYKtShPLFqtfz5d+8xz/25xFjIZEOsOOEjckK8vnjqs2ldSGbNpnQffReqhFIa8G3qqUOg5EMH38WuvBeX9rE6nX1fPHN2zjtst6WVNi8ESzec/153LzRcu2JEIQli3v/L1zGOzrqOom3dHr429et4NbBzeW3WexhP+BJ06z3ufhFQPGzAeXw8arS6zrxRd185FXbeeNFeYy7+jt4N+fGieaSM0rRljLES+r+1OXEGOwSe3zTLdv9LG9RNBnMThnbSvnrG1+lbAgnG1s6GjhNZdWf0pXSvHWKkkbpbr0NppUOsPhUyHuuGozH3z5hRX3dTvsvOO6cyruM9jbYQR4T4a4Ykv9Ad5a8viPA13Aa4BXA13mtmVJtE4fvyAIq5vF6ND57MQMsWSGgb7GGJmFswnqparwK6U+DHwDQ/zXAl9TSv3lvL5tEZDBJoIg1MNiuHqyOfu9jWkH0+Pz0N3unndmTy0K+VaMPvwxAKXU3wIHgI/P6xubzExcBpsIglA7/jYXs8n0ggqiqjE0GqTVZWdrA127RoC3SRY/cBLIdZq7gdFaPlwp9d+VUk8opQ4ppb6tlPIopc5RSu1WSj2rlPquUqphUdVUOmMONpEWCIIg1IbVnr2Z/XqGRoNc3NvR0A6gA30dPDcxk21TUw9lhV8p9Xml1OeAIPCEUurrSqmvAYeAqnlESqle4M+AK7TWOzDm9L4F+CTwGa31NozeP/+l7lWXISL97QVBqJPcZo3NIJXOcPhkqGJK6XwYsAK8p+rP56+kkHvNv/cB9+Vs//c6P79FKZUEvMAp4CbgDvP9bwAfAf6hjs8sSz0N2gRBEKC4PXs1IvEUvzoyRjJtdLKxKbjhgu7s5xTyzPgM8VQmrwtnI8it4L2yRGbP9/aeKHtspUEs31jIorTWo0qpTwHDwCzwC4ybyLTW2no2GQFK5mQppd4NvBugv7+4bLkUIvyCINTLXNuGeE3737P3BB/5yeG8be9+0VY+9IqLSu5vBWAbPU6x2+dhbZuLp06Hi96Lp9K8//sHyx7btF4BSqlOjBTQc4CNQCtwS63Ha63v1lpfobW+Yt26dTUdM+fqER+/IAi1MSf8tQ1jOXomQrvbwUPvv5GH3n8jW9e1cnQiUnb/odEgbW4H5zRhst+GjhbGwrGi7RPhyjexZprGNwPPa60nAJRSPwCuA9YopRym1d9HjYHiWsha/NLtUhCEGrHas9dq8Q8HovR3ednk9wJGA7gTgWjZ/Q+OBrl4o68pox17fG5Gp4uFf7yK8NeSx//GWraVYBi4RinlNbt8vhg4DPwaeIO5z53Aj2r4rJqYEVePIAh1YrRnd9Vs8Q8Homzu8mZf9/tbGQ5EMXpa5pNMZzhyqvGBXYtun4fxUAnhL7Etl1pcPXfVuC0PrfVu4PvAfmDI/K67gQ8A/0Mp9SxGUdhXalhDTYiPXxCE+eBvddZk8aczmpHAbNbaB+j3tzCbTHOmRNuHZ8ZmSKQyZbuDLpTudjeTkUR2VKPFWGierh6l1MuBVwC9ZlqnhQ+oKXFUa/3XwF8XbD4KXFXL8fUyJ/zi4xcEoXb8rS6marD4x0IxEukM/bnCb1r/w4FI0fB3q4Nmsyx+qzvpmZk4G9e05K3TXsG1VMniP4mRhRMz/7b+/Jhl2ritnkHrgiAIFv5WF5M1WPzDpi8/T/j9rXnv5WIFdrc0IbALc7O7xwpcO+PhON0V5hFXSud8HHhcKfUvOemXy5pIPIVNQYtTLH5BEGrH3+piKlrd4rfEfbN/Tsj7Og1Le3hytmj/oZEgO3qbE9iFuQlkha6dsVCM7gqzCipV7g4ppQ5iDGI5WPinMctuLJFEquZ5u4IgCBZ+r4upaIJ0xgjQHjsT4Y++tY9wLP9mMDwZxW5TbMhp/e5x2lnv8xRZ/Ml0hiOnw01z88DcBLKJgpTO8dA8LX6MQesrCunMKQjCfPC3utAapqMJutrcfOk/jnL/0Glev7OPm7fPDUsaDkTZuMaD055vM/f7vQwH8nP5nzodJpHKNLxwK5euVjd2myqy+MfDMa48p/wM8UqunmXbc78ckXgarwR2BUGoE3+bYR1PRRO4nXZ+fMAoLxoaDRYJf65/36K/y8tDz0zkbTtkds4c7GtMK+ZS2G2KtW2uPB9/PJVmKpqkp8Ig+lry+K9RSj2qlJpRSiWUUmml1Pyn/DYRacksCMJ8sDp0Ts4k+PGBk0QSaVpd9qx4WxjCXxyo7fd7GQvFiSXT2W1Do0HaPQ42l7hRNJIenyevYGvctP4rDaKvJY//C8DtwDNAC/BO4IsLWGfTqHfQuiAIAuR06Iwm2LXnOBeub+dlF6/nYI7wh2NJApFEaYvf3DYyNefnHxoNsmNjY1sxl6K73ZNn8Vs3gXkFd3PRWj8L2LXWaa3116ij585iMiM+fkEQ5oEl/L95eoJDoyHuuLqfgb4OJsLxrKieCBhZO+VcPQDHJw3hT6QyPHkq3LTCrVx6fO4Ci99YbyVXTy0qGTWHpRxQSv0dRmvlpjV3WwiRREoatAmCUDedrU4Avr9vhBanndfu7OVps+vl0EiQnu2ebPC2ksVvZfY8PRYmkc40NaPHorvdQ8Cs3nU5bNkbVc8CXT1vN/f7EyACbAJuW/hyG080nhaLXxCEunE77LS5HSTTmlddsgGfx8n2jT5siqy7J1u81VUs/F2tLrwue3YfayTiYgi/JfATM4bVPxaO4zD7D5WjqkrmZPfEzNYNm0zXz7JDXD2CIMwXf6uLmXiKO67eDIDX5WBbd1s2wDsciNLR4qSjxVl0rFKKfr8326UzG9gtcZNoNFbbhrFQjN41Ldkc/kqxhVqyev5d8b7oAgAADjlJREFUKeVTSvkxGq79k1Lq041adKNIpTPEUxkJ7gqCMC82rvGwo9fHJTl++R29HRwcCaK15vhk6VROi36/N+vjHxoJMtDbsSjFpFZ/ICubZzxcuWoXanP1dGitQ8DrgW9qra/G6LW/rIiYfXqkQZsgCPPhM2++lK/eeWWeWA/2dnBmJs5YKM6JMjn8FkYRV5R4Ks1TpxcnsAtzFv+4Wb07FopVrNqF2oTfoZTaALwJ+OnCltg8ZhJGOyHJ4xcEYT5s6GgpspQt8T5wYpqRqdmS/n2L/i4v8VSGh589s2iBXTDiC0b1riH84+F49mZQjlqE/2PAA8CzWutHlVJbMXL6lxXSi18QhEazfUMHNgW/OHyaVEZXtfgBfnrwFLA4gV0wBsmsa3MzbhaQTUeTFTN6oLbg7j3APTmvj7IMs3qs6Vti8QuC0ChaXHbO627nl4fHgNKpnBbWe788PIbP46i4b6Pp8bkZC8ezs3Yb4eNfEYjFLwhCM9jR20E4ZuhLJTHv7WxBKQjHUgz0LU5g18IawWi5exrh418RWMFdr0uCu4IgNI5B08/vsCk2dJS3pN0OOxs7jN78A73Na8xWiu52o3rXquBthI9/RRARV48gCE3Aaqvc29mCw15ZMjf5LeFfHP++RY/PqN616giqCX9VlVRK/aXW+uPmv91a6+rzyRaJe/ae4MAJY6blM2MzgLh6BEFoLNs3+LDbVE0++36/l0eOBrJPCYuFFcx94mQIp13R6S0uMsul0rD1DwD/AbwB+Li5+XfAZQ1ZaQP42E8Pk0xnslb+JZvWVP0PC4Ig1EOLy87rd/YyuKm6++amC3uYCMez4xgXC2sE46HRIN3tnqrxhUrm8ZPAG4GtSqmHzNddSqkLtNZPNWi98yaaSBGOpXj/LRfwRzdsW+rlCIJwFvP3b7ykpv1u2bGeW3asb/JqirF67x89E2Fnf/UbVCWH1TTwIeBZ4Abg/5jbP6iU+s8FrbIBZIcNVGg9KgiCsBrI9elXasdsUcnifxnwYeBc4NPAQSCitf6DhS2xMcxFryunLQmCIJzt+L0uHDZFKqMrTt6yKGvxa60/pLV+MXAM+GfADqxTSv1WKfWTRi14vsz1nBaLXxCE1Y3NprLN2mrRxFrSOR/QWu/VWt8NjGitXwgsudVfa6GCIAjCasCq1q1FE6sKv9b6/Tkv32FuOzO/pTWOiXAcl8NWsje2IAjCaqPHFPxq7RqgzgIurfXj81tS4xkLxejxuRe1LFoQBGG5Yvn2a4l7rtjK3bFQXDJ6BEEQTDaY7SIWmtWzrBkPx7hgfftSL0MQBGFZ8JYrN7Glq5XO1vKzdi1WrMU/Lha/IAhClq42N68c3FDTvitS+CPxFOF4qqZ8VUEQBCGfFSn82eItsfgFQRDqZmUKvxRvCYIgzJumBXeVUhcA383ZtBWjBcQa4F3AhLn9Q1rr++v57DFp1yAIgjBvmib8ZgfPSwGUUnZgFLgPo+r3M1rrT833s8ezVbti8QuCINTLYrl6Xgw8p7U+3ogPGw/HcTts+FpWbDaqIAjCkrFYwv8W4Ns5r/9EKXVQKfVVpVRnqQOUUu9WSu1VSu2dmJjIe8+o2q0+bEAQBEEopunCr5RyAa8G7jE3/QNGq+dLgVPA/y51nNb6bq31FVrrK9atW5f33lgoJs3ZBEEQ5sliWPwvB/ZrrccAtNZjWuu01joD/BNwVb0fOB6OS0aPIAjCPFkM4b+dHDePUiq3tOx1wKF6P3A8FJfiLUEQhHnS1OioUqoVeAnwnpzNf6eUuhTQGENe3lPi0LJE4ilm4inJ6BEEQZgnTRV+rXUE6CrY9vaFfKaMXBQEQVgYK65yV0YuCoIgLIwVK/yS1SMIgjA/VpzwT5iunlrGiwmCIAjFrDjhHwvF8Dht+DxStSsIgjAfVqDwx6VqVxAEYQGsOOEfD0vVriAIwkJYecIfiot/XxAEYQGsOOEfC8Vk8pYgCMICWFHCPxNPEUmkpV2DIAjCAlhRwj83clGEXxAEYb6sKOEfC8mQdUEQhIWyooR/PGxW7YrFLwiCMG9WlvCHpGpXEARhoawo4R8LxWhx2ml3S9WuIAjCfFlZwh+O0+NzS9WuIAjCAlhRwj8eiskAFkEQhAWysoQ/LCMXBUEQFsqKEX6ttVG1K4FdQRCEBbFihH8mniKaSEuDNkEQhAWyYoR/btauWPyCIAgLYcUIf3bkovj4BUEQFsSKEX6reEssfkEQhIWxcoQ/LEPWBUEQGsGKEf6xUByvy06bVO0KgiAsiBUk/DGZtSsIgtAAVozwj4firBM3jyAIwoJZOcIfluItQRCERrBihH8sFKdHLH5BEIQFsyKEP601s0mZtSsIgtAIVoTwp9IakBx+QRCERrAihD+ZzgBIS2ZBEIQGsCKEP2UKf4+4egRBEBbMihD+pOnqkVm7giAIC2dFCH8qo2mVql1BEISG0DThV0pdoJQ6kPMnpJT6b0opv1Lql0qpZ8y/O6t9VjKdkcCuIAhCg2ia8Gutn9JaX6q1vhS4HIgC9wEfBB7UWp8HPGi+rkgqnZGqXUEQhAaxWK6eFwPPaa2PA68BvmFu/wbw2moHJ9NaLH5BEIQGsVjC/xbg2+a/e7TWp8x/nwZ6Sh2glHq3UmqvUmqv4eoRi18QBKERNF34lVIu4NXAPYXvaa01oEsdp7W+W2t9hdb6Co0UbwmCIDSKxbD4Xw7s11qPma/HlFIbAMy/x2v5EPHxC4IgNIbFEP7bmXPzAPwYuNP8953Aj2r5ELH4BUEQGkNThV8p1Qq8BPhBzua/BV6ilHoGuNl8XRURfkEQhMbQ1IoorXUE6CrYNomR5VMXMmtXEAShMayIyl2bUrRK1a4gCEJDWBHC77TLnF1BEIRGsSKE32FbEcsUBEFYEawIRW33iJtHEAShUawI4ZccfkEQhMaxIoRfEARBaBwi/IIgCKsMEX5BEIRVhgi/IAjCKkOEXxAEYZUhwi8IgrDKEOEXBEFYZYjwC4IgrDKUMQRreaOUCgNPLfU6lhlrgTNLvYhlhpyTYuScFLOazslmrfW6wo0rpRfCU1rrK5Z6EcsJpdReOSf5yDkpRs5JMXJOxNUjCIKw6hDhFwRBWGWsFOG/e6kXsAyRc1KMnJNi5JwUs+rPyYoI7gqCIAiNY6VY/IIgCEKDEOEXBEFYZSxr4VdK3aKUekop9axS6oNLvZ6lQCm1SSn1a6XUYaXUE0qp/2pu9yulfqmUesb8u3Op17rYKKXsSqnHlFI/NV+fo5TabV4v31VKuZZ6jYuJUmqNUur7SqknlVJHlFLXrvbrRCn1383fm0NKqW8rpTyr/TqBZSz8Sik78EXg5cB24Hal1PalXdWSkAL+XGu9HbgG+GPzPHwQeFBrfR7woPl6tfFfgSM5rz8JfEZrvQ2YAv7Lkqxq6fg/wM+11hcCl2Ccm1V7nSileoE/A67QWu8A7MBbkOtk+Qo/cBXwrNb6qNY6AXwHeM0Sr2nR0Vqf0lrvN/8dxvhl7sU4F98wd/sG8NqlWeHSoJTqA14JfNl8rYCbgO+bu6yqc6KU6gBeBHwFQGud0FpPs8qvE4wi1RallAPwAqdYxdeJxXIW/l7gRM7rEXPbqkUptQXYCewGerTWp8y3TgM9S7SspeKzwPuBjPm6C5jWWqfM16vtejkHmAC+Zrq/vqyUamUVXyda61HgU8AwhuAHgX2s7usEWN7CL+SglGoD7gX+m9Y6lPueNnJyV01erlLqVmBca71vqdeyjHAAlwH/oLXeCUQocOuswuukE+OJ5xxgI9AK3LKki1omLGfhHwU25bzuM7etOpRSTgzR/5bW+gfm5jGl1Abz/Q3A+FKtbwm4Dni1UuoYhgvwJgz/9hrzkR5W3/UyAoxorXebr7+PcSNYzdfJzcDzWusJrXUS+AHGtbOarxNgeQv/o8B5ZgTehRGU+fESr2nRMX3XXwGOaK0/nfPWj4E7zX/fCfxosde2VGit79Ja92mtt2BcF/+mtX4r8GvgDeZuq+2cnAZOKKUuMDe9GDjMKr5OMFw81yilvObvkXVOVu11YrGsK3eVUq/A8OXaga9qrf9miZe06CilXgg8BAwx58/+EIaf/3tAP3AceJPWOrAki1xClFI3AH+htb5VKbUV4wnADzwGvE1rHV/K9S0mSqlLMYLdLuAo8AcYxt2qvU6UUh8F3oyRHfcY8E4Mn/6qvU5gmQu/IAiC0HiWs6tHEARBaAIi/IIgCKsMEf7/194dq1YRhGEYfj9QCBjQQisL09iKQbEUBLHUiKUg8QLsvINILkMszA3YicROkBSCHiIqQryCFKJYBPNb7IQcJMoJ5ETDvE+3u7MwW+zP8O/yjSR1xsIvSZ2x8EtSZ47KZuvS1CT5yfC77I6Fqvryj6YjTZ2/c6p7Sb5V1exfrh8by3aRjjxbPdIekiwmeZbkJbCaZDbJapI3SUZJbrVxcy3//kmST0lWklxP8qpl4F9p404keZxkrYWodZc0q/+HK35177dWz0ZV3U6yCDwCLlTV5k6sb1V9TXIaeA2cB84BnxlSU9cZokbeMmS83wTuV9VCkmXgfVU9TXIKWAPmq+r74T2pNLDHL8GPqrq4x/kXY/EGAZaTXGWIzjjLbsTxRlWNAJKsM2x8UklGwFwbc4MhWO5hO55hiFEY30hGOhQWfunPxlfjd4EzwKWq2mrJoDPt2njOy/bY8Ta771iAO1X1cXrTlSZjj1+azEmGPQC2klxjaPHsx3PgQUuJJMn8QU9QmpSFX5rMCnC5tW/uAR/2ef8ScBx419pBSwc8P2liftyVpM644pekzlj4JakzFn5J6oyFX5I6Y+GXpM5Y+CWpMxZ+SerMLwoQRdw54S+jAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "ca_df.plot(x='Frame')\n", "plt.ylabel('# salt bridges')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[1] Oliver Beckstein, Elizabeth J. Denning, Juan R. Perilla, and Thomas B. Woolf.\n", "Zipping and Unzipping of Adenylate Kinase: Atomistic Insights into the Ensemble of OpenClosed Transitions.\n", "Journal of Molecular Biology, 394(1):160–176, November 2009.\n", "00107.\n", "URL: https://linkinghub.elsevier.com/retrieve/pii/S0022283609011164, doi:10.1016/j.jmb.2009.09.009.\n", "\n", "[2] Richard J. Gowers, Max Linke, Jonathan Barnoud, Tyler J. E. Reddy, Manuel N. Melo, Sean L. Seyler, Jan Domański, David L. Dotson, Sébastien Buchoux, Ian M. Kenney, and Oliver Beckstein.\n", "MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations.\n", "Proceedings of the 15th Python in Science Conference, pages 98–105, 2016.\n", "00152.\n", "URL: https://conference.scipy.org/proceedings/scipy2016/oliver_beckstein.html, doi:10.25080/Majora-629e541a-00e.\n", "\n", "[3] Naveen Michaud-Agrawal, Elizabeth J. Denning, Thomas B. Woolf, and Oliver Beckstein.\n", "MDAnalysis: A toolkit for the analysis of molecular dynamics simulations.\n", "Journal of Computational Chemistry, 32(10):2319–2327, July 2011.\n", "00778.\n", "URL: http://doi.wiley.com/10.1002/jcc.21787, doi:10.1002/jcc.21787." ] } ], "metadata": { "kernelspec": { "display_name": "Python (mda0170)", "language": "python", "name": "mda0170" }, "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.6.7" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": true } }, "nbformat": 4, "nbformat_minor": 2 }