{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Write your own contacts analysis method\n", "\n", "The `contacts.Contacts` class has been designed to be extensible for your own analysis. Here we demonstrate how to define a new method to use to determine the fraction of native contacts.\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 molecular and data visualisation:**\n", " \n", "* [matplotlib](https://matplotlib.org)\n", "* [pandas](https://pandas.pydata.org)" ] }, { "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 salt bridges" ] }, { "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": [ "## Define your own function" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Any function you define *must* have `r` and `r0` as its first and second arguments respectively, even if you don't necessarily use them:\n", "\n", " - `r`: an array of distances between atoms at the current time\n", " - `r0`: an array of distances between atoms in the reference\n", "\n", "You can then define following arguments as keyword arguments.\n", "\n", "In the function below, we calculate the fraction of native contacts that are less than `radius`, but greater than `min_radius`." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "def fraction_contacts_between(r, r0, radius=3.4, min_radius=2.5):\n", " is_in_contact = (r < radius) & (r > min_radius) # array of bools\n", " fraction = is_in_contact.sum()/r.size\n", " return fraction" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then we pass `fraction_contacts_between` to the `contacts.Contacts` class. Keyword arguments for our custom function must be in the `kwargs` dictionary. Even though we define a `radius` keyword in my custom analysis function, it is *not* automatically passed from `contacts.Contacts`. We have to make sure that it is in `kwargs`. " ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "ca = contacts.Contacts(u, \n", " selection=(sel_acidic, sel_basic),\n", " refgroup=(acidic, basic),\n", " method=fraction_contacts_between,\n", " radius=5.0,\n", " kwargs={'radius': 5.0,\n", " 'min_radius': 2.4}).run()" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "ca_df = pd.DataFrame(ca.timeseries, \n", " columns=['Frame', \n", " 'Contacts from first frame'])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Plotting" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "scrolled": false }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEGCAYAAABsLkJ6AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOy9eXxb9Znv/34kW/K+SMpuO4ltICtJIKsNhVJKgU4LbefOQFf665R2utz2dnpvYZYyl7a3nZl22jLTZZi5wO1KKfObDneaDlAILbEDSSCEkI3Y2eyslrxvkmV97x86R5Es2ZYcL1H8vF8vvSJ9zzlffc+xcj7n+zzP93nEGIOiKIoy+3DM9AAURVGUmUEFQFEUZZaiAqAoijJLUQFQFEWZpagAKIqizFJyZnoAmeDz+cySJUtmehiKoihZxSuvvOI3xswZ2Z5VArBkyRJ2794908NQFEXJKkTkRKp2NQEpiqLMUlQAFEVRZikqAIqiKLOUrPIBKMqlzNDQEK2trQwODs70UJRZSl5eHhUVFeTm5qa1vwqAokwSra2tFBcXs2TJEkRkpoejzDKMMQQCAVpbW1m6dGlax6RlAhKRR0TkvIi8Mcp2EZGHRKRJRF4XkWvitn1ERI5Yr4/EtV8rIvusYx4S/R+jZDmDg4N4vV69+Sszgojg9XozmoGm6wN4DLh1jO23AVdYr3uBH1gD8gAPAJuAjcADIlJuHfMD4ONxx43Vv6JkBXrzV2aSTH9/aQmAMeb3QPsYu9wB/MhEeQkoE5EFwDuAZ40x7caYDuBZ4FZrW4kx5iUTzUf9I+DO8cbRMxhOZ7iKoihKGkxWFNAioCXuc6vVNlZ7a4r2JETkXhHZLSK7z3b0TNJwFeXy5OzZs9x1113U1NRw7bXXcvvtt/Pmm29OqK/vfOc79Pf3T+jY1157ja1bt2Z0zEMPPcTy5cv5wAc+MKHvTJe2tjY2bdrEunXrePHFF7n99tvp7OxM+/jHHnuM06dPp9x26NAh1q5dy7p162hubp6sIU8Zl3wYqDHmYWPMemPM+iGcRCJawEZRUmGM4T3veQ833ngjzc3NvPLKK3z961/n3LlzE+pvugXg+9//Ps8++yw//elPE9rD4cmd+T/33HOsXr2aPXv2cP3117N161bKysoS9jHGEIlEUh4/lgD86le/4g//8A/Zs2cPNTU1afU3oxhj0noBS4A3Rtn2T8DdcZ8PAwuAu4F/Grmfte1QXHvCfqO9XPNrzb7WTqMolyIHDhyY0e9/7rnnzPXXX59yWyQSMV/84hfNypUrzapVq8zjjz9ujDFm27Zt5oYbbjDve9/7zFVXXWXe//73m0gkYr773e+a3Nxcs2rVKnPjjTcaY4z55Cc/aa699lqzYsUK8+UvfznW986dO82WLVvM1VdfbTZs2GA6OztNZWWl8fl8Zs2aNebxxx83L7zwglmzZo1Zs2aNWbt2renu7k4Y3yc+8YnY9/393/+9eeCBB8wHP/hBU1dXZ+666y4zMDBg7rnnHrNq1Sqzdu1a8/zzzxtjjHn00UfNHXfcYW6++WazePFi8w//8A/mW9/6llm7dq3ZtGmTCQQCCd+zZ8+ehLH19/ebxYsXm7a2NnPs2DFz5ZVXmg996ENmxYoV5vjx4+YjH/lI7Jr9/d//vfnlL39pCgsLzZVXXhk73ubXv/61mTdvnlm4cKG58cYbU/Y32jVcvHixue+++8yaNWvMtddea1555RVzyy23mOrqavODH/wgtt/f/u3fmvXr15vVq1cnHB9Pqt8hsNukunenaky549gC8E7gN4AAm4GdVrsHOAaUW69jgMfattPaV6xjbx9vDK75teaffteU8qQVZaaJ/4/310+9Yf7oh42T+vrrp94Y8/u/+93vms9//vMptz355JPm5ptvNuFw2Jw9e9ZUVlaa06dPm23btpmSkhLT0tJihoeHzebNm82LL75ojDGxG6ONfTMNh8PmhhtuMHv37jXBYNAsXbrU7Ny50xhjTFdXlxkaGjKPPvqo+fSnPx079g/+4A/M9u3bjTHG9PT0mKGhoaQxxn/fAw88YK655prYDfab3/ym+ehHP2qMMebgwYOmsrLSDAwMmEcffdTU1NSY7u5uc/78eVNSUhK7YX7+85833/72t5O+Z+TY4gVARMyOHTuMMcbs3r3b3HzzzbH9Ojo6jDHG3HDDDWbXrl0pr/MDDzxg/u7v/s4YY5L6G+0a2mP4/ve/Hxv36tWrY+c0d+5cY4wxTz/9tPn4xz9uIpGIGR4eNu985zvN7373u6QxZCIA6YaB/hzYAVwlIq0i8jER+aSIfNLaZStwFGgC/hn4lDW7aAe+AuyyXg9abVj7/It1TLMlAmPiznHQ0BRIZ8iKosSxfft27r77bpxOJ/PmzeOGG25g165dAGzcuJGKigocDgdr167l+PHjKft44oknuOaaa1i3bh379+/nwIEDHD58mAULFrBhwwYASkpKyMlJXl5UX1/PF77wBR566CE6OztT7jOSd7/73eTn58fG/8EPfhCAZcuWsXjx4phv461vfSvFxcXMmTOH0tJS3vWudwGwevXqUc9lNBYvXszmzZsBqK6u5ujRo3z2s5/lP//zPykpKcmor5H9QeprGH++9rg3bdoUOye3201nZyfPPPMMzzzzDOvWreOaa67h0KFDHDlyJOMxxZPWQjBjzN3jbDfAp0fZ9gjwSIr23cCqdL7fpsidw85j7YTCEVw5l7z7QpnFPPCuldP+nStXruTJJ5/M+Di32x1773Q6U9rcjx07xje/+U127dpFeXk599xzT0bx5vfddx/vfOc72bp1K/X19Tz99NMsW7ZszGMKCwszHr/D4Yh9djgcGfsP4r+zvLycvXv38vTTT/PDH/6QJ554gkceSbqVpd3feNcwftwjzykcDmOM4f777+cTn/hERmMYi6y6ixa5cxgYGua1lvQ99ooyW7jpppsIBoM8/PDDsbbXX3+dF198keuvv55f/OIXDA8P09bWxu9//3s2btw4Zn/FxcX09EQj77q7uyksLKS0tJRz587xm99EJ+xXXXUVZ86cic0menp6CIfDCccCNDc3s3r1ar70pS+xYcMGDh06lNG5XX/99THn8JtvvsnJkye56qqrMuojU/x+P5FIhPe973189atf5dVXXwVIOrd0Ge0apss73vEOHnnkEXp7ewE4deoU58+fz3gc8WRVKohCdw5hgYYmPxuXemZ6OIpySSEi/Nu//Ruf//zn+Zu/+Rvy8vJYsmQJ3/nOd7juuuvYsWMHa9asQUT427/9W+bPnz/mjfjee+/l1ltvZeHChWzbto1169axbNkyKisrqa+vB8DlcvGLX/yCz372swwMDJCfn89vf/tb3vrWt/KNb3yDtWvXcv/997N9+3a2bduGw+Fg5cqV3HbbbRmd26c+9Sn+9E//lNWrV5OTk8Njjz2W8JQ8FZw6dYqPfvSjseidr3/96wDcc889fPKTnyQ/P58dO3bEzFTjsWbNmpTXMF1uueUWDh48yJYtWwAoKiriJz/5CXPnzs2on3gkar3JDtavX28W3vMd3DkOfvnJupkejqIkcPDgQZYvXz7Tw1BmOal+hyLyijFm/ch9s8oEBFBX42PPyU76groqWFEU5WLIOgGor/USjhh2HhsrM4WiKIoyHlknAOsXe3A5HTQ0+Wd6KIqSRDaZVJXLj0x/f1knAPkuJ9csLqOhWdcDXCoMDUd441TXTA9jxsnLyyMQCKgIKDOCseoB5OXlpX1MVkUB2Wxc6uWh544wODRMXq5zpocz69m67wyf/8VrbP/STSwqSy8i4nKkoqKC1tZW2traZnooyizFrgiWLlkpAHOKo+Ff3QNDKgCXAGe7BjEGjrX1zWoByM3NTbsSk6JcCmSdCQigND9a77JrYGiGR6LAhb/DyfaJZY5UFGVmyGoB6B5UAbgUsP8OLR0qAIqSTWSlAJTkRS1XOgO4NOgaiK7J0BmAomQXWSkAagK6tLD/Dq0qAIqSVWS3APSrAFwKqA9AUbKTrBSAkpgPQNNBXAr0WALQ0T9Ej/plFCVryEoByHU6KHQ51QR0idA1MER5QVSUW9oHZng0iqKkS1YKAERnASoAM48xhq6BIVYtKgU0EkhRsomsFYBSFYBLgv7QMOGIYeVCSwDUD6AoWUPWCkBJfi7dKgAzjr0GoMpTQHFejgqAomQRWZkKAqIzgIu92ZzrHkyoK7CwLD9laonBoWEcIinrEPcFw+TlOnE6JKPv7guGKXRP3eXv7A/R3heKffYWuim17PRj0TM4RHFe6v16g2GKRozZnoWV5udS5SkYNxKoPxQmL8eJI8X16g2GKXQ5EcnsWiqKMjHSmgGIyK0iclhEmkTkvhTbF4vIcyLyuoi8ICIVVvtbReS1uNegiNxpbXtMRI7FbVubycBL8i5uBnDkXA+b/tdz3PSt38Ven3t8T8p9//ifdvDgf+xPag8PR3jrN1/gH59vyui7z3cPsu7BZ3nxyNQkDRscGuaGv3sh4dxu+tYLDA1HxjzuZKCfdQ8+y/Yjyam2dx5rZ83/fIajbb0J7XYobml+LpXlBbR0jO4EDoaHqf/G8zy+qyVpm783yPqvPstzBy+uxqmiKOkzrgCIiBP4HnAbsAK4W0RWjNjtm8CPjDFXAw8CXwcwxmwzxqw1xqwFbgL6gWfijvvv9nZjzGuZDPxifQAHz0aLOv/F7cv57l1r2bCknENnkws9D0cM+09388z+c0lpfve2dnG+J8izB89m9N0n2vsJDUd4vXVqUii/erKDroEhPvPWWr5711o+8ZZqAn2hcb/v4NluwhHDsweSz+fZA2cZjhiO+fsS2u2/QUl+DlXeAlra+4lEUqdDbu0YoKN/iH0pUkc3ne9lcCiScpuiKFNDOjOAjUCTMeaoMSYEPA7cMWKfFcDz1vttKbYD/CHwG2PMpBiJS/Nz6QsNj/tUOxq2+egDm6u4Y+0iNizxcKpjgOERN68zXQOEI4bzPUGaRzz9NlpFafaf7qazP0S6+HuCCWOYbBqbAjgdwiduqOaOtYv4xA01VvvYRXTs8aSqtdDQFG3z9wYT2u21GNEZQD7BcIS2EfvY2OahVOc91jZFUaaGdARgERA/Z2+12uLZC7zXev8eoFhEvCP2uQv4+Yi2r1lmo2+LiDvNMQNQmh+1RfdMcDFYS3s/viIXBa5oP5WeAsIRw5mugRH7Xfhs3wRjn5v9FLqcGAM7MihQ47ds81MVMtnQ7GdNRWnMlu8pdLFiQQkNzekJQNP5Xs51D8ba2/tCHDjTDYC/N1Ho4n0AlZ6ChH5GYqeKSHXeY21TFGVqmKwooC8CN4jIHuAG4BQwbG8UkQXAauDpuGPuB5YBGwAP8KVUHYvIvSKyW0R2xxfasB2aEzUDnWzvj92wIBrFYrfHY9/MClzOhDKUA6FhXj3RyR9tqKTQ5Rz35hpPwHpCnorUCT2DQ7ze2kV9rS+hvb7Wy6snOhkIDY9yZHQ8Ba6oE7wx7nzixS0wigAU510QgNHOy25PNdOyt2k6CUWZPtIRgFNAZdznCqsthjHmtDHmvcaYdcBfWG2dcbv8EfBvxpihuGPOmChB4FGipqYkjDEPG2PWG2PWz5kzJ9ZekndxAtDS0U9l+QUBsN+3jljJ2tLRj9Mh3LZqAS8dDcRuXLtPtBMajvCWK+ewqdpLY1MGMwBLAE53DhKeoAlrNF4+2s5wxFBXkygAdbU+QsMRdp9oH/XYlo4B6mt9lBfkJsx2Gpr9FLlzWFSWn2wCGhiiOC8Hp0NYVJaPyOirge32lDMty3l8rjvI4NDoIqUoyuSRjgDsAq4QkaUi4iJqynkqfgcR8YmI3df9wCMj+ribEeYfa1aARGP+7gTeyGTgF5MRNDwc4XTnYOypH2BBWR5OhyQ9gZ5s72dBaR5vudJH92A4Vvu2oSlAjkPYuMRDXY2Xo/4+TnemlwbBfooejhjOdA2Os3dmNDT7yct1cM3isoT2jUs85DgkyYxlY4yhpb2fxZ4CttR4aWzyx5zejU1+Ni31MK/ETaAvWQBsMc7LdTKvOG/MGYA9w0h1ne1trWNEEimKMnmMKwDGmDDwGaLmm4PAE8aY/SLyoIi829rtRuCwiLwJzAO+Zh8vIkuIziB+N6Lrn4rIPmAf4AO+msnAY0VhJiAAZ7oGGY4YKj0XyhfmOh0sKM1LskG3tPdT5SmIPVFvt8xADU1+1lWVUejOiZlbGsZxstoEekPkOiXW/2TS2BRgwxIP7pzE9QyF7hzWVZUlmHbiaesJEgxHqPJGz/V01yDH/H20dvRzPNBPXa0PX5E7pQnI/ltA1JSWyo5vC8ympR4gcaY1EBqmrScY26Z+AEWZHtLyARhjthpjrjTG1Bhjvma1fdkY85T1/kljzBXWPn9imXXsY48bYxYZYyIj+rzJGLPaGLPKGPNBY0xiiM04XMwMwH76jPcBACkXMp1sH6CyvIA5xW6umldMY7Ofzv4Qb5zuionCVfOK8Ra6aEzTEezvDbLCSp0wmTbvtp4gh8/1JJl/bOpqfOw71ZUyjXb8NYkJWnMgZtqqr/XiLXInmYBGCkCFJz+lqHUNDNETDLNxqTdpptVq3fDt79VIIEWZHrI6FQRMTADsG0y8DwCsp9cRT6b+3iBV3uh+dbVedh/v4HdvtmHMhRuWwyFsqfHSEGc2GQt/b5DVi0pwOmRSn3btp/v62pEBWFjtvmjE0tFkobLHUVlewBJvAQtL82hs8tPQ7MdX5OKqecX4ily094USHLipZgBnuwcJhhPt+PZ1XeorZGFZ4kzLFoN1VeW4cxwqAIoyTWStAOTlOnHlOCZkAjrZ3k+OQ1hQmpfQXukpwN8bpD8UDS21b1IV5VFTUX2Nj2A4wve2NZGf62Rt5QU7e32tL+VagZGEwhG6B8PMK85jYVkeJycxfXJDk5+SvJxYYraRrK0sIz/XmdIMdDIQHUdFeT4iQl2tjx1HAzQ2B9hS40NE8BW5iRgS1jx0Dw5Rkn8hPURleQHGRCN9EvqPzTDyqSxPnGnZN/wqTwGVaaSTUBRlcshaAYCoGWgiheFbOgZYWJZPjjPx9G2TkO2EjL8xAWyq9uB0CG+e62XjUk9CbqD6GtsPMLYZyHaieovc1oxjcm52xhgamgJsqfGOmpfIleNg41JPSl9FS0c/80rcsVxI9bVeOvuHaOsJUl/jtcbsAhLXAiTNAKzZ0siUELEZhqcgaaZ1sn2A/FwnviKXZYZTJ7CiTAdZLwAT9QFUjbD/A1RaT/onA4kx6bYwFOflsqYi+nR93Yg4+ypvARXl+eM6gm0nqq/IFc2dM0kCcLK9n1OdA0njGsl1tT6a2/o4OyL6aOQ1qY/zI9imLm+h2zqHqIgFw8MMDkUSBMA2q6WK8ikvyKXEWi8wcqZV5SlARKjyFNDa3p+WKU1RlIsja7OBwsQFoLW9n1tWzktqt2+A9tNqS/sABS4n3kJXbJ/6Wh+vnuykLoWdvb7Gx2/eOIMxZtSMlrYT1VvkptJTQKAvlJAZtLmtl/d+v3HMBVupiFg3zLpxBMAed2Ozn/deUxFrb23vZ3P1hXOaW5LHFXOLCIYjMQGcU2zNAKyVzPGrgGPHFbtx5Tg4MSJnUEvcwrv4mdaV84qtbVHxrSjPpycYprN/iPK4664oyuST1QJQkpczat6Z0egLhgn0hZIigCCaMqHA5UxYlVpZXpBwM/9o/VIWleWzYkFJ0vG1c4voHgzTEwzHYuNH4o+fAcQJzrL50f6e2X+OroEh7n1LdcYppheW5lHtKxxzn+XzS/AUumhoCsQEIBge5kz3YNI1+cb7VhMevvAkbs8A7FxG3QPRJ/iSOAFwOISrF5Wy60RHQl8t7f2stKqGxVZdB/q5Ym4RJ9v72WKZmeJFWAVAUaaWrBaA0vxcmtv6xt8xjvhol5HYJgjbPt3a0Z90U/QUurhrY1XKvmM28p7gqAJgm098lg8AojMNWwAam/1cOa+IP799eUbnlS4Oh7Cl2ktjsz82UzndOYgxyWGx1y72JHwuzc/F6ZCYH+NCJtDEc62r9fGPzx+J+QeGI4ZTnQPcumoBcMHU1tLRT3tfiP7QcOzvEZ9O4uqKxMVsiqJMLrPOB2Db91P5AAAqLLu8McbKF5Sfcr9UeIssG3nf6JlBA30h8nIdFLicF3wO1owjGB5m1/H2UeP4J4u6Wi9nrIVe8d8/2jWxcTgEb6Er5sfoTmECAqiv8RIx8LIVbnq2e5ChYRPr31PootCaaY387so4UVQUZWrJegHoHhwaNf98KkZbBGZjr2QNWE+m490U4/FZM4DAGGYpf08Qb6EbEYndCG1H8KsnOhkciozryL1YYhFL1sK1kdFOYxG/GCyVDwCi8fzRcNPE/m0xFREqrZmWHS1k/z2K3Dl4Cl0aCqoo00BWC0BJfi7GQE8w/ZTQrR0DFLlzKB+lPGKlJ5/+0DB7W6K57DITgOgMoK139BmAvy+Erzi634UbYfRm19jsx+kQNlV7Rj1+MljsLWBRWX6sPkBLez+uHAdzi8fPyO0rcsX8GHYI7kgBcOU42BAXbppqhmGf90hxsLe1ajoIRZlysl4AILN8QCfb+2OLnVJh36TseP7RZgqp8BSOPwMI9AbxxTk3K+Ny5zQ0+bk6Lo//VCEi1NV42XE0QCRiaOmIXpNUdXpH4oufAVgpJVL5O+prvBw538v57kFa2/txSLTmsk20fGQ/JwOJdRmi2/J1BqAo00BWC8BE8gG1tKdeA2Bj3/Dt1bKpnMWjket0UFaQm5QwLR5/bzDmLIYL6Sd6BofY29qVEH8/ldTX+ujsH+LAme5YtFM6xPsAugaGyLdWZKfqH6CxOWBlVM0nN27hXZU902rtTJmTKVXNAEVRJpfLQgDSnQEYE33aHeup3r4RHjrbg6/ITb7LOeq+qfAWupISpsV/f6A3FDMVRb8vn4GhYbbuOxPN4z9KHp/Jps4Ku2xo8tPSPpC2qctX7GZgaJj+UDhpFXA8KxaUUFaQG+2/YyDJmW7/DQ6d7UkSn9GqsymKMrlktQBkWhSmrTfI4FBkzJtdvsvJHMsWXpVBBJBNqpTJNt0DYcIRE4sWggupE362swV3joNrqsoz/s6JYC/02vrGWboGhtKOdrIXxfl7QnQPji4AdrhpQ5OfE4HkWVf859G2aSSQokwtWS0AmZaFTOVwTIUdnpmJ/d/GV+TG35d6BtAWWwMQ5wOwnn73tnSyYYknlotnOqiv9WXs7LZnL/6+IF0DiYngUvV/umsQf28w6Sm/Ir4a28jZQfnYtYUVRZkcslsAbBNQmgnh7CfK8W529vZMIoBsvEWu2ErZkcQvArOJvxGOrOM71dhmIEhf7OyxB3pDdA2ER50BQOL52DMdm/iZ1sjvtquzaWEYRZlasloACl1OnA5JewZgR5ZUjOPwjOWsycABbOMrctM9GCYUTq71a4dPxjuB42+Eo+Xxnyo2VXuxA3/SFYALGUGD0XKQYwiAXVcAUl/z2ExrxDa7OptGAinK1JLVAiAiGa0GPtnez9xi97hmlpFJyzLBvkG2p1gNHEsFXZgYb1/lKRgzj/9UUZqfy+qKMsqsLJ3pEB/q2j2GExiI1RWA1LOpKk9ByroM9rbnDp7n3f+4nXf/43b+6Ic7Uq4NaGnv5xM/3j2hpIAT4Wcvn+T7LzRNy3cpylST1bmAIJoQrmsgvYVge052sDxFEreR3LRsLh/avJh1VZnnooklTOsNMn/Ejc3fG0Lkwk3U5k+uW0rnwFDGyd8mg8+9rTaWHiMd8nKdFOflcL4nOGbSO5uP1i9hTrE7we9hc/fGKpYtKEmqywDw4S2L+cWuFgAiBn73Zhu/2XeWj7+lOmG/X+05xdP7z3H76vPcsXZR2ucxEYwxfG9bE539IT5+fXVCWKuiZCNZLwCl+blphYGe7Rqkua2PP95QOe6+viI3X7lz1YTGE0uZnCIU1N8bxFPgSrrR37Z6wYS+azK4aVlyWuzx8BW5Y3mExpoBAKxcWDrqzGZTtZdN1anNXreuWhBLHgdw07deoKHZnyQADdZ6jYYm/5QLgF1zAeD11s6kZHmKkm1k/SNMSZomIHth11QnWrtQNCWFCWjEIrBsxVvoovl8tPTleAIwWdTX+Nh5rD3BtzIQGubVE9EopoamwJQXkYmv9jZe5TdFyQbSEgARuVVEDotIk4jcl2L7YhF5TkReF5EXRKQibtuwiLxmvZ6Ka18qIi9bff5CRCZ0Z0x3BtDQFKC8IDdlHv/JxM7zk2oGMHIRWLbiK3Jz2qooNm0CUOuNrRy22X2indBwhJuXz+VU58CUO40bmv3MK3GzcmHJuJXfFCUbGFcARMQJfA+4DVgB3C0iK0bs9k3gR8aYq4EHga/HbRswxqy1Xu+Oa/8b4NvGmFqgA/jYRE4gnRmAMYbGZj9barxp5bu5GApdTtw5jpQpoaNpILJfAOJnMaWjJNWbbDZXexEh4cbb0BQg1yl87m1XArB9Cm/KkYhhR3OA+hof9bU+9pzszLhqm6JcaqQzA9gINBljjhpjQsDjwB0j9lkBPG+935ZiewISzcR2E/Ck1fR/gDvTHXQ8dkrosab/x/x9nOkanHLzD0QjX+ITpsUT6A0llJfMVuJFLN3ooYulrMDFqoWlNMaZXhqb/ayrLGfVohLmlbgTtk02h8720N4Xoq7WR12Nl9BwhF3H26fs+xRlOkhHABYBLXGfW622ePYC77XevwcoFhHbu5cnIrtF5CURsW/yXqDTGGOH76TqEwARudc6fndbW1vS9tL8XIaGDQNDoz+N2Xnvp2uhVXzKZJvBoWF6guFYzH82Myd+BjBNJiCIFrJ59WQHfcEwnf0h9p3qoq7Wi4hQX+OjsdmfUW2ITLBnHvW1XjYu9ZDrFDUDKVnPZDmBvwjcICJ7gBuAU4B9R15sjFkPvB/4jojUZNKxMeZhY8x6Y8z6OXPmJG1PJyNowxE/C0vzWOLNPK5/IniL3EkpoW2T0OU2A5hOAbiu1kc4Yth5vJ2XjgYw5oKo19X66Ogf4uDZ7in57oZmP9W+QhaU5lPgymFdZXksAklRspV0BOAUEB87WWG1xcPSXo8AACAASURBVDDGnDbGvNcYsw74C6ut0/r3lPXvUeAFYB0QAMpEJGe0PtNlPAEYjhh2HA1QV+sbtQbAZBOfMtnGFoTLwgdgiZjL6SAvd/oCydYv9uByOmhs8tPQFKDA5WSNVTfYXkU9FWagUDjCzmPtCZla62q97D/dTWf/6Km/FeVSJ53/vbuAK6yoHRdwF/BU/A4i4hMRu6/7gUes9nIRcdv7APXAARM12G8D/tA65iPAv0/kBGIZQftTC8CB0910DQxNeZnFeHzFbgJ9wQS/hC0IqRZEZRt2pFNJfs60iSpE02Zcs7iMhqYADc1+Ni71xGoRLCjNp3pO4ZQ8le9t7aQ/NJzwG7qu1ocxsKNZw0GV7GVcAbDs9J8BngYOAk8YY/aLyIMiYkf13AgcFpE3gXnA16z25cBuEdlL9Ib/DWPMAWvbl4AviEgTUZ/A/57ICVxICJd6NXBDLP5/+vLseAtdDA0buuNWKLelSASXrfgKbQGYPvOPTX2NjwNnujna1pdUPCfVWoHJoKHJj0g0EslmTWUZhS6nmoGUrCatlcDGmK3A1hFtX457/yQXInri92kEVo/S51GiEUYXxXgmoIYmP1fMLWJuSXK+makiPmWyHSYZSJEILlspyc8h1ynTav+3qav18a1n37TeJ4p6fa2XH790gr2t0dTak0VjU4BVC0spK7jwt8t1Oti41DOlkUeKMtVcFqkgAH768glePpr8n3HnsXbu3lg1rWOKCUBPkJo5RUDUB1DgcibUvs1WRARvoXtGBGBNRSlF7hxcOQ6Wz09c1Be/ViCVAPz4pRNcW1XOioWJxxlj+PZvj3CmM3UBmj0tHfx/1y1Naq+v9bHt8EG+8MRrOC1T2O2rF/DWZXMnenqKMq1k/d2oOC+HLdVeTgT6OGetTo1nfmke71qzcFrHZD/lxy8GO3yuZ0L1BS5V3r12YSyd83SS43TwkbrFuJzOpEV98WsFPn9z4nFtPUH+6ldv8M6rF/C991+TsO3Q2R4eeu4I3kIX7hT1jReV5fOuq5N/Q+9YOZ+f7zzJS5YfoKN/iH2nulQAlKwh6wXA4RB+fu/mmR5GAheKpkTt/sHwMLuOT/9MZCr589uXz9h3//d3LBt1W12tl0e2H6M/FE6Ybdm5oHY0B4hETIJ42PH8//ez17GwLH1Rq/QU8Nyf3Rj7/L1tTfzd04fx9wYvC1+PcvmT9cngLkXKC3IRgTbL7v/qiU4GhyJJTktl8qmv8TE0bNh5LHGVrm2rb+8LcehsT+K25gBLfYUZ3fxTfrcVJdSokUFKlqACMAXkOB2UF7hiM4DGZj9Oh7CpWtMHTzUbllhrBUbchBua/aypiKalboyL3BkajvDy0cCkRImtXlRKcV4OjbpCWMkSVACmiGg6iKgANDT5ubqilOJpypszm8l3OVlXVZaQpuFkoJ/WjgHee00F1b7ChG2vt3bSFxqelDQhToewudqroaFK1qACMEV4C90EekP0DA6xt7VLzT/TSH2tj/2nu2NlOe0bcn2tl7paLy/HrRXYfiSACGwZpTBNplxX66OlfYAWrWesZAEqAFOEt8hFoC/EzmPtDEdMUsy6MnXYaSHsVboNTX7mFrupmVNEfY0voa5AQ7OflQtLKJ+kHE32d2uiOCUbUAGYInxFbvw9QbY3+XHnOLimqnymhzRruLriwirdWB5/KxfUlpoLawX6Q2H2nOyY1NlZzZwi5ha7p7Q2gaJMFioAU4SvyEVPMMwLh9vYsMRDXq5zpoc0a8h1OthU7aWxyc/hcz0E+kIxG3/8WoFdxzsYGjbUTWKeKBGhvtYXCzdVlEsZFYApwo4DP+bvU/PPDFBX4+V4oJ9f7m4FLphmILpWYE9LB789cI5cp7BhyeTOzupqvAT6Qhw+1zP+zooyg6gATBHxaZ/VATz9XHdF9Jr/5OUTsTz+NvZagV/sbmFdVfmkp+ewZxvqB1AudVQApgg7HURJXg6rFpXO8GhmH1fNK8ZX5CIUjiTNwOy1AqHw1CzOW1iWz1JfoS4IUy55sj4VxKXKHGsGsLnai3OKC9EryUQdvj7+797TSTd5e63Ay8faE0xDk0ldjZdf7TnF/95+LNb29uXzqEpRle63B85xYpSw0RyHcOe6RTOSeC8baO3o51x3kGsXa5DFRFABmCLmFLuZW+zmnVcvmOmhzFreuXo+24+0UZfiKf+dVy+gtWOANZVlU/Ldt6ycz892nuQr/3Eg1rbrWDs//NC1Cfu194X4+I93Y8bwF3cPDPHZt10xJePMdr7yHwfYfbyDV/7q7TM9lKxEzFi/vEuM9evXm927d8/0MNLGGDOtFbOUS4u+YJiwFQn0P5/az3OHzvPqX709YUb469fP8OmfvcpPPraJ1RXJpsK7Hn6J0vwcHr93y7SNO1sYjhjWPvgMPYNh9v31LbrSfgxE5BWrNnsC6gOYQvTmP7spdOdQmp9LaX4ub7lyDl0DQxw4nVi0vqHZT5E7h83Vnti+8a/rar28eqKTgdDwDJ3Fpcu+U130WJUAW9pT13JQxkYFQFGmATvZ3MgFYg1NfjYt9ZDjTP1fsa7WR2g4wu4T7Sm3z2YS8j1p6o0JoQKgKNPA3JI8rpxXlJCJtLWjnxOB/jET0W1c4iHXKbqyOAWNzX4WWSm8WztUACaCCoCiTBN1NT52HW8nGI6ac+waBWMJQKE7h3WV5Vp7eASDQ8PsPt7BLSvnUZyXozOACZKWAIjIrSJyWESaROS+FNsXi8hzIvK6iLwgIhVW+1oR2SEi+61tfxx3zGMickxEXrNeayfvtBTl0qO+1sfgUIRXT1xIROcrcnPlvKIxj6ur9fLG6S46+0Nj7jebePVEB8FwhOtqfVR5ClQAJsi4AiAiTuB7wG3ACuBuEVkxYrdvAj8yxlwNPAh83WrvBz5sjFkJ3Ap8R0Ti4+7+uzFmrfV67SLPRVEuaTZVe3BI1HRhjKGxOVqIZrxggfpaH8bAS0d1FmDTYBVZ2rjUQ5WnQNNvT5B0ZgAbgSZjzFFjTAh4HLhjxD4rgOet99vs7caYN40xR6z3p4HzwJzJGLiiZBsleblcXREtVnPkfC9tPcG0FqKtqSijwOWkQc1AMRqaAqyxiixVegpo6RjQ5HsTIB0BWAS0xH1utdri2Qu813r/HqBYRBJ+2SKyEXABzXHNX7NMQ98WkZRVtEXkXhHZLSK729ra0hiuoly61Nd62dvaxTP7zwKkXKQ2EleOg41LPVppzKJ7cIjXWztjvpNKTwGhcIQ2qwKfkj6T5QT+InCDiOwBbgBOAbHAZRFZAPwY+KgxJmI13w8sAzYAHuBLqTo2xjxsjFlvjFk/Z45OHpTspr7Gx3DE8M8vHqPKU0ClJzk1xGjHHW3r42zX4BSP8NLn5aPtRMwF8awsj0YCqR8gc9IRgFNAZdznCqsthjHmtDHmvcaYdcBfWG2dACJSAvwa+AtjzEtxx5wxUYLAo0RNTYpyWXPN4nLcOQ66BoYyykNUp5XGYjQ0+cnLdXDN4qg7scoSUfUDZE46uYB2AVeIyFKiN/67gPfH7yAiPqDderq/H3jEancB/0bUQfzkiGMWGGPOSNQDdifwxsWejKJc6uTlOlm/pJyGpkBa5h+b5fNL8BS62LrvDAtK8zL+3pULSyktuDxSJTQ2+9mwxIM7J1pkaVF5PiI6A5gI4wqAMSYsIp8BngacwCPGmP0i8iCw2xjzFHAj8HURMcDvgU9bh/8R8BbAKyL3WG33WBE/PxWROYAArwGfnLzTUpRLl7deNZddxztiq4PTweEQrr/Cx7+/dprnDp3P+DvfsXIe//ShpFQwWUfP4BBvnuvl3WsWxtrcOU7ml+RpOogJkFY2UGPMVmDriLYvx71/EngyxXE/AX4ySp83ZTRSRblMuKduCbeump9QNCgdvnrnKt6/sSrj73us8Tjbj/gJD0dGTTmRLZzviTp6K8oTfSeV5RoKOhE0HbSiTDM5TkfSDSwdivNy2VSdef2Ctt4gv3njLK+f6uKaquzOmx/ojS6Gswsu2VR6ChLSbCjpkd2PA4qijMsWSzQaLwMHcsAK9fSNmD1VevI52z0YS7OhpIcKgKJc5niL3KxYUHJZLCTzWwIwcgZQ5SnAGDjVoX6ATFABUJRZQH2tl1dOdjA4lN1PyP7eECLgKUg2AYFGAmWKCoCizALqan2EwhF2H++Y6aFcFP7eIOUFriRndmwtgM4AMkIFQFFmARuXeMhxSNankwj0hvAWupLa5xS5ceU4NBIoQ1QAFGUWUOjOYV1VWdY7ggN9wST7P0TXSVSW53MyoAKQCSoAijJLqKvxse9UF10DQzM9lAnj7w0lRQDZRLOCqgBkggqAoswS6mt9RLK8roC/NziqAGhhmMxRAVCUWcLayjLyc51ZawYKhofpGQyn9AFAdDVwz2CYrv7sneFMN7oSWFFmCRfqCmQ+AzDG0NE/hDHRois5Tgel+Zknl4tEDA7H2BXQRsNeBewrHt0EBPDG6S6WzS8GoDQ/N+vTX0wlKgCKMouor/Xyv7Yeoq0nyJxRbqSp+PZvj/DQc0cS2h69ZwNvXTY3oe24v49bv/t7fvXpepbNL0nY9lpLJ+//55f4z8+9hSpv5qkwYmkgRpkBLPFF+/zAv7wca7v+Ch8//timjL9rtqACoCiziCvnRZ+Mjwf6MhKAZ/afZfmCEu7eWIkx8MBT+3m9tStJAPad6mJwKMKhMz1JAnDwTDf9oWGeP3SOe+qXZjz2C6uAU4/7qnnFfPeutTEn94tH/Dx/6Dw9g0MU510eqbAnG50bKcosYiLFU/y9QQ6d7eFdaxbw4S1L+Ejdkmj65RQRN3abP0V5RjuPz0RMUPF9zhlFAESEO9Yu4sNblvDhLUv4aP0ShiOGncfaJ/R9swEVAEWZRUykeEqjdcOujytgU+nJT9mHLSyBvlDSNr9lwnnpaIDwcCRp+3jYfaZaB5CKa6qi1de2Z6nTezpQAVCUWcREiqc0Nvkpycth1aLSWFulp4DWlAIQ7dffkzwDsJ/gewbDvHG6O9Oh4+8Jkp/rpNCdnuU6L9fJhiUeGi+DJHhThQqAoswyMi2esr3Jz+ZqL8646J3K8gLOpEi/fHKMGUCgN0S1rxCYWG3jQF8o7ad/m7paL4fP9dCWQpAUFQBFmXVksmL2ZKCf1o4B6msT6xfb6ZdPdw7G2sLDEU53RmcAgRQ+AH9vkCvnFbNsfvGEBMDfG8y4ipptttJiMalRAVCUWYZdPCWd1NB28rj62sRKZKnSL5/pGiQcMTgdErP3x2M/wdfX+th9IvPU1P7eEHMynAGsWlRKSV6OmoFGQQVAUWYZseIpneP7ARqa/MwtdlMzpyipD0iMJrJnFcsXFOPvDcYWjUF0dtDRH8Jb5Ka+1ksoHOGVE5mlpg70BvEWZjYDcDqEzdXerM+COlWkJQAicquIHBaRJhG5L8X2xSLynIi8LiIviEhF3LaPiMgR6/WRuPZrRWSf1edDIjKx5YGKomREuqGgkYhhR3OA+lofI/97zi1OTr9sv19XWU4wHKEvdOEJv70/hDEwp8jFxqXeaGrqDMxAkYiZkA8AojmQWjsGNFNoCsYVABFxAt8DbgNWAHeLyIoRu30T+JEx5mrgQeDr1rEe4AFgE7AReEBE7KrUPwA+DlxhvW696LNRFGVcKtMUgMPnegj0hairSS5E73AIFeX5Cb6ElvYBnA5htRUtFB8JdKGYu5sidw5rKssyWg/QNTDEcMSMmghuLGz/hc4CkklnBrARaDLGHDXGhIDHgTtG7LMCeN56vy1u+zuAZ40x7caYDuBZ4FYRWQCUGGNeMtF54o+AOy/yXBRFSYM5RW7cOY5xq2fZT+gjHcA2leWJ2TdPtvezsCyPeaV5QDR3v00sj491A6+v8bKvtTPt1NR2XxOZAdTMKWReiXtCjufLnXQEYBHQEve51WqLZy/wXuv9e4BiEfGOcewi6/1YfSqKMgXYT+/jmUQamwMs9RWysCw/5fYqT0FCHyfb+6ksL4jl6ol3BI8s5l5npaZ+Oc3U1G09iQKSCSJCfY2PHc0BIhEz/gGziMnKBfRF4B9F5B7g98ApYFKqT4vIvcC9AFVVVZPRpaLMeqpShII+s/8sX/73/QxbzttAb5C7N47+f67Sk0+3lX65tCCX1o5+bl4+L3aTjk8HYb/3WU7cdVVl5OU6+G+/eI2CFAu7BPjz25dz57roc6E9A5iIAEBUcP7/Pac4dLaHFQsTcxS9caqLv/zVGzz20Q2UFWQ+w8hm0hGAU0Bl3OcKqy2GMeY01gxARIqA9xljOkXkFHDjiGNfsI6vGNGe0Gdc3w8DDwOsX79e5VtRJoFKTwG7R0Th/OurrQTDw9y6agEATgd8dIykbRcKsfeT4yzE3xui0lOAx5oBBBJmACFynUJJfvSW485x8pU7VvHqyc6UfT938By/fKXlggD0ZpYGYiR2GGtjsz9JAH615xSvtXTyuzfbuGPt7DJEpCMAu4ArRGQp0Zv0XcD743cQER/QboyJAPcDj1ibngb+V5zj9xbgfmNMu4h0i8hm4GXgw8A/XPTZKIqSFlWeC8VTSgtyGbYifm5dNZ+vv3d1Wn3EO5NzrZz7lZ4CXDnRWgHxi8HsEM74aKL/sr6S/7K+klR85T+c/OSlEwwODZOX68TfG8QhUD7BJ/QFpflU+wppaPLzJ9dXJ2yzndGNTYFZJwDj+gCMMWHgM0Rv5geBJ4wx+0XkQRF5t7XbjcBhEXkTmAd8zTq2HfgKURHZBTxotQF8CvgXoAloBn4zWSelKMrYVJQnLuTaf7qL7sHwqA7fVMQvBrP7sWcF3iJXgg8g0xDO+lovwXCEV61Zir83hKfQlZCOIlPqar3sPNbOUFwiukBvkINnuhGZnVFCaa0DMMZsNcZcaYypMcbYN/cvG2Oest4/aYy5wtrnT4wxwbhjHzHG1FqvR+PadxtjVll9fsbErxpRFGVKiTffADRYK2W3pAj5HI2SvFzKCnJp6eiPhZRWlkcdxr4id5IPIBP7fWytgHVT9k9gEdhI6mt89IWG2dtywey0w3JC37Fm4axcK6ArgRVlFlLpid6o7Sf3xmY/V84rYm5xXmb9lBdwsn2Ak+39FLqcMfu/r8iVkBAu0JvZDCC2VsASpkBvcML2f5stNd7ok35cWoiGpgDF7hw+cUNN9PMsmwWoACjKLKQ4L5fyglxa2vsJhofZdbydupr0zT82VVZa6NaOfio9BTEbv7fwwgzAGJPxDACiawVeb+2ke3CIQF9owhFANmUFLlYuLEm4yTc2+9lU7WHZ/GLmFs++tQIqAIoyS6n0RBdyvXqik8GhSEb2f5sKTz6tHQMcD/THfAIQNQF19g8xNByhNxgmGI7gyziVs71WoB1/z8XPACBqBtpzsoP+UJjWjn5OBPqpq4mmuqivnX1rBVQAFGWWUukpoLVjgMZmPw6BTdWejPuo8hQQGo7Q3NYb8yvAhXDNjr5QXDH3zJ7g7bUCzx86R19o+KJnABAVlaFhw67jHbEMobbw1dV4CfSFOHyu56K/J1tQAVCUWUpleQGtHf28eMTP1RVllEygcHqlFU1kzAUHMBB72m/rDSatAk4Xd060otfWfWcT+rwYNiwpJ9cpNDb5aWj24ytyc+W8aKbTWM6gWWQGUgFQlFlKlaeAoWHDay2dSfn+M+kj9t6baAKCqPPX3zvxNA71tb5YvqCLjQICKHDlsK6qnO1NfhqbA9TVeGN+i4Vl+Sz1FcZqIM8GVAAUZZZiRwJBYsH3TFhYFi0yDxdmA0CsclegL3hRaRzixzUZPgC7z/2nu2nrCSYJX12Nl5ePBhLWClzOqAAoyizFfnp35zi4ZnH5OHunxpXjYGFpVEgqypN9AP6eEH4rkZsdIpoJKxaWUJofNU1Nhg8AEqubjYx8qq+NrhV4vTV1iorLjclKBqcoSpaxsCwfh8D6JeXk5Ton3E9FeT6h4Qj5rgt9FLtzcOU48PcFGQgNU5qfiysn8+dNp0PYUu3lP/efnbQZwJrKMgpdTrxF7oTIJYAt1dG1AtuPBLh28dhO8ZePBmhoDvCFt1+Z0fc/vvMkv953Jva5yJ3DN953dUzophMVAEWZpeQ6HXz8LdVsXjox+7/N+zdVcbZrMKFNRPAVugj0hhgIDV/Uzfue+iXMK3FT4Jqc21Wu08F/fdsVKfMKlRe6WLEgulbgczdfMWY/D//+KM8dOs+nbqxJW0CNMXznt0cYNoaK8nzCw4YXj/h5+4p5vPeaivE7mGRUABRlFnP/bcsvuo/REqh5rXQQ/aHhWBroibC52svm6osTqZHYK39TUV/r49GGY/SHwqOKTng4wsvHomnNTnUOJNVMHo3mtj7Odg/ytfes4gObFhOJGK796rNsb/LPiACoD0BRlCnBVxSdAQR6g/iKsyfPfl2NN7ZWYDT2tnbRGwwDJFRFG49GaxXydVbIqcMh1NX4aGwKMBPp0FQAFEWZErxFbgK9wWgm0EkI4ZwuNi71xNYKjEb8ttYMBKChyc+isvyE8Nm6Wi9nuwc56u+b2IAvAhUARVGmBG+Ri7beIJ39Q5PmwJ0OClw5rKssHzMxXEOzn+ULSnDnONKeAdg1F+prvQl1EexQ17EEZ6pQAVAUZUqYU+RmaDhq1pisEM7por42ulagsz+UtG0gNMyrJzq5rtYby6eUDqPVXFjsLWBRWX5CltLpQgVAUZQpIf6pfzLSOEwn9bVejIEdKVYF7z7RTmg4Ql2tL1pbuX0grT5Hq7kgItTVeNlxNMDwNCeiUwFQFGVKiLf7e7NsBmCvFUhlBmpoCpDjEDYu8VBZnk9Le39aDtyxai7YKS8OnO6elPGniwqAoihTQrzZJ9tMQLlOBxuXemIZQ+NpbPazrqqMQncOlZ4CeoLhWL6i0Riv5kKdNSuY7oI0KgCKokwJ8WafbHIC29TX+jjq7+NM1wUTT1f/EPtOdcVu5PF1kcdivJoLc0vyuGJu0bRnIlUBUBRlSii3cv+4nA6K3dm35tS+ycc7Z3ccDWDMhdTRsdrK4/gB0qm5UF/rY9fxdoLh4YsdetqoACiKMiXkOh2UF+TiK3IlhD1mC8vmF+MpdCWEZzY2+8nPdbK2sgxIfwbQ0DR+zYW6Gi+DQxH2nJy+RHRpCYCI3Coih0WkSUTuS7G9SkS2icgeEXldRG632j8gIq/FvSIistba9oLVp71t7uSemqIoM423yJ11DmAbh0PYUuNl2+HzfO3XB/jarw/w9P6zbFzqiSW2K3Ln4Cl00dIxugD0DA6xt7Vr3JoLm6q9OGR61wOMOy8TESfwPeDtQCuwS0SeMsYciNvtL4EnjDE/EJEVwFZgiTHmp8BPrX5WA78yxrwWd9wHjDG7J+lcFEW5xLjxyjm4c7PX0HDHmoX8/s02fvrySQAcIty5bmHCPnYk0GjsPNbOcMSMW3OhND+X1RVl0QyjFz/0tEjHMLcRaDLGHAUQkceBO4B4ATBAifW+FDidop+7gccnPlRFUbKNv/yDFTM9hIvilpXz2bdy/pj7VHoKeONU16jbG5oCaddcqK/x8vDvj9IbDFM0DX6TdKR5EdAS97nVaovnr4EPikgr0af/z6bo54+Bn49oe9Qy//yVjGIkFJF7RWS3iOxua2tLY7iKoijTR6WngFOdA6Mu4mps9qddc6G+1kc4Yth5bHpWBU/W3Oxu4DFjTAVwO/BjEYn1LSKbgH5jzBtxx3zAGLMauN56fShVx8aYh40x640x6+fMmTNJw1UURZkc7NrKZ7sHk7b5e4McOtszavz/SK5dXI4rxzFtaSHSEYBTQGXc5wqrLZ6PAU8AGGN2AHlA/BnfxYinf2PMKevfHuBnRE1NiqIoWYVdC/lkINkPYBeYHy3+fyR5uU7WLy6ftvUA6QjALuAKEVkqIi6iN/OnRuxzEngbgIgsJyoAbdZnB/BHxNn/RSRHRHzW+1zgD4A3UBRFyTIurAVIIQBNforzcli9qDTt/uprfRw624O/NzhpYxyNcQXAGBMGPgM8DRwkGu2zX0QeFJF3W7v9GfBxEdlL9En/HnMhOcZbgBbbiWzhBp4WkdeB14jOKP55Us5IURRlGllQlodDSBkK2tDsZ3O1F6cj/XUQdlqIVInoJpu03MzGmK1EnbvxbV+Oe38AqB/l2BeAzSPa+oBrMxyroijKJUeu08HCsvykxWAt7f20tA/wsfqlGfW3elEpxe4cGpv9vGvNwvEPuAiyN0BXURTlEqGyvCDJBGTb8dO1/9vkOB1sqvZOiyNYBUBRFOUiqfIUcHJEPqCG5gBzi93Uzk2vYHw89bVeTrb3j7nAbDJQAVAURblIKj35+HuDDISiidyMMexo9lNX451QHiR71tA4xemhVQAURVEuEjspnO0IPnyuB39viLoMzT82V8wtYk6xe8rNQNmXo1VRFOUSwxaAJ19p5ap5xbx8LLP4/5HYZSK3H/Hzr6+0ptxnUXk+m6vHTjA3HioAiqIoF0mNrwh3joOHf38h2n3Z/GIWleVPuM+3LZ/Hv792mj/75d6U250O4fUHbqHwInIGqQAoiqJcJKUFubz852+jeyAca/MVX1wVtHddvYD1i8sJDyfnGNp2+DwPPLWf1o4BrppfPOHvUAFQFEWZBMoKXJQVTF7pSxFh4SgziDVWQZqT7f0XJQDqBFYURckyqtKsRDYeKgCKoihZRnlBLoUu50WvE1ABUBRFyTJEhEpP8urjTFEBUBRFyUIqPQVj1iJOBxUARVGULKTKU0BL+wAXEi9njgqAoihKFlJZns/A0DD+3tCE+1ABUBRFyUKqvBcfCaQCoCiKkoXYpShbL8IPoAKgKIqShVSMUYs4XVQAFEVRspB8l5M5xe6LigRSAVAURclSooVoVAAURVFmHZXl+bSMqESWCWkJgIjcKiKHRaRJRO5Lsb1KRLaJyB4RhPAI5QAACmNJREFUeV1Ebrfal4jIgIi8Zr1+GHfMtSKyz+rzIZlI2RxFUZRZTJWngDNdAwwNRyZ0/LgCICJO4HvAbcAK4G4RWTFit78EnjDGrAPuAr4ft63ZGLPWen0yrv0HwMeBK6zXrRM6A0VRlFlKhaeAiIHTnRObBaQzA9gINBljjhpjQsDjwB0j9jFAifW+FDg9VocisgAoMca8ZKLL2H4E3JnRyBVFUWY5F5sVNB0BWAS0xH1utdri+WvggyLSCmwFPhu3ballGvqdiFwf12d8nbNUfQIgIveKyG4R2d3W1pbGcBVFUWYHsVrEE/QDTJYT+G7gMWNMBXA78GMRcQBngCrLNPQF4GciUjJGP0kYYx42xqw3xqyfM2fOJA1XURQl+5lfkkeuUyY8A0inItgpoDLuc4XVFs/HsGz4xpgdIpIH+Iwx54Gg1f6KiDQDV1rHV4zTp6IoijIGToewqCx/wmsB0pkB7AKuEJGlIuIi6uR9asQ+J4G3AYjIciAPaBOROZYTGRGpJursPWqMOQN0i8hmK/rnw8C/T+gMFEVRZjEXUxdgXAEwxoSBzwBPAweJRvvsF5EHReTd1m5/BnxcRPYCPwfusZy7bwFeF5HXgCeBTxpj2q1jPgX8C9AENAO/mdAZKIqizGIuRgDSKgpvjNlK1Lkb3/bluPcHgPoUx/0r8K+j9LkbWJXJYBVFUZREqjwFdPQP0T04RElebkbH6kpgRVGULMbOCjqRWYAKgKIoShZT6ckHJhYKqgKgKIqSxdiLwY4H+jI+VgVAURQliykrcLHYW8CuY+3j7zwCFQBFUZQsp67Gx8vH2glnmBROBUBRFCXLqa/10hsMs7e1K6PjVAAURVGynLoaHwCNTf6MjlMBUBRFyXI8hS5WLCihoVkFQFEUZdZRX+vl1ROdDISG0z5GBUBRFOUyoK7WR2g4wu4T6UcDqQAoiqJcBmxc4iHHITQ0BdI+RgVAURTlMqDQncO6qjIaM/ADqAAoiqJcJtTV+Nh3qouu/qG09lcBUBRFuUyor/VhDOw4mp4ZSAVAURTlMmFtZRn5uc60zUBp1QNQFEVRLn1cOQ42LvXwy92t7GgefxagAqAoinIZ8ac31lDkzsFgYm2/HWVfFQBFUZTLiM3VXjZXexPafvDB1PuqD0BRFGWWogKgKIoyS0lLAETkVhE5LCJNInJfiu1VIrJNRPaIyOsicrvV/nYReUVE9ln/3hR3zAtWn69Zr7mTd1qKoijKeIzrAxARJ/A94O1AK7BLRJ4yxhyI2+0vgSeMMT8QkRXAVmAJ4AfeZYw5LSKrgKeBRXHHfcAYs3tyTkVRFEXJhHRmABuBJmPMUWNMCHgcuGPEPgYosd6XAqcBjDF7jDGnrfb9QL6IuC9+2IqiKMrFko4ALAJa4j63kvgUD/DXwAdFpJXo0/9nU/TzPuBVY0wwru1Ry/zzVyIiqb5cRO4Vkd0isrutrS2N4SqKoijpMFlO4LuBx4wxFcDtwI9FJNa3iKwE/gb4RNwxHzDGrAaut14fStWxMeZhY8x6Y8z6OXPmTNJwFUVRlHQE4BRQGfe5wmqL52PAEwDGmB1AHuADEJEK4N+ADxtjmu0DjDGnrH97gJ8RNTUpiqIo00Q6C8F2AVeIyFKiN/67gPeP2Ock8DbgMRFZTlQA2kSkDPg1cJ8xpsHeWURygDJjjF9EcoE/YPTFajFeeeWVXhE5nMaYZxM+os525QJ6TZLRa5LMbLomi1M1ijEmVXviTtGwzu8ATuARY8zXRORBYLcx5ikr8uefgSKiDuH/YYx5RkT+ErgfOBLX3S1AH/B7INfq87fAF4wxY9YyE5Hdxpj14w54FqHXJBm9JsnoNUlGr0maAnCpoH+wZPSaJKPXJBm9JsnoNdGVwIqiKLOWbBOAh2d6AJcgek2S0WuSjF6TZGb9NckqE5CiKIoyeWTbDEBRFEWZJFQAFEVRZilZIQDjZSOdDYhIpZVx9YCI7BeRz1ntHhF5VkSOWP+Wz/RYpxsRcVqZaP/D+rxURF62fi+/EBHXTI9xOhGRMhF5UkQOichBEdky238nIvLfrP83b4jIz0Ukb7b/TiALBCAuG+ltwArgbmvdwWwjDPyZMWYFsBn4tHUd7gOeM8ZcATxnfZ5tfA44GPf5b4BvG2NqgQ6iK9VnE98F/tMYswxYQ/TazNrfiYgsAv4rsN4Ys4ro2qO70N/JpS8ApJeN9LLHGHPGGPOq9b6H6H/qRUSvxf+xdvs/wJ0zM8KZwUo18k7gX6zPAtwEPGntMquuiYiUAm8B/jeAMSZkjOlklv9OiGY9yLeyEBQAZ5jFvxObbBCAdLKRzipEZAmwDngZmGeMOWNtOgvMm6FhzRTfAf4HELE+e4FOY0zY+jzbfi9LgTaimXb3iMi/iEghs/h3YuUd+ybRlDVngC7gFWb37wTIDgFQ4hCRIv5fe/cTYmUVh3H8+8AUokIWujDDpiDahaMhQiH9kXAhZhQVKP6BtoILN7YTxWXbVomLZhMpNDsJdSFCCSo1aH+IxshFKrgQxcWQT4tzprmJ6AzM3Hdu5/ms7vvnwpnLufzu+zvvPC8cB/bZvt17zOWe3mbu65W0Bbhh+0LXY1lAhoC1wOe2RyixK/9p9zQ4T56mXAG9ADwLLAE2dzqoBWIQCsBM0kibUIPzjgOjtk/U3dclrazHVwI3uhpfB14Dtkq6SmkNvkXpfy+rl/rQ3ny5Blyz/X3d/ppSEFqeJ5uACds3bU8CJyhzp+V5AgxGAfg3jbSu0n8MjHU8pr6rve0vgJ9sf9ZzaAzYVV/vAr7p99i6YvuA7edsD1PmxWnb24EzwAf1tNY+k7+APyW9XHe9DVyh4XlCaf1skLS4fo+mPpNm58mUgfhP4IelkXY8pL6T9DpwFhhnut/9KWUd4CtgNfAH8KHtW50MskOS3gD2294i6UXKFcEzwCVgxwNPovtfk7SGsij+JPA7sIfyY6/ZeSLpIPAR5W66S8AnlJ5/s/MEBqQARETE3BuEFlBERMyDFICIiEalAERENCoFICKiUSkAERGNGnr8KRFtkPQ35TbbKdtsX+1oOBHzLreBRlSS7the+ojjQz3ZMREDLy2giEeQtFvSmKTTwClJSyWdknRR0rikd+t5wzV//5ikXyWNStok6VzN4F9fz1si6aik8zWsrblk21g4cgUQUT3QApqw/Z6k3cBh4BXbt6bihG3flrQc+A54CXge+I2S0nqZEmHyAyVjfiuwx/Y2SUeAK7a/lLQMOA+M2L7bv780osgaQMS0e7bXPGT/tz2xCQKOSNpIieRYxXS08oTtcQBJlykPYLGkcWC4nvMOJcBuf91eRIln6H2gTURfpABEPF7vr/PtwApgne3JmkS6qB7rzZG537N9n+nvmoD3bf8yf8ONmJmsAUTMzlOUZxBMSnqT0vqZjZPA3ppKiaSRuR5gxEylAETMzijwam3r7AR+nuX7DwFPAD/WNtGhOR5fxIxlETgiolG5AoiIaFQKQEREo1IAIiIalQIQEdGoFICIiEalAERENCoFICKiUf8AZp8IoEbhVCIAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "ca_df.plot(x='Frame')" ] }, { "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 }