{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Getting Started with Projections" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here, we are showing some basic functionalities for analysing MD trajectories using `Projection` classes. First, let's import HTMD:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "2024-06-11 18:08:42,913 - numexpr.utils - INFO - Note: NumExpr detected 20 cores but \"NUMEXPR_MAX_THREADS\" not set, so enforcing safe limit of 8.\n", "2024-06-11 18:08:42,914 - numexpr.utils - INFO - NumExpr defaulting to 8 threads.\n", "2024-06-11 18:08:43,027 - rdkit - INFO - Enabling RDKit 2022.09.1 jupyter extensions\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "\n", "Please cite HTMD: Doerr et al.(2016)JCTC,12,1845. https://dx.doi.org/10.1021/acs.jctc.6b00049\n", "HTMD Documentation at: https://software.acellera.com/htmd/\n", "\n", "You are on the latest HTMD version (2.3.28+7.g46a4ca83e.dirty).\n", "\n" ] } ], "source": [ "from htmd.ui import *\n", "from moleculekit.config import config\n", "from matplotlib import pylab as plt\n", "\n", "config(viewer='webgl')" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Projecting simulations" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "-" } }, "source": [ "Get the data for this tutorial [here](http://pub.htmd.org/tutorials/projections/data.tar.gz). Alternatively, you can download the data using `wget`:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import os\n", "assert os.system('wget -rcN -np -nH -q --cut-dirs=2 -R index.html* http://pub.htmd.org/tutorials/projections/data/') == 0" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "-" } }, "source": [ "We will be doing two types of data projections:\n", "* Projecting single Molecule objects\n", "* Projecting lists of simulations" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Projecting single Molecule objects" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One can load a previously simulated single MD trajectory (in this case, NTL9 protein) into a `Molecule` object and visualize it using:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "2024-06-11 18:09:06,033 - moleculekit.molecule - INFO - Removed 17 atoms. 631 atoms remaining in the molecule.\n" ] }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "28b3ab2842774d8cafec010edd0afd56", "version_major": 2, "version_minor": 0 }, "text/plain": [] }, "metadata": {}, "output_type": "display_data" }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "218b94aa4cd94b84960cca840bc8873d", "version_major": 2, "version_minor": 0 }, "text/plain": [ "NGLWidget(max_frame=9)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "mol = Molecule('data/ntl9_structure.psf')\n", "mol.read('data/ntl9_trajectory.xtc')\n", "mol.filter('protein')\n", "mol.view()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Getting the RMSD along time" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "HTMD has many projection classes. One of those classes is `MetricRmsd`." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using the `MetricRmsd` class, one may calculate the RMSD of a `Molecule` object against a reference structure (`refmol`, e.g. the crystal structure) for a given atom selection (`trajrmsdstr`):" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "scrolled": true, "slideshow": { "slide_type": "-" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[16.360521 15.793981 14.122061 14.7720175 13.908776 12.928344\n", " 11.541591 10.691926 10.778459 12.866117 ]\n" ] } ], "source": [ "crystal = Molecule('data/ntl9_crystal.pdb')\n", "metr = MetricRmsd(refmol=crystal, trajrmsdstr='protein and name CA')\n", "proj = metr.project(mol)\n", "print(proj)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The values are in Angstroms. They are quite high values, because the sample trajectory is of the unfolded protein (see the visualization above)." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Plot the RMSD along time" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The previously obtained RMSD values can be easily plotted using matplotlib:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "slideshow": { "slide_type": "-" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAGwCAYAAABcnuQpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAABStUlEQVR4nO3dd1hTZ8MG8PskgbCjbJmCgog4EFHqqHvV7efe2mqtW1tb7bB2ql221rZq6+hw1q1ttW7c4sC9UKaIKMqGMHK+PxD6UkcJJJwk3L/rynWZnENy43n75vY5z3mOIIqiCCIiIiIjJZM6ABEREVFFsMwQERGRUWOZISIiIqPGMkNERERGjWWGiIiIjBrLDBERERk1lhkiIiIyagqpA+ibRqNBYmIibG1tIQiC1HGIiIioDERRREZGBtzc3CCTPX/sxeTLTGJiIjw9PaWOQUREROUQHx8PDw+P5+5j8mXG1tYWQNFfhp2dncRpiIiIqCzS09Ph6elZ8j3+PCZfZopPLdnZ2bHMEBERGZmyTBHhBGAiIiIyaiwzREREZNRYZoiIiMioscwQERGRUWOZISIiIqPGMkNERERGjWWGiIiIjBrLDBERERk1lhkiIiIyaiwzREREZNRYZoiIiMioscwQERGRUWOZKSdRFLH/2j2Ioih1FCIioiqNZaac1p6Kx5hVpzF6VQSS0nKljkNERFRlscyUU36hBuYKGQ5ev4+OCw9h45kEjtIQERFJgGWmnEY2r4k/p7REQ89qyMgtwBu/n8fLP5/GvXSO0hAREVUmlpkKqO1si03jX8CbXerAXC7D/mvJ6PjVIWw+y1EaIiKiysIyU0EKuQwT2tTGziktUd9dhfTcAszYcB5jfzmD5AyO0hAREekby4yO+LvYYsuE5pjZuQ7M5AL2Xr2HTgvDsS3yDkdpiIiI9IhlRocUchkmtq2NHZNbIsjdDqnZ+Zi6LhLjfzuD+xlqqeMRERGZJJYZPQhwtcOWCS0wo6M/FDIBuy/fQ6eFh7DjfCJHaYiIiHSMZUZPzOQyTGnvh+2TWiKwhh0eZedj8tpzmLD6LB5kcpSGiIhIV1hm9CzQzQ5bJ7bA1PZ+UMgE/HUpCZ0WhuOPC3eljkZERGQSWGYqgblChukd/bF1YgsEuNriYVYeJq45i4lrzuJhVp7U8YiIiIway0wlCnJXYfuklpjSrjbkMgF/XLiLTgsPYdcljtIQERGVF8tMJTNXyDCjUx1sndAC/i42eJCZh/G/ncWUtefwiKM0REREWmOZkUh9DxV2TG6JiW1rQSYA288nouPCcPx9OUnqaEREREaFZUZCSoUcMzsHYMuEFvBztsGDTDXG/XoG09dHIjWbozRERERlwTJjABp6VsOOyS0xvnXRKM2Wc3fQaWE49l29J3U0IiIig8cyYyAszOSY1TUAm15rDl8nayRnqPHyz6fx+obzSMvJlzoeERGRwWKZMTDBXtXx55RWGPeiLwQB2HQ2AZ0WHsKBa8lSRyMiIjJILDMGyMJMjrdfqouN41+Aj6M17qWrMXpVBGb+fh7puRylISIi+l+Slpnw8HD06NEDbm5uEAQBW7dufWKfq1evomfPnlCpVLC1tUVYWBji4uIqP6wEQrzt8eeUVni5pQ8EAfj9TAI6LwzHoRv3pY5GRERkMCQtM1lZWWjYsCEWL1781O23bt1Cy5YtERAQgIMHD+L8+fN47733YGFhUclJpWNpLsd73QOx4dUXUNPBCnfTcjFyxSnM2nQBGRylISIigiAayG2cBUHAli1b0Lt375LXBg0aBDMzM/z6669lfh+1Wg21+p8bOaanp8PT0xNpaWmws7PTZeRKl5NXiM92X8PKozEAADeVBRb0a4BWfk7SBiMiItKx9PR0qFSqMn1/G+ycGY1Ggz/++AP+/v7o3LkznJ2d0axZs6eeivpf8+bNg0qlKnl4enpWTuBKYGkux/s96mHduDB42VshMS0Xw5efwttbLiJTXSB1PCIiIkkYbJlJTk5GZmYm5s+fjy5duuDvv/9Gnz590LdvXxw6dOiZPzd79mykpaWVPOLj4ysxdeUI83XArmmtMPIFbwDAmpNx6LwwHEejHkicjIiIqPIppA7wLBqNBgDQq1cvTJ8+HQDQqFEjHDt2DEuWLEHr1q2f+nNKpRJKpbLSckrFylyBD3oFoXOQK97ceAEJj3Iw9KeTGB7mjVldA2CtNNhDS0REpFMGOzLj6OgIhUKBwMDAUq/XrVu3ylzNVBbNazli97QXMSzMCwDw64lYdPkmHMdvpUicjIiIqHIYbJkxNzdHaGgorl+/Xur1GzduwNvbW6JUhslaqcDHvetj9SvN4F7NEvEPczD4xxN4f9slZOdxLg0REZk2Sc9FZGZmIioqquR5dHQ0IiMjYW9vDy8vL8ycORMDBw7Eiy++iLZt22LXrl3YsWMHDh48KF1oA9aitiN2TWuFT/+8hrWn4vDz8VgcuH4fn/drgGa+DlLHIyIi0gtJL80+ePAg2rZt+8TrI0eOxKpVqwAAK1aswLx585CQkIA6derggw8+QK9evcr8Gdpc2mVKwm/cx6xNF5CYlgtBAEY1r4k3OwfA0lwudTQiIqL/pM33t8GsM6MvVbXMAEB6bj4+/eMq1kUUXdFV08EKn/dviNCa9hInIyIiej6TWGeGKs7Owgzz/68BVo0OhaudBWJSsjFg6XF8vPMKcvMLpY5HRESkEywzVUCbOs7YPf1F9A/xgCgCPx2JxkvfHMaZ2EdSRyMiIqowlpkqQmVphs/7N8TKUaFwsVPi9oMs9F9yDPP+vMpRGiIiMmosM1VM2wBn/D2tNfo2dodGBJaG30a3RYdxLSld6mhERETlwjJTBamszPDVgEb4aUQTONkqcet+FiavOYdCjUnPBSciIhPFMlOFdQh0wa6prWBnocDN5Ez8cfGu1JGIiIi0xjJTxTnYKDG2lS8A4Ju9Nzg6Q0RERodlhjCqRU1UszLDrftZ2HE+Ueo4REREWmGZIdhamP0zOrPvJgoKNRInIiIiKjuWGQIAjGxeE9WtzBD9IAvbIjk6Q0RExoNlhgAANkoFxr1YCwCwaD9HZ4iIyHiwzFCJES94w8HaHLEp2dh87o7UcYiIiMqEZYZKWCsVeLV10dyZb/ffRD5HZ4iIyAiwzFApw8K84WhjjviHOdh0JkHqOERERP+JZYZKsTJXYHzrorkz3+6PQl4BR2eIiMiwsczQE4aFecPJVok7qTn4/Uy81HGIiIiei2WGnmBhJseENkWjM9/tj4K6gHfVJiIiw8UyQ081uKkXXOyUSEzLxYbTnDtDRESGi2WGnsrCTI6JbWsDKBqdyc3n6AwRERkmlhl6poGhnqihskBSei7WR3DuDBERGSaWGXompUKOCcWjMwc4OkNERIaJZYaea0ATD7hXs0RyhhprTsZJHYeIiOgJLDP0XErFP3Nnfjh0Czl5HJ0hIiLDwjJD/6lfiAc8qlvifoYaq0/GSh2HiIioFJYZ+k/mChkmtysanVly6Bay8wokTkRERPQPlhkqk76NPeBlb4UHmXn49ThHZ4iIyHCwzFCZmMn/GZ1ZGn4bWWqOzhARkWFgmaEy6xPsjpoOVniYlYefj8dIHYeIiAgAywxpQSGXYUp7PwDAsvDbyOToDBERGQCWGdJKz4Zu8HW0Rmp2Pn4+FiN1HCIiIpYZ0o5CLsPUDv+MzqTn5kuciIiIqjqWGdJa9wZuqO1sg7ScfKw6GiN1HCIiquJYZkhrcplQMnfmx8O3kZbD0RkiIpIOywyVS7f6NeDvYoOM3AKsOBItdRwiIqrCWGaoXOQyAVPb+wMAVhyJRlo2R2eIiEgaLDNUbl2DXBHgaosMdQF+OnJb6jhERFRFscxQuclkAqY9vrJp5dEYPMrKkzgRERFVRSwzVCGdAl1Rt4YdMtUF+PEwR2eIiKjyscxQhchkAqY/Hp1ZdSwGDzk6Q0RElYxlhiqsY6ALgtztkJ1XiKXht6SOQ0REVQzLDFWYIAiY3qHoyqZfjsXiQaZa4kRERFSVsMyQTrQLcEZDDxVy8guxLJxzZ4iIqPKwzJBOCIKAaR0fj84cj0FyRq7EiYiIqKpgmSGdaePvhEae1ZCbr8HSQxydISKiysEyQzojCAKmPx6d+e1ELJLTOTpDRET6xzJDOvWinyNCvKtDXaDB9wd5ZRMREekfywzp1P9e2bTmVByS0jg6Q0RE+sUyQzrXorYDmta0R16BBt8fjJI6DhERmThJy0x4eDh69OgBNzc3CIKArVu3lto+atQoCIJQ6hEWFiZNWCqzoiubilYFXncqHompORInIiIiUyZpmcnKykLDhg2xePHiZ+7TpUsX3L17t+Tx559/VmJCKq/mtRzRzMceeYUafHeAozNERKQ/Cik/vGvXrujatetz91EqlXB1da2kRKRL0zv6Y9CyE9hwOh6vtakFj+pWUkciIiITZPBzZg4ePAhnZ2f4+/tj7NixSE5Ofu7+arUa6enppR4kjTBfBzSv5YD8QpGjM0REpDcGXWa6du2K1atXY//+/fjyyy8RERGBdu3aQa1+9r1/5s2bB5VKVfLw9PSsxMT0b8Xrzvx+OgHxD7MlTkNERKZIEEVRlDoEUDRpdMuWLejdu/cz97l79y68vb2xbt069O3b96n7qNXqUmUnPT0dnp6eSEtLg52dna5jUxkMX34Sh28+wMAmnljQr4HUcYiIyAikp6dDpVKV6fvboEdm/q1GjRrw9vbGzZs3n7mPUqmEnZ1dqQdJa9rjdWc2nk1AbEqWxGmIiMjUGFWZSUlJQXx8PGrUqCF1FNJCiHd1tPZ3QqFGxLf7OXeGiIh0S9Iyk5mZicjISERGRgIAoqOjERkZibi4OGRmZuKNN97A8ePHERMTg4MHD6JHjx5wdHREnz59pIxN5VA8d2bz2QREP+DoDBER6Y6kZeb06dMIDg5GcHAwAGDGjBkIDg7GnDlzIJfLcfHiRfTq1Qv+/v4YOXIk/P39cfz4cdja2koZm8qhkWc1tAtwhkYEvt337NOERERE2jKYCcD6os0EItKvCwmp6Ln4KGQCsGdGa9RyspE6EhERGSiTnQBMxq2BRzV0qOsCjQgs4ugMERHpCMsMVappHYru2bT9fCKikjMkTkNERKaAZYYqVZC7Cp0CXSCKwNd7OTpDREQVxzJDla543Zk/Lt7F9SSOzhARUcWwzFClC3SzQ9cgV4gi8M2+G1LHISIiI8cyQ5KY1sEfggD8eTEJV+/yZqDFzsenous3hzFh9Rnk5hdKHYeIyCiwzJAk6rja4qX6RSs5f8O5MxBFEb+eiEX/Jcdx9W46/ryYhHG/stAQEZUFywxJZlp7PwgCsOtyEi4npkkdRzLZeQWYvj4S7229hLxCDVr5OcLSTI7wG/cxYfVZ5BVopI5IRGTQWGZIMn4utujRwA1A1b2yKSo5E70WH8XWyETIZQLeeakufhnTFMtHNYFSIcP+a8mYtOYs8gtZaIiInoVlhiQ1pb1f0YrAV+7hYkLVGp3ZeSERvRYfwc3kTDjbKrF2bBjGvugLQRDQvJYjfhrZBOYKGf6+cg9T151DAQsNEdFTscyQpGo726BXI3cAwNd7q8aVTXkFGszdfhmT1pxDVl4hwnztsXNKSzT1sS+1Xys/JywdHgJzuQx/XkzC9A3nUagx6buPEBGVC8sMSW5yu9qQCcC+a8mIjE+VOo5eJabmYOCy41h1LAYAMKFNLfz2cjM421o8df+2dZzx/dDGUMgE7DifiJm/s9AQEf0bywxJztfJBn2CPQCY9ujM4Zv30f3bIzgXlwo7CwV+GtEEb3YJgEL+/P8MOwS6YPGQYMhlAjafu4NZmy5Aw0JDRFSCZYYMwpT2tSGXCTh4/T7Oxj2SOo5OaTQiFu27iRErTuFhVh6C3O2wc3IrdAh0KfN7dAmqgW8GNYJMAH4/k4B3tl6Cid/wnoiozFhmyCB4O1ijb3DR3JmFe0xndOZRVh5Gr4rAV3tuQBSBwU09sXF8c3g5WGn9Xt0buGHhwEYQBGDtqTi8v/0yCw0REVhmyIBMbucHhUzA4ZsPcDrmodRxKiwyPhXdvz2CQzfuw8JMhi/6N8S8vg1gYSYv93v2auSOz/s1hCAAvxyPxUc7r7LQEFGVxzJDBsPLwQr9Qormziw04rkzoiji1+Mx6L/kGO6k5sDH0RpbJrQo+d0qql+IB+b3rQ8AWHE0GvP/usZCQ0RVGssMGZSJbWvDTC7gaFQKTt5OkTqO1rLUBZi2PhLvbbuM/EIRXeq5YtukFqhbw06nnzMw1Asf9w4CACwNv40v/77BQkNEVRbLDBkUT3sr9G/iCcD4RmeikjPQ+7uj2PZ4Nd93u9XFD8Maw87CTC+fNyzMG3N7BAIAFh+Iwjf7quYqykRELDNkcIpHZ07cfojjt4xjdGbH+UT0XHy0ZDXfdePC8EqrotV89WlUCx+8260ugKJbQnx3IEqvn0dEZIhYZsjguFezxKBQLwBFozOGfPqkeDXfyWvPITuvEC/4OuCPKa0QWtP+v39YR15p5Yu3ugQAAD7ffR3Lwm9V2mcTERkClhkySBPa1oK5XIZT0Q9xzEBHZxJTczBgaenVfH99uSmcbJWVnuW1NrUwo6M/AODTP69h+ZHoSs9ARCQVlhkySDVUlhjS7PHozB7DG50Jv3Ef3RYdRmR80Wq+y0eWbTVffZrS3g9T2tUGAHy08wp+OR4jWRYiosrEMkMG67U2taBUyHA69hEO33wgdRwARav5fr33BkauPIVH2fkIcrfDH1NaoX3dsq/mq0/TO/rjtTa1AABztl3GmpNxEiciItI/lhkyWC52FhjazBuAYcydeZiVh1GrIvD13puPV/P1wsbxzeFpr/1qvvoiCALe7FwHr7T0AQC8veUiNpyOlzgVEZF+scyQQRvfxhcWZjKci0vFoRv3JctxLu4Rui86jPDHq/l+2b8h5vWtX6HVfPVFEAS8060uRjWvCQB4a9MFbDmXIG0oIiI9Ypkhg+Zsa4FhxaMzEsydEUURPx+LwYClx5GYlgsfR2tsndgC/6ej1Xz1RRAEvN8jEMPCvCCKwOsbzmP7+USpYxER6YWiLDstWrSozG84ZcqUcocheppXW9fC6pNxOJ+QhgPXk9EuoHLmp2SpCzBr80XseFwCuga54rN+DWCrp0XwdE0QBHzYMwgFhSLWRcRj+vpIKGQCXqpfQ+poREQ6JYhl+Keuj49P2d5MEHD79u0Kh9Kl9PR0qFQqpKWlwc5Ot0vKU+WZ9+dVLA2/jfruKmyf1ELvi9FFJWdg/G9nEZWcCYVMwKyuAXi5pY/eP1cfNBoRMzdewKazCVDIBHw/tDE61XOVOhYR0XNp8/1dpjJjzFhmTENKphqtPjuA7LxC/DiiCToG6m90Zvv5RMzadAHZeYVwsVPiuyGN0aQSF8HTh0KNiBkbIrEtMhFmcgFLh4dU2ggXEVF5aPP9zTkzZBQcbJQY+XhCq77mzuQVaPD+tkuY8ng13+a1ilbzNfYiAwBymYAv+zdEtwY1kF8oYvyvZyWdUE1EpEtlmjPzbwkJCdi+fTvi4uKQl5dXattXX32lk2BE/zaulS9+ORaDK3fTsfvyPXQJ0t2pkjupOZi4+iwi41MBAJPa1sb0jv6Qy4zvtNKzKOQyfD2wEQoLRey6nIRxv5zGilGhaFHbUepoREQVonWZ2bdvH3r27AkfHx9cv34dQUFBiImJgSiKaNy4sT4yEgEAqlubY3QLHyw+EIWv995Ap0AXyHRQNg7duI9p687hUXY+7CwUWDiwkcEsgqdrZnIZFg0OxoTVZ7D3ajJe/jkCq0Y3RZivg9TRiIjKTevTTLNnz8brr7+OS5cuwcLCAps2bUJ8fDxat26N/v376yMjUYlXWvnAVqnAtaQM7LqcVKH3KtSIWLjnBkYZ6Gq++mKukOG7oY3Rpo4TcvM1GLMqAhExD6WORURUblqXmatXr2LkyJEAAIVCgZycHNjY2ODDDz/EggULdB6Q6H9VszLH6Mer23699wY0mvLNnXmYlYdRK0/hm31Fq/kOaWZ4q/nqk1Ihx5JhIWjl54jsvEKMXhmBs3GPpI5FRFQuWpcZa2trqNVqAICbmxtu3bpVsu3BA8O4fw6Ztpdb+sDWQoEb9zLxx8W7Wv988Wq+h28+gIWZDF8NaIhP+xjmar76ZGEmx7LhTfCCrwMy1QUYufwULiSkSh2LiEhrWpeZsLAwHD16FADQrVs3vP766/jkk08wZswYhIWF6Twg0b+pLM3wSktfAMA3+26isIyjM6IoYtXR6CdW8+3b2LBX89UnS3M5lo9qgqY17ZGhLsCwn07i0p00qWMREWlF6zLz1VdfoVmzZgCAuXPnomPHjli/fj28vb2xfPlynQckeprRLWvCzkKBqORM7Lzw38v0Z6kLMGVdJObuuIL8QhEv1XfF9kktEODKtYeszBVYMToUId7VkZ5bgGHLT+Lq3XSpYxERlRkXzSOj9e2+m/hyzw34Olljz/TWz7yM+ua9DIz/7Qxu3c+CQiZg9kt1MaZFTaNczVef0nPzMXz5KZyPT4WDtTnWjQuDn4ut1LGIqIrS66J5vr6+SElJeeL11NRU+Pr6avt2ROU2qkVNVLMyw+37Wdh+/s5T99kWeQe9vjuKW/ez4GKnxLpxYUZ7WwJ9s7Mwwy9jmiLI3Q4pWXkY/ONJRCVnSh2LiOg/aV1mYmJiUFhY+MTrarUad+48/QuFSB9sLcwwttXjuTN7b6KgUFOyTV1QiPe2XsLUdZHIzitEi9qms5qvPqkszfDby81Qt4YdHmSqMeTHE4h+kCV1LCKi5yrzonnbt28v+fPu3buhUqlKnhcWFmLfvn2oWbOmTsMR/ZeRzWti+ZFoxKRkY2tkIvqFeCDhUTYmrjmH849X853crjamdTCt1Xz1qZqVOVa/0gyDl53A9XsZGPLjCawf9wK8HKrGZetEZHzKPGdGJisaxBEE4Yn74piZmaFmzZr48ssv0b17d92nrADOmTF9Sw7dwvy/rsHL3grv9wjE67+fR2p2PlSWZvh6YCO0DXCWOqJRepCpxqBlJxCVnAn3apZY/2oYPKqz0BBR5dDrXbN9fHwQEREBR0fjuJ8Ly4zpy84rQKsFB5CS9c99whp4qPDdkMZVZhE8fUlOz8WgZSdw+0EWPO0tseHVF1BDZSl1LCKqAvQ6ATg6OvqJIpOamqrt2xDpjJW5AuNb1yp5PrSZF34f/wKLjA4421lgzdgweDtYIf5hDgYvO4F76blSxyIiKkXrMrNgwQKsX7++5Hn//v1hb28Pd3d3nD9/XqfhiMpqZPOaeL2jP5YMC8EnfepDqahaq/nqk6uqqNB4VLdETEo2Bv94AskZLDREZDi0LjNLly6Fp6cnAGDPnj3Yu3cvdu3aha5du2LmzJk6D0hUFuYKGSa390OXIFepo5gk92qWWDs2DG4qC9y+n4WhP55ESqZa6lhERADKUWbu3r1bUmZ27tyJAQMGoFOnTnjzzTcRERGh84BEZBg87a2wdlwYXO0scDM5E0N/OolH/zNPiYhIKlqXmerVqyM+Ph4AsGvXLnTo0AFA0X1vnrb+zPOEh4ejR48ecHNzgyAI2Lp16zP3ffXVVyEIAr7++mttIxORjng7WGPN2GZwslXiWlIGhi0/ibTsfKljEVEVp3WZ6du3L4YMGYKOHTsiJSUFXbt2BQBERkaidu3aWr1XVlYWGjZsiMWLFz93v61bt+LkyZNwc3PTNi4R6Zivkw3Wjm0GRxtzXE5Mx/AVJ5Gey0JDRNLRuswsXLgQkyZNQmBgIPbs2QMbGxsARaefJkyYoNV7de3aFR9//DH69u37zH3u3LmDSZMmYfXq1TAzM9M2LhHpQW1nW6x+JQz21ua4kJCGkStOIYOFhogkUuYVgIvl5eXhjTfeeOL1adOm6SJPKRqNBsOHD8fMmTNRr169Mv2MWq2GWv3PxMT0dN79l0gf6rja4reXm2HwjydwLi4Vo1dG4OcxTWGt1Pr/VoiIKkTrkRkXFxeMGTMGR44c0UeeUhYsWACFQoEpU6aU+WfmzZsHlUpV8iierExEuhfoZoffXm4GWwsFTsc+wphVEcjJ027uHBFRRWldZtauXYu0tDS0b98e/v7+mD9/PhITE3Ue7MyZM/jmm2+watUqre5wPHv2bKSlpZU8iicrE5F+1PdQ4deXm8FWqcDJ6Id45ZcI5Oaz0BBR5dG6zPTo0QObNm1CYmIiXnvtNaxduxbe3t7o3r07Nm/ejIKCAp0EO3z4MJKTk+Hl5QWFQgGFQoHY2Fi8/vrrz72hpVKphJ2dXakHEelXI89qWDUmFNbmchyNSsG4X8+w0BBRpdG6zBRzcHDA9OnTcf78eXz11VfYu3cv+vXrBzc3N8yZMwfZ2dkVCjZ8+HBcuHABkZGRJQ83NzfMnDkTu3fvrtB7E5HuhXjbY+XoprA0kyP8xn1MWH0WeQUaqWMRURVQ7pl6SUlJ+OWXX7By5UrExcWhX79+ePnll5GYmIj58+fjxIkT+Pvvv5/7HpmZmYiKiip5Hh0djcjISNjb28PLywsODg6l9jczM4Orqyvq1KlT3thEpEdNfeyxfFQTjFkVgf3XkjF57Vl8PzQEclnZTxUTEWlL6zKzefNmrFy5Ert370ZgYCAmTpyIYcOGoVq1aiX7NGrUCMHBwf/5XqdPn0bbtm1Lns+YMQMAMHLkSKxatUrbaERkAJrXcsRPI0Ix5ucI7L58Dx/tvIK5Pct2NSIRUXloXWZGjx6NwYMH4+jRowgNDX3qPr6+vnjnnXf+873atGkDURTL/NkxMTFl3peIpNPSzxFfD2yECavPYtWxGHjZW2FMSx+pYxGRiRJELdpEQUEBli1bhr59+8LV1Thu6Jeeng6VSoW0tDROBiaqZEsP3cK8v65BEIClw0LQqZ5x/P8GEUlPm+9vrSYAKxQKvPHGG6UWpSMiepZxL/picFMviCIwdV0kLiakSR2JiEyQ1lczNWvWDOfOndNHFiIyMYIg4KNe9fCivxNy8gsx5ucIJDyq2JWORET/pvWcmQkTJuD1119HQkICQkJCYG1tXWp7gwYNdBaOiIyfQi7Dd0OC0X/JcVxLysCYVRHY+Fpz2FnwXmtEpBtazZkBAJnsycEcQRAgiiIEQUBhoWEtlMU5M0SGITE1B72/O4rkDDVa+TlixahQmMnLvdQVEZk4bb6/tR6ZiY6OLncwIqq63KpZYsWoUAxYehyHbz7Ae1svYV7f+lrdroSI6Gm0LjPe3t76yEFEVUCQuwrfDg7G2F9OY11EPLwcrDChTW2pYxGRkdN6jHfevHlYsWLFE6+vWLECCxYs0EkoIjJd7eu6YE73QADAZ7uuY+cF3d+oloiqFq3LzNKlSxEQEPDE6/Xq1cOSJUt0EoqITNuoFj4Y3aImAGDGhvM4E/tQ2kBEZNS0LjNJSUmoUaPGE687OTnh7t27OglFRKbv3W6B6BjogrwCDcb+cgaxKVlSRyIiI6V1mfH09MTRo0efeP3o0aNwc3PTSSgiMn1ymYBvBjVCfXcVHmblYfTKCDzKypM6FhEZIa3LzCuvvIJp06Zh5cqViI2NRWxsLFasWIHp06dj7Nix+shIRCbKylyB5SObwL2aJW4/yMKrv52BusCwlncgIsOn9Tozoihi1qxZWLRoEfLyiv4VZWFhgbfeegtz5szRS8iK4DozRIbvelIG+v1wDBnqAvQJdsdXAxrykm2iKk6b72+ty0yxzMxMXL16FZaWlvDz84NSqSxXWH1jmSEyDodv3seolREo1IiY0t4PMzr6Sx2JiCSktxtN/i8bGxuEhoYiKCjIYIsMERmPVn5O+KR3EABg0b6b2HQmQeJERGQsuJY4ERmMQU298FqbWgCAWZsv4PitFIkTEZExYJkhIoMys1MddG9QA/mFIl799TSikjOkjkREBo5lhogMikwm4Iv+DRHiXR3puQUYvSoCDzLVUsciIgPGMkNEBsfCTI5lw0PgZW+F+Ic5GPvLaeTm85JtInq6cpWZW7duYfLkyejQoQM6duyIKVOm4NatW7rORkRVmIONEitHh0JlaYZzcamYvj4SGk25Lr4kIhOndZnZvXs3AgMDcerUKTRo0ABBQUE4efIk6tWrhz179ugjIxFVUbWcbLBseAjM5AL+upSEBbuuSR2JiAyQ1uvMBAcHo3Pnzpg/f36p12fNmoW///4bZ8+e1WnAiuI6M0TGb8u5BExffx4A8EmfIAxt5i1xIiLSN72uM3P16lW8/PLLT7w+ZswYXLlyRdu3IyL6T32CPTC9Q9EienO2XcbB68kSJyIiQ6J1mXFyckJkZOQTr0dGRsLZ2VkXmYiInjClfW30beyOQo2ISWvO4erddKkjEZGBUGj7A2PHjsW4ceNw+/ZtNG/eHIIg4MiRI1iwYAFef/11fWQkIoIgCJjftwESU3Nw4vZDjFkVga0TW8DFzkLqaEQksXLdaPLrr7/Gl19+icTERACAm5sbZs6ciSlTphjczeE4Z4bItKRl56PvD0dx634W6rnZYcOrL8BaqfW/y4jIwFXKjSYBICOjaGVOW1vb8r6F3rHMEJmeuJRs9Pn+KFKy8tA+wBnLRjSBXGZY/5AioorR6wTgdu3aITU1FUBRiSkuMunp6WjXrp32aYmItOTlYIUfRzaBUiHDvmvJ+GgnLz4gqsq0LjMHDx5EXl7eE6/n5ubi8OHDOglFRPRfGntVx8KBjQAAq47FYMWRaGkDEZFkynyi+cKFCyV/vnLlCpKSkkqeFxYWYteuXXB3d9dtOiKi53ipfg3M7hqAeX9dw0d/XIFHdUt0qucqdSwiqmRlLjONGjWCIAgQBOGpp5MsLS3x7bff6jQcEdF/GfeiL2IfZmPNyThMXReJ9a+GoYFHNaljEVElKnOZiY6OhiiK8PX1xalTp+Dk5FSyzdzcHM7OzpDL5XoJSUT0LIIg4MOe9ZDwKAfhN+7j5Z9PY8uE5vCobiV1NCKqJBW6mskY8GomoqohIzcf/Zccx7WkDPi72GDja81hZ2EmdSwiKie9Xs1ERGSIbC3MsGJUKJxtlbhxLxMTfjuL/EKN1LGIqBKwzBCRyXCrZokVo0JhZS7HkagHeHfLJZj44DMRgWWGiExMkLsK3w4OhkwA1p+Oxw+HbkkdiYj0jGWGiExO+7oueL9HPQDAZ7uuY8f5RIkTEZE+aX1DE1EUcebMGcTExEAQBPj4+CA4ONjg7slERFXbyOY1EZuSjRVHo/H67+dRQ2WBJjXtpY5FRHqg1cjMgQMHUKtWLTRr1gwDBgxA//79ERoaCj8/P4SHh+srIxFRubzTrS46Brogr0CDsb+cRsyDLKkjEZEelLnMREVFoXv37qhZsyY2b96Mq1ev4sqVK/j999/h4eGBl156Cbdv39ZnViIirchlAr4Z1AgNPFR4lJ2PMasi8CjryduxEJFxK/M6M5MmTcLVq1exb9++J7aJoogOHTogMDDQ4FYB5jozRJSckYs+3x3DndQcNK1pj19faQqlgot8Ehkyvawzc/DgQUybNu2p2wRBwLRp03DgwAGtghIRVQZnWwusGBUKW6UCp2Ie4q2NF3jJNpEJKXOZiYuLQ/369Z+5PSgoCLGxsToJRUSka3VcbfH9sMZQyARsjUzEwr03pY5ERDpS5jKTmZkJK6tn3+vEysoK2dnZOglFRKQPrfyc8EmfIADAon03sfFMgsSJiEgXtLo0+8qVK0hKSnrqtgcPHugkEBGRPg0M9UJsSja+P3gLszdfgFs1CzSv5Sh1LCKqgDJPAJbJZBAE4annmYtfFwQBhYWFOg9ZEZwATET/ptGImLLuHHZeuAs7CwU2T2iO2s62Usciov+hzfd3mUdmoqOjKxyMiMgQyGQCvujfEHfTcnEm9hFGr4rAlgkt4GijlDoaEZVDmUdmjBVHZojoWVIy1ej7wzHEpmSjkWc1rBsXBgszXrJNZAj0cmn2w4cPkZBQerLc5cuXMXr0aAwYMABr1qzROmh4eDh69OgBNzc3CIKArVu3lto+d+5cBAQEwNraGtWrV0eHDh1w8uRJrT+HiOhpHGyUWDEqFCpLM0TGp2L6+khoNCb97zsik1TmMjNx4kR89dVXJc+Tk5PRqlUrREREQK1WY9SoUfj111+1+vCsrCw0bNgQixcvfup2f39/LF68GBcvXsSRI0dQs2ZNdOrUCffv39fqc4iInqWWkw2WDQ+BuVyGvy4lYcGua1JHIiItlfk0k4+PD1auXIk2bdoAAL744gssWbIE165dg0KhwBdffIGNGzfixIkT5QsiCNiyZQt69+79zH2Kh5z27t2L9u3bl+l9eZqJiMpi67k7mLY+EgDwSZ8gDG3mLW0goipOL6eZkpKS4OPjU/J8//796NOnDxSKojnEPXv2xM2b+luEKi8vD8uWLYNKpULDhg2fuZ9arUZ6enqpBxHRf+kd7I7pHfwBAHO2XcbB68kSJyKisipzmbGzs0NqamrJ81OnTiEsLKzkuSAIUKvVOg0HADt37oSNjQ0sLCywcOFC7NmzB46Oz14TYt68eVCpVCUPT09PnWciItM0pX1t/F9jDxRqRExacw5X7/IfQ0TGoMxlpmnTpli0aBE0Gg02btyIjIwMtGvXrmT7jRs39FIc2rZti8jISBw7dgxdunTBgAEDkJz87H8xzZ49G2lpaSWP+Ph4nWciItMkCALm9a2PF3wdkKkuwJhVEbiXnit1LCL6D2UuMx999BG2bdsGS0tLDBw4EG+++SaqV69esn3dunVo3bq1zgNaW1ujdu3aCAsLw/Lly6FQKLB8+fJn7q9UKmFnZ1fqQURUVuYKGZYMC0EtJ2vcTcvFmFURyFIXSB2LiJ6jzIvmNWrUCFevXsWxY8fg6uqKZs2aldo+aNAgBAYG6jzgv4miqJfTWURExVRWZlg5qin6fH8UlxPTMWXtOfw4oglkMkHqaET0FFrdm8nJyQm9evV66rZu3bpp/eGZmZmIiooqeR4dHY3IyEjY29vDwcEBn3zyCXr27IkaNWogJSUF33//PRISEtC/f3+tP4uISBteDlb4aWQTDFp2AvuuJWPjmQQMCOUcPCJDVOYy88svv5RpvxEjRpT5w0+fPo22bduWPJ8xYwYAYOTIkSWXff/888948OABHBwcEBoaisOHD6NevXpl/gwiovIK9qqONzrVwSd/XsW8v66iY6ALqlubSx2LiP5FqxtN2tjYQKFQPPVmk0DR5LmHDx/qNGBFcZ0ZIqqI/EINui86guv3MjC4qSfm9W0gdSQig5GlLkBWXgGcbS10/t56WWembt26MDc3x4gRI3Do0CE8evToiYehFRkioooyk8vwcZ8gAMDaU/E4G/dI4kREhmPjmQQ0n7cf8/+SduXsMpeZy5cv448//kBOTg5efPFFNGnSBD/88AMXpSMikxda0x79QjwAAO9uuYSCQo3EiYikJ4oiVp+MRYFGhKudtHecL3OZAYBmzZph6dKluHv3LqZMmYINGzagRo0aGDp0KK8wIiKTNrtrAFSWZrhyNx2/noiVOg6R5E7HPsKNe5mwMJOhT2MPSbNoVWaKWVpaYsSIEfjggw/QtGlTrFu3DtnZ2brORkRkMBxslHizSx0AwJd/30AyF9OjKm7NyTgAQM+GblBZmkmaResyc+fOHXz66afw8/PDoEGDEBoaisuXL5daQI+IyBQNCvVCQ89qyFQX4OM/rkodh0gyD7Py8MfFuwCAIQZwU9Yyl5kNGzaga9eu8PPzQ0REBL788kvEx8fjs88+Q0BAgD4zEhEZBLlMwMe9giATgO3nE3E06oHUkYgkselMAvIKNKjnZoeGHiqp42h3abaXlxeGDh0KFxeXZ+43ZcoUnYXTBV6aTUS69v62S/j5eCx8nazx19RWUCrkUkciqjSiKKLdl4cQ/SALn/apjyHNvPTyOdp8f5d50TwvLy8IgoA1a9Y8cx9BEAyuzBAR6dqMTnXwx8Uk3L6fhZ8OR2Ni29pSRyKqNMdvpSD6QRaszeXo2chN6jgAtCgzMTExeoxBRGQ8VJZmeLdbXUxbH4lF+26iZ0M3eNpbSR2LqFKsfjzxt3ewO2yUWt0VSW/KdTXTs9y5c0eXb0dEZLB6NXLDC74OUBdo8MGOy1LHIaoU9zPU2H05CQAw1AAm/hbTSZlJSkrC5MmTUbs2h1qJqGoQBAEf9a4HM7mAvVeTsefKPakjEendhtPxKNCIaORZDYFuhjMPtcxlJjU1FUOHDoWTkxPc3NywaNEiaDQazJkzB76+vjhx4gRWrFihz6xERAaltrMtxrbyBQDM3X4Z2XkFEici0h+NRsTaU0WnmIbqadJveZW5zLz99tsIDw/HyJEjYW9vj+nTp6N79+44cuQI/vrrL0RERGDw4MH6zEpEZHAmt/ODezVL3EnNweL9UVLHIdKb8Jv3kfAoB3YWCnRvYBgTf4uVucz88ccfWLlyJb744gts374doijC398f+/fvR+vWrfWZkYjIYFmayzG3Zz0AwI+HbyMqOUPiRET6UTzxt29jD1iaG9ZyBGUuM4mJiQgMDAQA+Pr6wsLCAq+88oreghERGYuOgS7oUNcZ+YUi3tt6GWVcvovIaNxNy8H+a8kADO8UE6BFmdFoNDAz++feC3K5HNbW1noJRURkbN7vUQ8WZjIcv52C7ecTpY5DpFPrI+JRqBHR1Mcefi62Usd5QpkvEBdFEaNGjYJSWXSb79zcXIwfP/6JQrN582bdJiQiMgKe9laY3M4Pn+++jo92XkXbAGfYWUh78z0iXSgo1GDdqXgAhjkqA2hRZkaOHFnq+bBhw3QehojImL3Sygebzibg9v0sfPX3jZK5NETG7MD1+0hKz4W9tTm6BLlKHeepylxmVq5cqc8cRERGT6mQ46NeQRj600n8cjwG/UI8EOQu/U34iCpi9clYAEC/EA+DvQ+ZTlcAJiKq6lrUdkTPhm7QiMA7Wy9Bo+FkYDJe8Q+zcejGfQDA4KaGeYoJYJkhItK5d7vVhY1SgfPxqVgXES91HKJyWxcRB1EEWtZ2hI+j4V70wzJDRKRjznYWeL2TPwBgwa5rSMlUS5yISHt5BRqsj0gAAAwx0Im/xVhmiIj0YHiYNwJr2CEtJx/z/7omdRwire25cg8PMtVwslWiY6CL1HGei2WGiEgPFHIZPu4TBAD4/UwCImIeSpyISDtrThVN/B3YxBNmcsOuC4adjojIiDX2qo7BTT0BAO9uuYT8Qo3EiYjK5vb9TByNSoEgAIMe/2/YkLHMEBHp0ZudA1DdygzX72Vg1dEYqeMQlUnx3bHb+DvBo7qVxGn+G8sMEZEeVbc2x+yudQEAC/fewN20HIkTET1fbn4hNp4pmvg7tJm3xGnKhmWGiEjP+oV4IMS7OrLzCvHRzitSxyF6rl2XkvAoOx81VBZoU8dJ6jhlwjJDRKRnMpmAj3sHQS4T8OfFJBy8nix1JKJnWnOy6BTToFAvKAx84m8x40hJRGTk6taww6jmNQEA72+/jNz8QmkDET3FjXsZOBXzEHKZgIGhhj/xtxjLDBFRJZnWwQ8udkrEpmRjyaFbUschekLxqEz7AGe4qiwkTlN2LDNERJXE1sIM73UPBAB8f/AWYh5kSZyI6B85eYXYdPbxxN8w45j4W4xlhoioEnWrXwOt/ByRV6DBnO2XIYq8ESUZhh0XEpGRWwBPe0u0qu0odRytsMwQEVUiQRDwQc96MJfLEH7jPnZdSpI6EhEAYPXjU0yDm3pBJhMkTqMdlhkiokrm62SD8a19AQAf7LiCTHWBxImoqrt0Jw3n41NhJhfQP8R4Jv4WY5khIpLAhLa14WlviaT0XCzad1PqOFTFrXm84m/neq5wslVKnEZ7LDNERBKwMJPjw55FN6JcfiQa15MyJE5EVVWmugDbzt0BAAxp5iVxmvJhmSEikkjbAGd0rueCQo2Id7de5GRgksS2yDvIyiuEr6M1XvB1kDpOubDMEBFJaE6PerA0kyMi5hE2nb0jdRyqYkRRxG8nik4xDWnmBUEwrom/xVhmiIgk5F7NElM7+AEA5v15FanZeRInoqokMj4VV++mw1whw/819pA6TrmxzBARSWxMCx/4OdsgJSsPn+++LnUcqkKKV/ztXr8GqlubS5ym/FhmiIgkZq6Q4aPeRZOB15yKQ2R8qrSBqEpIy8nHjguJAIx34m8xlhkiIgMQ5uuAvsHuEEXg3a0XUajhZGDSry1nE5Cbr0EdF1uEeFeXOk6FsMwQERmI2S/VhZ2FApfupGP1yVip45AJE0WxZMXfoWHGO/G3GMsMEZGBcLJVYmaXAADA57uvIzkjV+JEZKoiYh7hZnImLM3k6B3sLnWcCmOZISIyIEOaeqGBhwoZuQWY9+c1qeOQiVrzeOSvZ0M32FmYSZym4lhmiIgMiFwm4OPeQRAEYMu5Ozh+K0XqSGRiHmbl4c+LRTc4HRpm3BN/i7HMEBEZmAYe1TCsmTcA4L1tl5BXoJE4EZmSjWfikVeoQZC7HRp4VJM6jk6wzBARGaA3OtWBo405opIzsfxItNRxyERoNCLWnooHAAx9XJhNgaRlJjw8HD169ICbmxsEQcDWrVtLtuXn5+Ott95C/fr1YW1tDTc3N4wYMQKJiYnSBSYiqiQqKzO8/VJdAMCifTeR8Chb4kRkCo7fTkH0gyzYKBXo2dBN6jg6I2mZycrKQsOGDbF48eIntmVnZ+Ps2bN47733cPbsWWzevBk3btxAz549JUhKRFT5+gS7o6mPPXLyC/HBjitSxyETUHzJf+9gN1grFRKn0R1Jf5OuXbuia9euT92mUqmwZ8+eUq99++23aNq0KeLi4uDl9fRJS2q1Gmq1uuR5enq67gITEVUiQSiaDPzSN4ex58o97L1yDx0CXaSORUYqOSMXf1++BwAY0tR0TjEBRjZnJi0tDYIgoFq1as/cZ968eVCpVCUPT0/PygtIRKRj/i62eLmVDwBg7o7LyMkrlDgRGavfTyegQCOisVc1BLrZSR1Hp4ymzOTm5mLWrFkYMmQI7OyefRBmz56NtLS0kkd8fHwlpiQi0r0p7fzgprJAwqMcfHcgSuo4ZIQKNWLJTSWHmNDE32JGUWby8/MxaNAgaDQafP/998/dV6lUws7OrtSDiMiYWSsVmNOjHgBgafgt3LqfKXEiMjbhN+/jTmoO7CwU6N6ghtRxdM7gy0x+fj4GDBiA6Oho7Nmzh+WEiKqkzvVc0LaOE/ILRczZdgmiyBtRUtmtPlE0KtMvxBMWZnKJ0+ieQZeZ4iJz8+ZN7N27Fw4ODlJHIiKShCAI+KBnEJQKGY5GpWDHhbtSRyIjcTctB/uvPZ7428w055FKWmYyMzMRGRmJyMhIAEB0dDQiIyMRFxeHgoIC9OvXD6dPn8bq1atRWFiIpKQkJCUlIS8vT8rYRESS8HKwwsS2tQEAH+28gvTcfIkTkTFYdyoeGhFo5mOP2s62UsfRC0nLzOnTpxEcHIzg4GAAwIwZMxAcHIw5c+YgISEB27dvR0JCAho1aoQaNWqUPI4dOyZlbCIiyYx70Rc+jta4n6HGwj03pI5DBq6gUIN1EUWnmIaGmd7E32KSrjPTpk2b55735TlhIqLSLMzk+KBnPYxYcQo/H4tBvxAP1HNTSR2LDNT+a8m4l66GvbU5Otcz3TWKDHrODBERPelFfyd0a1ADGhF4d+slaDT8hx893erHl2P3b+IBpcL0Jv4WY5khIjJC73ULhLW5HOfiUrHhNNfToifFP8xG+M37AIDBoU9fNd9UsMwQERkhV5UFpnf0BwDM33UND7N4YQSVtvZUHEQRaOXniJqO1lLH0SuWGSIiIzWqeU0EuNoiNTsfC/66JnUcMiB5BZqSEbuhzUx7VAZgmSEiMloKuQwf9w4CAKw/HY8zsQ8lTkSG4u8rSXiQmQcnWyXa1zXdib/FWGaIiIxYk5r2GNDEAwDwzpZLKCjUSJyIDEHxfZgGhXrCTG76X/Wm/xsSEZm4WV3ropqVGa4lZeDn47FSxyGJ3bqfiWO3UiATgEFNTf8UE8AyQ0Rk9OytzfFWlwAAwFd/X0dSWq7EiUhKax+PyrSp4wz3apYSp6kcLDNERCZgYBNPNPKshqy8Qnz8xxWp45BEcvMLsfFsAoCqMfG3GMsMEZEJkMkEfNw7CDIB2HnhLg4/Xl+Eqpa/Lt1FanY+3FQWaFPHWeo4lYZlhojIRAS5qzDihZoAgDnbLkNdUChtIKp0JRN/m3pBLhMkTlN5WGaIiEzIjE7+cLJVIvpBFpYdui11HKpE15MyEBHzCHKZgIGhnlLHqVQsM0REJsTOwgzvdqsLAFh8IApxKdkSJ6LKsuZk0ZVsHeu6wMXOQuI0lYtlhojIxPRs6IbmtRygLtDg/e2XIIq8EaWpy84rwOZzdwAAQ6rQxN9iLDNERCZGEAR82CsIZnIBB67fx99X7kkdifRs5/m7yMgtgJe9FVrWdpQ6TqVjmSEiMkG1nW3w6ou1AAAfbL+M7LwCiRORPq1+fIppSDMvyKrQxN9iLDNERCZqYtva8KhuicS0XCzaFyV1HNKTS3fScD4hDWZyAf1CPKSOIwmWGSIiE2VpLscHPesBAH46fBs372VInIj0YfXjy7G7BNWAo41S4jTSYJkhIjJh7eu6oGOgCwo0It7degkaDScDm5KM3Hxsiyya+FuVVvz9N5YZIiIT936PQFiYyXAy+iHm8Oomk7ItMhHZeYWo5WSNZj72UseRDMsMEZGJ86huhfl9G0AQgN9OxOG9bSw0pkAUxZJTTEOaeUMQqt7E32IsM0REVUDvYHd83q9hSaF5f/tlFhojdy4+FVfvpsNcIcP/NXaXOo6kWGaIiKqIfiEe+Oz/ikZofjkei7ksNEat+D5M3RvUQDUrc4nTSItlhoioCunfxBMLHhean4/H4oMdV1hojFBadj52nE8EAAxt5i1xGumxzBARVTEDmnhiQd8GAIBVx2Lw4U4WGmOz+VwC1AUaBLjaorFXNanjSI5lhoioChoQ6okF/1cfALDyaAw+2nmVhcZI/O/E36HNvKr0xN9iLDNERFXUwFAvzO9bVGhWHI3Gx3+w0BiDU9EPEZWcCStzOXoHV+2Jv8VYZoiIqrBBTb3waZ+iQrP8SDQ+YaExeGtOFY3K9GzoBlsLM4nTGAaWGSKiKm5IMy980icIAPDTkWjM++saC42BSslU46+LSQA48fd/scwQERGGNvPGx72LCs2y8NuYz0JjkDaeSUBeoQYNPFSo76GSOo7BYJkhIiIAwLAwb3z0uNAsDb+N+btYaAyJRiNi7eNTTEOaVt37MD0NywwREZUYHuaND3sV3Wl76aHb+Gz3dRYaA3HsVgpiUrJhq1SgR0M3qeMYFJYZIiIqZcQLNUsKzQ8Hb+FzFhqDsPpkLACgT2N3WCsVEqcxLCwzRET0hBEv1MTcHoEAgO8P3sIXf7PQSCk5PRd7rtwDUDRhm0pjmSEioqca1cIH7z8uNN8duIWv9txgoZHIhtPxKNCICPGujgBXO6njGByWGSIieqbRLXwwp3tRofl2fxQWstBUukKNiLWn4gEUrfhLT2KZISKi5xrT0gfvdqsLAFi0Pwpf770pcaKqJfzGfdxJzYHK0gwv1a8hdRyDxDJDRET/6ZVWviWF5pt9N/H13hsSJ6o6iif+9gvxgIWZXOI0hollhoiIyuSVVr5456WiQvP13pv4hiM0epeYmoP915IBcOLv87DMEBFRmY190RdvvxQAAFi49wYW7WOh0ad1EfHQiECYrz1qOdlIHcdgscwQEZFWxr1YC7O6FhWar/bcwOL9LDT6UFCowfqIohV/eR+m52OZISIirY1vXQtvdSkqNF/8fQPfHYiSOJHp2XctGffS1XCwNkfneq5SxzFoLDNERFQur7WphTe71AEAfL77OguNjq0+WTQq07+JJ8wV/Lp+Hv7tEBFRuU1oUxszO/9TaH44eEviRKYhLiUb4TfuA+BNJcuCZYaIiCpkYtvaeKOTPwBgwa5rWHKIhaai1j6eK9PKzxFeDlYSpzF8LDNERFRhk9r54fWORYVm/l/XsCychaa88go02BBRvOIvJ/6WBcsMERHpxOT2fpjeoajQfPrnNfwYflviRMZp9+UkpGTlwcVOifZ1naWOYxRYZoiISGemdvDDtA5+AIBP/ryKnw6z0GhrzeOJvwObeMJMzq/pspD0byk8PBw9evSAm5sbBEHA1q1bS23fvHkzOnfuDEdHRwiCgMjISElyEhFR2U3r4I+p7YsKzcd/sNBo49b9TBy/nQKZAAzkxN8yk7TMZGVloWHDhli8ePEzt7do0QLz58+v5GRERFQR0zr4YUq72gCKCs3yI9ESJzIOax+PyrQLcIZ7NUuJ0xgPhZQf3rVrV3Tt2vWZ24cPHw4AiImJqaRERESkC4IgYHpHf4gAvt0fhY92XoFMAEa38JE6msHKzS/ExrMJAHgfJm1JWmb0Qa1WQ61WlzxPT0+XMA0RUdUlCAJmdPSHKAKLD0Thgx1XIAAYxULzVH9evIvU7Hy4V7NEa39O/NWGyc0smjdvHlQqVcnD09NT6khERFWWIAh4vZM/JrSpBQCYu+MKfj4WI20oA1U88XdwU0/IZYLEaYyLyZWZ2bNnIy0treQRHx8vdSQioipNEATM7FwHrz0uNO9vv4xfjsdIG8rAXEtKx+nYR1DIBAxown+Ea8vkTjMplUoolUqpYxAR0f8QBAFvdq4DUQSWHLqFOdsuQxAEDA/jonDAP6MyHQNd4GxnIXEa42NyIzNERGSYBEHAW13q4NUXfQEA7229hN9OxEqcSnrZeQXYcvYOAK74W16SjsxkZmYiKuqfu6xGR0cjMjIS9vb28PLywsOHDxEXF4fExEQAwPXr1wEArq6ucHXl7dCJiIyNIAiY1TUAIoBl4bfx7tZLEISq/SW+43wiMtQF8HawQvNaDlLHMUqSjsycPn0awcHBCA4OBgDMmDEDwcHBmDNnDgBg+/btCA4ORrdu3QAAgwYNQnBwMJYsWSJZZiIiqhhBEDC7awDGtiq6qumdLZdKTrNURasf/+5DmnpBxom/5SKIoihKHUKf0tPToVKpkJaWBjs7O6njEBHRY6IollpQb17f+hhcxVa9vZiQhh6Lj8BcLsPx2e3gYMM5n8W0+f7mnBkiIpKEIAh4t1tdjHm87szszRexPqJqjdCsOVU0Z6hLkCuLTAWwzBARkWQEQcB73etidIuaAIBZmy9iQ0TVWFIjIzcf2yKL5oQO5Yq/FcIyQ0REkhIEAXO6B2JU85oQReCtzRew4bTpF5qtkYnIzitEbWcbNPWxlzqOUTO5dWaIiMj4CIKA93sEAgBWHYvBW5suQADQ30QXkBNFEasfX5Y+pKkXBIETfyuCZYaIiAxCcaHRiCJ+OR6LNzddgCAI6BfiIXU0nRFFEem5BTh5OwXXkjKgVMjwf41N5/eTCssMEREZDEEQ8EHPehBF4NcTsZi58TwEAP9nBIUmO68A99LVuJeei3vpuUh+/Oek4j9nFL2em68p+ZnuDdygsjKTMLVpYJkhIiKDIggCPuxVDyJE/HYiDm9sPA+ZDOgTLE2hySvQ4H6mGklpuUh+XFTuZfxTWooLTEZuQZnfs5qVGbztrTChbS09Jq86WGaIiMjgCIKAD3sGQRSLFpV7fcN5CBDQO9hdZ59RqBGRkqVGcnpRUSkaOVEj+fFoSvGfU7LyyvyeVuZyuNpZwNlOCRc7C7jYWcDZVglXVdGfXWyLtlmYyXX2exDLDBERGSiZTMBHvYKgEYG1p+IwY0MkBAHo1ej5hUYURaTl5ONeuvpxKSkeUSk9mnI/U41CTdnWjTWXy/6noCjhbFtUTlxVyscFpeh1G6WCk3klwDJDREQGSyYT8EnvIAAi1p6Kx/T1kcjOK0RNB+vSp3kyckuNqOQVaP7zvQFAJgBOtsrHIyhFhcTFzuKJ0ZXqVmYsKQaMZYaIiAxaUaGpD1EE1kXEY/bmi2X6uepWZiVlpLikOD8uKsXPHazNoZBzyTVjxzJDREQGTyYT8Gmf+rC1UGDLuTtQWf5TVJztik71FM1LKToF5GynhFLBeSlVBW80SURERAaHN5okIiKiKoNlhoiIiIwaywwREREZNZYZIiIiMmosM0RERGTUWGaIiIjIqLHMEBERkVFjmSEiIiKjxjJDRERERo1lhoiIiIwaywwREREZNZYZIiIiMmosM0RERGTUWGaIiIjIqCmkDqBvoigCKLqVOBERERmH4u/t4u/x5zH5MpORkQEA8PT0lDgJERERaSsjIwMqleq5+whiWSqPEdNoNEhMTIStrS0EQdDpe6enp8PT0xPx8fGws7PT6XuT9ng8DAuPh2Hh8TAsPB7/TRRFZGRkwM3NDTLZ82fFmPzIjEwmg4eHh14/w87Ojv9jNCA8HoaFx8Ow8HgYFh6P5/uvEZlinABMRERERo1lhoiIiIway0wFKJVKvP/++1AqlVJHIfB4GBoeD8PC42FYeDx0y+QnABMREZFp48gMERERGTWWGSIiIjJqLDNERERk1FhmiIiIyKixzJTT999/Dx8fH1hYWCAkJASHDx+WOpLJmTdvHkJDQ2FrawtnZ2f07t0b169fL7WPKIqYO3cu3NzcYGlpiTZt2uDy5cul9lGr1Zg8eTIcHR1hbW2Nnj17IiEhoTJ/FZM0b948CIKAadOmlbzG41H57ty5g2HDhsHBwQFWVlZo1KgRzpw5U7Kdx6TyFBQU4N1334WPjw8sLS3h6+uLDz/8EBqNpmQfHg89EUlr69atE83MzMQff/xRvHLlijh16lTR2tpajI2NlTqaSencubO4cuVK8dKlS2JkZKTYrVs30cvLS8zMzCzZZ/78+aKtra24adMm8eLFi+LAgQPFGjVqiOnp6SX7jB8/XnR3dxf37Nkjnj17Vmzbtq3YsGFDsaCgQIpfyyScOnVKrFmzptigQQNx6tSpJa/zeFSuhw8fit7e3uKoUaPEkydPitHR0eLevXvFqKiokn14TCrPxx9/LDo4OIg7d+4Uo6Ojxd9//120sbERv/7665J9eDz0g2WmHJo2bSqOHz++1GsBAQHirFmzJEpUNSQnJ4sAxEOHDomiKIoajUZ0dXUV58+fX7JPbm6uqFKpxCVLloiiKIqpqamimZmZuG7dupJ97ty5I8pkMnHXrl2V+wuYiIyMDNHPz0/cs2eP2Lp165Iyw+NR+d566y2xZcuWz9zOY1K5unXrJo4ZM6bUa3379hWHDRsmiiKPhz7xNJOW8vLycObMGXTq1KnU6506dcKxY8ckSlU1pKWlAQDs7e0BANHR0UhKSip1LJRKJVq3bl1yLM6cOYP8/PxS+7i5uSEoKIjHq5wmTpyIbt26oUOHDqVe5/GofNu3b0eTJk3Qv39/ODs7Izg4GD/++GPJdh6TytWyZUvs27cPN27cAACcP38eR44cwUsvvQSAx0OfTP5Gk7r24MEDFBYWwsXFpdTrLi4uSEpKkiiV6RNFETNmzEDLli0RFBQEACV/3087FrGxsSX7mJubo3r16k/sw+OlvXXr1uHs2bOIiIh4YhuPR+W7ffs2fvjhB8yYMQNvv/02Tp06hSlTpkCpVGLEiBE8JpXsrbfeQlpaGgICAiCXy1FYWIhPPvkEgwcPBsD/RvSJZaacBEEo9VwUxSdeI92ZNGkSLly4gCNHjjyxrTzHgsdLe/Hx8Zg6dSr+/vtvWFhYPHM/Ho/Ko9Fo0KRJE3z66acAgODgYFy+fBk//PADRowYUbIfj0nlWL9+PX777TesWbMG9erVQ2RkJKZNmwY3NzeMHDmyZD8eD93jaSYtOTo6Qi6XP9GQk5OTn2jbpBuTJ0/G9u3bceDAAXh4eJS87urqCgDPPRaurq7Iy8vDo0ePnrkPlc2ZM2eQnJyMkJAQKBQKKBQKHDp0CIsWLYJCoSj5++TxqDw1atRAYGBgqdfq1q2LuLg4APxvpLLNnDkTs2bNwqBBg1C/fn0MHz4c06dPx7x58wDweOgTy4yWzM3NERISgj179pR6fc+ePWjevLlEqUyTKIqYNGkSNm/ejP3798PHx6fUdh8fH7i6upY6Fnl5eTh06FDJsQgJCYGZmVmpfe7evYtLly7xeGmpffv2uHjxIiIjI0seTZo0wdChQxEZGQlfX18ej0rWokWLJ5YruHHjBry9vQHwv5HKlp2dDZms9NeqXC4vuTSbx0OPJJp4bNSKL81evny5eOXKFXHatGmitbW1GBMTI3U0k/Laa6+JKpVKPHjwoHj37t2SR3Z2dsk+8+fPF1Uqlbh582bx4sWL4uDBg596maOHh4e4d+9e8ezZs2K7du14maOO/O/VTKLI41HZTp06JSoUCvGTTz4Rb968Ka5evVq0srISf/vtt5J9eEwqz8iRI0V3d/eSS7M3b94sOjo6im+++WbJPjwe+sEyU07fffed6O3tLZqbm4uNGzcuuVyYdAfAUx8rV64s2Uej0Yjvv/++6OrqKiqVSvHFF18UL168WOp9cnJyxEmTJon29vaipaWl2L17dzEuLq6SfxvT9O8yw+NR+Xbs2CEGBQWJSqVSDAgIEJctW1ZqO49J5UlPTxenTp0qenl5iRYWFqKvr6/4zjvviGq1umQfHg/9EERRFKUcGSIiIiKqCM6ZISIiIqPGMkNERERGjWWGiIiIjBrLDBERERk1lhkiIiIyaiwzREREZNRYZoiIiMioscwQERGRUWOZISK9mjt3Lho1aiTZ57/33nsYN26czt4vNDQUmzdv1tn7EVHFcQVgIio3QRCeu33kyJFYvHgx1Go1HBwcKinVP+7duwc/Pz9cuHABNWvW1Ml7bt++HW+88QauXbv2xE0FiUgaLDNEVG5JSUklf16/fj3mzJlT6i7OlpaWUKlUUkQDAHz66ac4dOgQdu/erbP3LCwshJubG1atWoWuXbvq7H2JqPz4zwoiKjdXV9eSh0qlgiAIT7z279NMo0aNQu/evfHpp5/CxcUF1apVwwcffICCggLMnDkT9vb28PDwwIoVK0p91p07dzBw4EBUr14dDg4O6NWrF2JiYp6bb926dejZs2ep19q0aYMpU6bgzTffhL29PVxdXTF37txS+8ydOxdeXl5QKpVwc3PDlClTSrbJ5XK89NJLWLt2bbn+zohI91hmiKjS7d+/H4mJiQgPD8dXX32FuXPnonv37qhevTpOnjyJ8ePHY/z48YiPjwcAZGdno23btrCxsUF4eDiOHDkCGxsbdOnSBXl5eU/9jEePHuHSpUto0qTJE9t+/vlnWFtb4+TJk/jss8/w4YcfYs+ePQCAjRs3YuHChVi6dClu3ryJrVu3on79+qV+vmnTpjh8+LCO/1aIqLxYZoio0tnb22PRokWoU6cOxowZgzp16iA7Oxtvv/02/Pz8MHv2bJibm+Po0aMAikZYZDIZfvrpJ9SvXx9169bFypUrERcXh4MHDz71M2JjYyGKItzc3J7Y1qBBA7z//vvw8/PDiBEj0KRJE+zbtw8AEBcXB1dXV3To0AFeXl5o2rQpxo4dW+rn3d3dERcXB41Go9u/GCIqF5YZIqp09erVKzV51sXFpdToh1wuh4ODA5KTkwEAZ86cQVRUFGxtbWFjYwMbGxvY29sjNzcXt27deupn5OTkAAAsLCye2NagQYNSz2vUqFHyWf3790dOTg58fX0xduxYbNmyBQUFBaX2t7S0hEajgVqtLsdvT0S6ppA6ABFVPWZmZqWeC4Lw1NeKRz40Gg1CQkKwevXqJ97LycnpqZ/h6OgIoOh007/3ed5neXp64vr169izZw/27t2LCRMm4PPPP8ehQ4dKfu7hw4ewsrKCpaVlWX9lItIjlhkiMniNGzfG+vXr4ezsDDs7uzL9TK1atWBnZ4crV67A399fq8+ztLREz5490bNnT0ycOBEBAQG4ePEiGjduDAC4dOlSyZ+JSHo8zUREBm/o0KFwdHREr169cPjwYURHR+PQoUOYOnUqEhISnvozMpkMHTp0wJEjR7T6rFWrVmH58uW4dOkSbt++jV9//RWWlpbw9vYu2efw4cPo1KlThX4nItIdlhkiMnhWVlYIDw+Hl5cX+vbti7p162LMmDHIycl57kjNuHHjsG7dOq0m6larVg0//vgjWrRogQYNGmDfvn3YsWNHyaJ/d+7cwbFjxzB69OgK/15EpBtcNI+ITJYoiggLC8O0adMwePBgnbznzJkzkZaWhmXLlunk/Yio4jgyQ0QmSxAELFu27ImrkSrC2dkZH330kc7ej4gqjiMzREREZNQ4MkNERERGjWWGiIiIjBrLDBERERk1lhkiIiIyaiwzREREZNRYZoiIiMioscwQERGRUWOZISIiIqPGMkNERERG7f8BfGGiReP51BgAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "time_range = mol.fstep*np.arange(len(proj))\n", "plt.plot(time_range, proj)\n", "_ = plt.xlabel('Time (ns)')\n", "_ = plt.ylabel('RMSD to crystal')" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Projecting a list of MD simulations" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "With the arrival of Markov state models, the paradigm of running one single long simulation has shifted to running hundreds or thousands of shorter simulations." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "With the analysis of a large amount of MD trajectories in mind, in HTMD, the `simlist` function can be used to bundle all the trajectories together:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Creating simlist: 100%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 3/3 [00:00<00:00, 2150.56it/s]\n" ] } ], "source": [ "sims = simlist(glob('data/1/filtered/*/'), 'data/1/filtered/filtered.pdb')" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Use Metric class to project a list of simulations" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using the `Metric` class, one can calculate all the projections of an entire simlist at once. In this case, the `MetricDistance` projection class is used, which calculates the matrix of distances between two selections (`sel1` and `sel2`):" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Projecting trajectories: 100%|███████████████████████████████████████████████████████████████████████████████████████████████████████████| 3/3 [00:00<00:00, 12.45it/s]\n", "2024-06-11 18:09:18,666 - htmd.projections.metric - INFO - Frame step 0.001ns was read from the trajectories. If it looks wrong, redefine it by manually setting the MetricData.fstep property.\n" ] } ], "source": [ "metr = Metric(sims)\n", "metr.set(MetricDistance(sel1='protein and name CA', sel2='resname MOL and noh', periodic='selections', metric='distances'))\n", "data = metr.project()\n", "data.fstep = 0.1" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here, the projection is the matrix of distances between each carbon alpha of the protein and each heavy atom of the MOL molecule (which in the case of these trajectories, it's a ligand binding to the protein)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "The `Metric.project` method returns a `MetricData` object which contains the trajectory and projection information:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "MetricData object with 3 trajectories of 60.0ns aggregate simulation time\n", "Centers: None\n", "K: None\n", "N: None\n", "description: at 0x78decc2d44f0\n", "fstep: 0.1\n", "parent: at 0x55ab3b8027e0\n", "trajectories shape: (3,)\n" ] } ], "source": [ "print(data)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "## Working with MetricData objects" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `MetricData` object contains the information on each trajectory, including the corresponding projection. For example, for simulation 2:" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "sim: simid = 2\n", "projection: np.array(shape=(200, 2493))\n", "reference: np.array(shape=(200, 2))\n", "cluster: None\n", "\n" ] } ], "source": [ "print(data.trajectories[2])" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "The projection data can be easily accessed as a whole in the `MetricData.dat` property, which is an array of size equal to the number of trajectories. We can check the projection data of simulation 14 with:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[46.903584, 48.459248, 48.862976, ..., 52.245136, 51.477478,\n", " 52.081757],\n", " [42.60737 , 42.6003 , 43.467403, ..., 46.8291 , 47.731323,\n", " 45.771404],\n", " [45.674603, 44.244358, 44.22837 , ..., 53.92176 , 54.58503 ,\n", " 53.017155],\n", " ...,\n", " [37.174904, 38.622143, 39.415646, ..., 37.552174, 37.359898,\n", " 36.812588],\n", " [41.236233, 42.631603, 43.36949 , ..., 41.442356, 41.54479 ,\n", " 41.618164],\n", " [45.716003, 45.90263 , 46.41945 , ..., 44.92147 , 44.872147,\n", " 44.011333]], dtype=float32)" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.dat[2] # equivalent to data.trajectories[2].projection" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(200, 2493)" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.dat[2].shape" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(600, 2493)\n" ] } ], "source": [ "bundledprojs = np.vstack(data.dat) # np.concatenate can be used instead as well\n", "print(bundledprojs.shape)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Operations and calculations" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One can use numpy functions like `np.min`, `np.max`, `np.average` to do calculations over the projected data:\n", "* get the global minimum over all data" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3.2580273\n" ] } ], "source": [ "global_minimum = np.min(bundledprojs)\n", "print(global_minimum)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "* get the maximum distance on each frame:" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(600,)\n" ] }, { "data": { "text/plain": [ "array([54.39825 , 53.998096, 55.228672, 51.56171 , 47.534378, 48.122726,\n", " 54.234592, 52.287495, 56.084633, 59.717773, 58.638866, 59.020008,\n", " 58.540077, 61.697826, 62.762337, 64.88375 , 57.164936, 61.46919 ,\n", " 61.785923, 51.648315, 50.93795 , 49.167656, 51.03158 , 49.82826 ,\n", " 50.78649 , 53.217636, 53.40074 , 54.59295 , 54.21768 , 49.040817,\n", " 49.194637, 51.379505, 53.340816, 52.982063, 53.229916, 52.800793,\n", " 51.121067, 50.52078 , 50.51849 , 53.35154 , 54.840214, 50.428257,\n", " 46.97517 , 48.814857, 46.214386, 47.91875 , 46.25085 , 46.454384,\n", " 48.536522, 47.89751 , 48.45051 , 47.181534, 50.3011 , 47.37468 ,\n", " 44.809082, 44.76846 , 44.396427, 49.58187 , 49.452393, 44.699467,\n", " 49.425236, 49.30479 , 48.831734, 48.04636 , 46.891056, 46.507465,\n", " 53.23867 , 52.209763, 51.905357, 52.760147, 53.15323 , 52.926044,\n", " 54.556217, 55.611897, 52.046055, 51.92997 , 50.818188, 51.01256 ,\n", " 49.769115, 49.608753, 50.520557, 52.815094, 51.304066, 49.99397 ,\n", " 50.81137 , 51.18169 , 49.10415 , 48.43612 , 52.17957 , 52.73203 ,\n", " 54.345753, 55.31595 , 55.39871 , 54.839466, 53.79522 , 53.284203,\n", " 47.44944 , 47.25443 , 46.913074, 48.135517, 46.80235 , 44.899788,\n", " 44.704002, 45.82361 , 46.498 , 46.07886 , 46.025085, 44.451492,\n", " 43.400375, 44.925617, 42.78121 , 42.454422, 44.038147, 44.447018,\n", " 43.457703, 44.942665, 48.482765, 50.187347, 51.073917, 49.428154,\n", " 48.50524 , 53.109505, 53.980877, 56.261894, 58.69202 , 58.06872 ,\n", " 58.765503, 60.141937, 56.34336 , 58.222004, 59.827217, 57.84348 ,\n", " 57.62419 , 59.43873 , 63.62348 , 64.050095, 62.193985, 59.517525,\n", " 58.832714, 61.67116 , 60.90755 , 58.634247, 60.395714, 60.144405,\n", " 49.775978, 48.653896, 49.283756, 56.077488, 49.26251 , 49.198658,\n", " 46.681004, 45.703247, 44.092865, 46.880257, 47.840847, 48.930866,\n", " 49.26276 , 48.424484, 51.754177, 53.38187 , 50.517876, 50.88899 ,\n", " 48.957874, 50.86398 , 49.52347 , 51.52822 , 48.15141 , 49.787586,\n", " 51.295555, 51.341866, 51.10505 , 49.70342 , 49.799862, 51.583054,\n", " 50.491463, 52.18646 , 51.79193 , 52.027943, 48.697643, 45.54545 ,\n", " 47.359547, 49.45975 , 48.112698, 52.51411 , 55.656136, 56.812218,\n", " 57.211945, 61.67162 , 60.73994 , 55.838486, 60.59322 , 52.707558,\n", " 59.200516, 55.77005 , 49.959007, 52.27606 , 53.4685 , 50.195663,\n", " 51.83982 , 46.972984, 58.02898 , 55.06521 , 56.92774 , 58.603672,\n", " 63.162193, 65.49683 , 55.11986 , 53.408955, 56.06107 , 51.011993,\n", " 54.8826 , 56.91488 , 58.846996, 59.789986, 60.955803, 62.344776,\n", " 61.62976 , 64.16956 , 58.623585, 61.526814, 58.12016 , 59.657608,\n", " 54.12708 , 51.00011 , 55.123035, 65.394104, 58.940674, 58.358627,\n", " 55.95053 , 56.449535, 52.173878, 55.778507, 57.703293, 57.1664 ,\n", " 51.501755, 48.076046, 50.35523 , 49.17984 , 47.47704 , 49.505035,\n", " 51.633377, 49.627056, 46.418175, 47.777798, 47.576702, 46.133038,\n", " 45.152557, 45.07051 , 47.970413, 48.24823 , 45.96343 , 45.115925,\n", " 45.666515, 46.490425, 45.549976, 46.01706 , 45.453964, 45.515446,\n", " 45.9844 , 47.13223 , 48.836105, 47.906868, 48.961525, 48.519596,\n", " 48.76764 , 51.081234, 46.477844, 46.278347, 44.705566, 44.669548,\n", " 44.54288 , 45.592293, 46.39919 , 47.019867, 45.584335, 45.207325,\n", " 46.80933 , 47.135063, 45.75573 , 45.835392, 47.149216, 46.10952 ,\n", " 47.89244 , 45.813976, 47.050255, 48.913574, 44.89156 , 47.54856 ,\n", " 45.633278, 48.73144 , 46.68292 , 46.73932 , 46.027905, 45.327396,\n", " 44.535587, 45.38553 , 47.58102 , 45.6692 , 45.515438, 43.811703,\n", " 44.08468 , 43.837803, 43.91667 , 44.465282, 44.301987, 44.174892,\n", " 46.657036, 48.582195, 47.261246, 50.408527, 49.713585, 47.64656 ,\n", " 51.807453, 48.009163, 54.441814, 52.480442, 59.65425 , 58.270264,\n", " 55.087585, 49.165928, 48.586834, 45.99378 , 48.621094, 49.504105,\n", " 50.48702 , 51.915627, 50.50157 , 51.82167 , 53.708004, 56.92613 ,\n", " 57.959114, 55.430756, 55.34785 , 49.523216, 48.99247 , 53.743862,\n", " 58.072803, 54.772854, 55.80709 , 55.54166 , 58.58327 , 60.215855,\n", " 57.643593, 58.334824, 58.690002, 58.64841 , 59.40299 , 58.352924,\n", " 55.96094 , 56.374454, 56.63217 , 62.827335, 63.720257, 62.8327 ,\n", " 65.87671 , 66.41088 , 64.80266 , 64.03413 , 67.14985 , 66.259865,\n", " 66.994804, 66.42933 , 65.72711 , 67.06555 , 65.66709 , 67.08758 ,\n", " 66.65906 , 65.61816 , 66.19741 , 67.06138 , 64.4134 , 57.828926,\n", " 55.1021 , 56.112644, 53.194065, 53.170574, 51.061966, 49.55039 ,\n", " 53.800117, 48.612778, 45.44221 , 47.68182 , 45.53881 , 44.353558,\n", " 52.942116, 54.64368 , 53.352097, 57.189117, 50.115234, 54.437584,\n", " 53.27835 , 51.984848, 45.55794 , 48.964207, 44.674267, 46.13538 ,\n", " 48.80334 , 45.493282, 52.49222 , 49.369164, 58.97841 , 56.06166 ,\n", " 59.335297, 54.887905, 54.104248, 51.059814, 48.11115 , 54.588306,\n", " 53.761787, 53.98229 , 52.87779 , 54.162224, 59.98114 , 61.48215 ,\n", " 60.23783 , 58.839363, 54.243195, 56.94387 , 55.312115, 59.235455,\n", " 55.783768, 51.72166 , 48.023804, 49.09402 , 55.95665 , 52.98509 ,\n", " 53.49044 , 55.68321 , 56.54712 , 54.37112 , 57.356087, 53.08747 ,\n", " 56.73367 , 57.498043, 53.299854, 53.063465, 49.98735 , 47.216633,\n", " 46.92734 , 57.599716, 57.037083, 55.949203, 57.110626, 65.150894,\n", " 63.827003, 64.43011 , 63.66299 , 66.6226 , 65.895065, 63.64519 ,\n", " 64.096016, 66.76584 , 64.84866 , 62.557846, 64.06469 , 62.43578 ,\n", " 64.79868 , 61.621136, 67.271355, 66.31896 , 67.36595 , 63.282326,\n", " 59.729267, 62.05453 , 58.38107 , 51.13966 , 53.981583, 48.66814 ,\n", " 54.832336, 54.26662 , 52.33797 , 51.29319 , 48.323696, 49.479324,\n", " 48.924763, 48.606667, 51.652077, 50.38462 , 47.688026, 48.65353 ,\n", " 48.816105, 46.393612, 46.657574, 44.796963, 45.534355, 48.16918 ,\n", " 48.30007 , 46.139935, 46.76868 , 45.927097, 45.686913, 45.812035,\n", " 48.05082 , 45.105473, 44.635887, 45.34747 , 46.082073, 44.79861 ,\n", " 45.04778 , 46.79814 , 46.277466, 46.9441 , 45.258972, 45.755707,\n", " 45.79995 , 45.970932, 44.646935, 48.003036, 46.15737 , 46.15521 ,\n", " 47.198174, 46.622158, 45.844547, 47.209805, 48.396137, 48.344986,\n", " 44.753998, 48.834183, 49.35003 , 51.633755, 49.904186, 49.072594,\n", " 49.796757, 50.38063 , 48.495564, 48.504204, 50.45795 , 49.709827,\n", " 49.515602, 48.371593, 50.092323, 51.048428, 50.520065, 49.9667 ,\n", " 52.122448, 49.67022 , 50.168163, 48.1803 , 49.34701 , 47.397907,\n", " 54.972626, 55.23825 , 54.018005, 51.48479 , 48.516758, 47.862297,\n", " 47.40332 , 58.441353, 51.811367, 52.852455, 59.88181 , 59.017944,\n", " 54.427086, 56.746006, 56.548817, 59.612827, 59.378536, 58.37478 ,\n", " 56.771206, 53.756557, 55.679623, 56.572567, 54.465122, 47.28183 ,\n", " 49.553474, 49.845093, 50.886894, 50.45243 , 53.55668 , 51.049065,\n", " 53.551304, 58.681187, 56.133167, 55.18882 , 56.56295 , 55.322258,\n", " 54.262825, 56.165165, 57.029556, 66.84226 , 60.031006, 60.488953,\n", " 61.970337, 55.85079 , 54.538193, 51.219833, 54.06704 , 53.40585 ,\n", " 53.578785, 47.85516 , 53.82843 , 53.730743, 54.80984 , 55.104183,\n", " 49.593353, 53.362053, 46.985603, 49.279034, 52.47066 , 48.883205],\n", " dtype=float32)" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "max_eachframe = np.max(bundledprojs,axis=1)\n", "print(max_eachframe.shape)\n", "max_eachframe" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "* or get the average of each pair across all frames:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(2493,)\n" ] }, { "data": { "text/plain": [ "array([34.71057 , 34.667583, 34.678036, ..., 39.51156 , 39.480606,\n", " 39.56731 ], dtype=float32)" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "avg_acrossframes = np.average(bundledprojs,axis=0)\n", "print(avg_acrossframes.shape)\n", "avg_acrossframes" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `axis` property is important to make the operation over the desired dimension:\n", "\n", "* `axis=1`, in this case, operates over all pairs, to give a result for each frame\n", "* `axis=0`, in this case, operates over all frames, to give a result for each pair\n", "\n", "But how can we know to what particular atoms of the system corresponds a given indexed pair?" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Playing with the projection mapping" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The mapping (meaning of each index) of the projection can be obtained from the `MetricData.map` property:\n", "\n", "* the mapping is a pandas `DataFrame` object\n", "\n", "The `head` method of `DataFrame` can be used to print the top entries of the mapping:" ] }, { "cell_type": "code", "execution_count": 18, "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", " \n", " \n", " \n", " \n", " \n", " \n", "
typeatomIndexesdescription
0distance[4, 4480]distance between GLU 1 CA and MOL 278 C3
1distance[4, 4483]distance between GLU 1 CA and MOL 278 C2
2distance[4, 4486]distance between GLU 1 CA and MOL 278 C1
3distance[4, 4489]distance between GLU 1 CA and MOL 278 C6
4distance[4, 4492]distance between GLU 1 CA and MOL 278 C5
\n", "
" ], "text/plain": [ " type atomIndexes description\n", "0 distance [4, 4480] distance between GLU 1 CA and MOL 278 C3\n", "1 distance [4, 4483] distance between GLU 1 CA and MOL 278 C2\n", "2 distance [4, 4486] distance between GLU 1 CA and MOL 278 C1\n", "3 distance [4, 4489] distance between GLU 1 CA and MOL 278 C6\n", "4 distance [4, 4492] distance between GLU 1 CA and MOL 278 C5" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.map.head()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "One can assess particular elements of the mapping through:\n", "\n", "* `iloc` - i(nteger)loc(ation):" ] }, { "cell_type": "code", "execution_count": 19, "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", "
typeatomIndexesdescription
20distance[29, 4486]distance between ASP 3 CA and MOL 278 C1
21distance[29, 4489]distance between ASP 3 CA and MOL 278 C6
\n", "
" ], "text/plain": [ " type atomIndexes description\n", "20 distance [29, 4486] distance between ASP 3 CA and MOL 278 C1\n", "21 distance [29, 4489] distance between ASP 3 CA and MOL 278 C6" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.map.iloc[20:22]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ " * `loc` (using the DataFrame index):" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "type distance\n", "atomIndexes [29, 4486]\n", "description distance between ASP 3 CA and MOL 278 C1\n", "Name: 20, dtype: object" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.map.loc[20]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "* or direct assess to the array elements:" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "'distance between ASP 3 CA and MOL 278 C1'" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.map['description'][20]" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "One can use the mapping to find out to what particular pair of elements the minimum distance found in the trajectories corresponds to." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `np.where` functionality gives us the array coordinates where a given condition is matched (in this case, when the distance is equal to the minimum):" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3.2580273\n", "[26]\n", "[3.2580273]\n" ] } ], "source": [ "print(global_minimum)\n", "frame, index = np.where(bundledprojs == global_minimum)\n", "print(index)\n", "print(bundledprojs[frame, index]) # confirm that indeed these coordinates correspond the inquired value" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then, using the index in the mapping shown before:" ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
typeatomIndexesdescription
26distance[29, 4500]distance between ASP 3 CA and MOL 278 N2
\n", "
" ], "text/plain": [ " type atomIndexes description\n", "26 distance [29, 4500] distance between ASP 3 CA and MOL 278 N2" ] }, "execution_count": 23, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data.map.iloc[index]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "gives us the atom indexes and description of the pair." ] } ], "metadata": { "anaconda-cloud": {}, "celltoolbar": "Slideshow", "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.10.13" } }, "nbformat": 4, "nbformat_minor": 4 }