{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Evaluating convergence\n", "\n", "Here we evaluate the convergence of a trajectory using the clustering ensemble similarity method and the dimensionality reduction ensemble similarity methods.\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.20.1\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n", "* [scikit-learn](https://scikit-learn.org/stable/)\n", " \n", "**Optional packages for visualisation:**\n", "\n", "* [matplotlib](https://matplotlib.org)\n", "\n", "\n", "
\n", " \n", "**Note**\n", "\n", "The metrics and methods in the `encore` module are from (Tiberti *et al.*, 2015). Please cite them when using the ``MDAnalysis.analysis.encore`` module in published work.\n", "\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 encore\n", "from MDAnalysis.analysis.encore.clustering import ClusteringMethod as clm\n", "from MDAnalysis.analysis.encore.dimensionality_reduction import DimensionalityReductionMethod as drm\n", "\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The test files we will be working with here feature adenylate kinase (AdK), a phosophotransferase enzyme. (Beckstein *et al.*, 2009)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "u = mda.Universe(PSF, DCD)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Evaluating convergence with similarity measures" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The convergence of the trajectory is evaluated by the similarity of the conformation ensembles in windows of the trajectory. The trajectory is divided into windows that increase by `window_size` frames. For example, if your trajectory had 13 frames and you specified a `window_size=3`, your windows would be:\n", "\n", " - Window 1: ---\n", " - Window 2: ------\n", " - Window 3: ---------\n", " - Window 4: -------------\n", " \n", "Where `-` represents 1 frame.\n", "\n", "These are compared using either the similarity of their clusters (`ces_convergence`) or their reduced dimension coordinates (`dres_convergence`). The rate at which the similarity values drop to 0 is indicative of how much the trajectory keeps on resampling the same regions of the conformational space, and therefore is the rate of convergence. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Using default arguments with clustering ensemble similarity" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "See [clustering_ensemble_similarity.ipynb](clustering_ensemble_similarity.ipynb#Calculating-clustering-similarity-with-default-settings) for an introduction to comparing trajectories via clustering." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "ces_conv = encore.ces_convergence(u, # universe\n", " 10, # window size\n", " selection='name CA') # default" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The output is an array of similarity values, with the shape (`number_of_windows`, `number_of_clustering_methods`)." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0.48194205],\n", " [0.40284672],\n", " [0.31699026],\n", " [0.25220447],\n", " [0.19829817],\n", " [0.14642725],\n", " [0.09911411],\n", " [0.05667391],\n", " [0. ]])" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ces_conv" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This can be easily plotted as a line." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, 'Jensen-Shannon divergence')" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEICAYAAABS0fM3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3hUddr/8fedRgkdAkrvSEeJNKkqKmsB14a9LqKoqKvPus9v3V133ebaXRQQKxbsHbvSQekdNIQWQQHB0CEJ9++PGdY8MZBJmZwk83ldVy4yZ86c+egFc8/5nu+5v+buiIhI7IoLOoCIiARLhUBEJMapEIiIxDgVAhGRGKdCICIS41QIRERiXFQLgZmdYWarzSzNzO7K5/mBZpZpZovCP3+MZh4REfmlhGgd2MzigTHAYCADmGtm77r7ijy7Tnf3syI9br169bx58+YlF1REJAbMnz9/m7un5Pdc1AoB0ANIc/d0ADObBAwF8haCQmnevDnz5s0rgXgiIrHDzNYf6bloDg01AjbmepwR3pZXbzNbbGYfmlnHKOYREZF8RPOMwPLZlrefxQKgmbvvNrNfAW8DbX5xILMRwAiApk2blnROEZGYFs0zggygSa7HjYFNuXdw953uvjv8+2Qg0czq5T2Qu49391R3T01JyXeIS0REiiiahWAu0MbMWphZEjAceDf3DmZ2jJlZ+Pce4Tw/RjGTiIjkEbWhIXfPNrObgI+BeOBpd19uZiPDz48FzgduMLNsYB8w3NUOVUSkVFl5+9xNTU11zRoSESkcM5vv7qn5Pac7i0VEYlzMFILtew5yz3vL2Z+VE3QUEZEyJWYKwYy0bTw7ax0XjZ/Dll37g44jIlJmxEwhOKdrQ8Ze1p1vvt/FuWNmsfr7XUFHEhEpE2KmEACc3vEYXr2+N1k5hzjviVlM/WZr0JFERAIXU4UAoHPjmrxz00k0qVOVa56dy8Q5R2y/ISISE2KuEAAcW7MKr43szYC2Kdz99jL++v4Kcg6Vr2m0IiIlJSYLAUC1Sgk8eUUqV/VpzlMz1nL9xHnsOZAddCwRkVIXs4UAID7O+PM5HbnnnI58sWoLF4ydzebMfUHHEhEpVTFdCA67sk9znrrqRNb/uIdhY2ay7LvMoCOJiJQaFYKwQe3q8/oNfYg344Kxs/l0xQ9BRxIRKRUqBLm0P7YGb486iTYNqjFi4jwmTE+nvPViEhEpLBWCPOrXqMwrI3pzeodjuPeDldz9zjKycw4FHUtEJGpUCPJRJSmexy89gZEDWvHCnA1c/excdu7PCjqWiEhUqBAcQVyccdeQ4/jnrzsze82PnP/ELDZu3xt0LBGREqdCUIDhPZry3DU92Jy5n3Mfn8nCDTuCjiQiUqJUCCJwUut6vHVjH6okxTN8/BwmL90cdCQRkRKjQhCh1vWr8/aNJ9GpUU1ufHEBY75M04wiEakQVAgKoW61Srx4XU/O6dqQf3+8mv95fQkHszWjSETKt6gtXl9RVU6M55Hh3WheL5lHP/+WjTv2Mvay7tSqmhR0NBGRItEZQRGYGbcPbstDF3Vlwfqf+PXjs1i3bU/QsUREikSFoBjOPb4xL1zXkx17D3Lu4zP5eu32oCOJiBSaCkEx9WhRh7duPInaVZO4bMJXvLUwI+hIIiKFokJQAprXS+bNG/twQrNa3PbKYh789BvNKBKRckOFoITUqprE89f05PzujXn0828ZPWkR+7Nygo4lIlIgzRoqQUkJcfz7/C60TEnmvo9Ws+mnfYy7vDt1q1UKOpqIyBHpjKCEmRk3DmzNmEtOYOl3mZz7+CzStuwOOpaIyBGpEETJmV2OZdKIXuw9mM25j89kZtq2oCOJiORLhSCKjm9am7duPIlja1bmyqe/5pW5G4KOJCLyCyoEUdakTlVev6EPvVvV5XdvLOWfH67i0CHNKBKRskOFoBTUqJzIM1edyKU9mzJ26hpGvbSAfQc1o0hEygYVglKSEB/HvcM6cfdZHfho+fcMHz+bLTv3Bx1LRCTyQmBmyYU9uJmdYWarzSzNzO46yn4nmlmOmZ1f2PcoT8yMa/u2YPzlqXzzw26GjZnJqu93Bh1LRGJcgYXAzPqY2QpgZfhxVzN7PILXxQNjgCFAB+BiM+twhP3+BXxcyOzl1uAODXhtZG9y3Dn/idl8uXpL0JFEJIZFckbwEHA68COAuy8G+kfwuh5Amrunu/tBYBIwNJ/9bgbeAGLq07BTo5q8M6ovzepW5dpn5/Lq3I1BRxKRGBXR0JC75/2UiuRKZyMg9+sywtv+y8waAecCYyPJUdEcU7Myr17fm5Na1+N/3ljC+Glrgo4kIjEokkKw0cz6AG5mSWZ2B+FhogJYPtvyzpt8GPidux+1sJjZCDObZ2bztm7dGsFblx/JlRJ46soTObPLsfx98ir+9dEqNawTkVIVSa+hkcAjhL7NZwCfAKMieF0G0CTX48bApjz7pAKTzAygHvArM8t297dz7+Tu44HxAKmpqRXuUzIpIY5Hhx9PzSqJPDFlDT/tzeLeYZ2Ij8uvloqIlKwCC4G7bwMuLcKx5wJtzKwF8B0wHLgkz7FbHP7dzJ4F3s9bBGJFfJzxt2GdqF01kTFfrmHnviwevKgrlRLig44mIhVcJLOGnjOzWrke1zazpwt6nbtnAzcRmg20EnjV3Zeb2UgzG1mc0BWVmXHn6cfxhzPb88HSzVz33Dz2HMgOOpaIVHBW0Hi0mS109+ML2lZaUlNTfd68eUG8dal6dd5G7npjCV2b1OKZq06kVtWkoCOJSDlmZvPdPTW/5yK5WBxnZrVzHawOWscg6i5MbcITl3Vn+aadXDhuNt9n6i5kEYmOSArBA8AsM/urmf0VmAXcF91YAnB6x2N49uoT+W7HPs4fO4t12/YEHUlEKqACC4G7Pw+cD/xA6KavX7v7xGgHk5A+rerx8ohe7D2Yw/ljZ7Nik1pSiEjJirTX0CrgTeAdYLeZNY1eJMmrS+NavHp9bxLjjYvGz2buuu1BRxKRCiSSWUM3Ezob+BR4H/gg/KeUotb1q/H6DX1IqV6Jy5/6ii9XxVRHDhGJokjOCEYD7dy9o7t3cffO7t4l2sHklxrVqsJr1/emTf3q/Ob5ebyz6LugI4lIBRBRiwkgM9pBJDJ1q1Xipd/0JLV5bW59ZRHPz14XdCQRKecimQaaDkwxsw+AA4c3uvuDUUslR1W9ciLPXt2Dm19eyB/fWc6OPVncckprwq06REQKJZIzgg2Erg8kAdVz/UiAKifG88SlJ3B+98Y89Nk33PPeCq2FLCJFEkmvoXsgtEKZu2siexmSEB/Hfed1oWaVRJ6asZbMfVncd34XEuO1AqmIRC6SWUO9i7JCmZSOuDjjD2e2587T2/HWwu8YOXE++7MiWS5CRCQkkq+OD1O0FcqklJgZowa15t5hnfhi9RaueOprdu7PCjqWiJQT0VyhTErZZb2a8ejw41m4cQfDx81h664DBb9IRGJeNFcokwCc3bUhE648kbXb9nDhuNls3L436EgiUsZFUghGElqR7PAKZd2IbIUyCciAtim8cF0Pftx9gAvGzubbH3YFHUlEyrBIms5tc/dL3b2Bu9d398vc/cfSCCdF171ZHV4d2ZtD7lwwbjaLNv4UdCQRKaMiWZjm0Xw2ZwLz3P2dqKQ6ilhZmKakbPhxL5c99RXbdh9g/OWp9G1TL+hIIhKA4i5MU5nQcNC34Z8uQB3gWjN7uMRSSlQ0rVuV10f2pmmdqlzz7Fw+XLo56EgiUsZEUghaAye7+2Pu/hhwKtAeOBc4LZrhpGTUr1GZV0b0pnPjmox6aQGTvt4QdCQRKUMiKQSNgORcj5OBhu6eQ67eQ1K21ayayMRre9CvTQp3vbmUsVPXBB1JRMqISArBfcAiM3vGzJ4FFgL3m1ky8Fk0w0nJqpqUwJNXpHJ214b888NV/OPDlRR0jUhEKr6j9hqyUDvLT4DJQA/AgP91903hXe6MbjwpaUkJcTx8UTdqVklg3NR0Mvdm8bdzOxMfp86lIrHqqIXA3d3M3nb37oSWqZQKID7O+OvQTtSpmsSjX6SRuS+Lh4d3o1JCfNDRRCQAkQwNzTGzE6OeREqVmXH7ae24+6wOfLjse659dh57DmQHHUtEAhBJIRhEqBisMbMlZrbUzJZEO5iUjmv7tuCBC7oyO/1HLpnwFTv2HAw6koiUskhWKBsS9RQSqPO6N6ZGlURGvbSAC8fNZuK1PTmmZuWgY4lIKYmkxcR6oAmhewnWA3sjeZ2UL4M7NOD5a3qwOXM/5z0xi7XbtAaRSKyIZGGaPwG/A34f3pQIvBDNUBKMXi3rMmlEL/Zn5XDB2Fks+y4z6EgiUgoi+WZ/LnAOsAcgPHVUaxZXUJ0a1eS1kb2plBDPxePnqCWFSAyIpBAc9NBdRw6htYujG0mC1jKlGq+N7E2LlGRueHEBv311Mbu04plIhRVJIXjVzMYBtczsN4TuJn4yurEkaA1rVeGNG/pwy8mteWthBkMemc7Xa7cHHUtEoiCSi8X3A68DbwDtgD+Gm89JBZcYH8ftp7XjtZF9iI8zLho/m399tIqD2YeCjiYiJSiSi8W3ASvd/U53v8PdP4304GZ2hpmtNrM0M7srn+eHhu9NWGRm88ysbyHzSyno3qw2k2/px0WpTXhiyhqGjZnJN1r1TKTCiGRoqAbwsZlNN7NRZtYgkgObWTwwhtB9CB2Ai82sQ57dPge6uns34BpgQuTRpTQlV0rgn+d14ckrUvlh537OemwGT89Yy6FDalonUt5FMjR0j7t3JLROcUNgqplF0nW0B5Dm7unufhCYBAzNc+zd/nP7y2TCF6Sl7BrcoQEf3dqffq3r8Zf3V3DF01/zfeb+oGOJSDEU5sawLcD3wI9A/Qj2bwRszPU4I7zt/zCzc81sFfABobMCKeNSqldiwpWp/P3czsxfv4PTH57Ge4s3FfxCESmTIrlGcIOZTSE0jFMP+I27d4ng2Pn1Nf7FN353f8vdjwOGAX89QoYR4WsI87Zu3RrBW0u0mRmX9GzK5NH9aFEvmZtfXsitkxaSuU/TTEXKm0jOCJoBt7p7R3f/k7uviPDYGYRaUxzWGDji10Z3nwa0MrNfrK7u7uPdPdXdU1NSUiJ8eykNLeol8/rI3tx2alveW7KZIQ9PY/aaH4OOJSKFcMRCYGY1wr/eB2wwszq5fyI49lygjZm1MLMkYDjwbp73aB1e/AYzOwFIIjT0JOVIQnwco09twxs39KFSYjyXTJjD3z5YwYHsnKCjiUgEjtZ99CXgLGA+oSGd3EM9DrQ82oHdPdvMbgI+BuKBp919uZmNDD8/FjgPuMLMsoB9wEW5Lh5LOdOtSS0+uKUvf5+8kienr2X6t9t4eHg3jjumRsEvFpHAWHn73E1NTfV58+YFHUMK8OWqLdz5+hJ27sviztPbcW3fFsRpOUyRwJjZfHdPzfe5IxWC8FDNEbn7ghLIVmgqBOXHj7sP8Ps3l/LJih/o1bIOD1zYjUa1qgQdSyQmFbUQfBn+tTKQCiwmNDzUBfjK3QO5C1iFoHxxd16bl8E97y0nLrxW8tBuDQlfGhKRUnK0QnDEi8XuPsjdBwHrgRPCs3a6A8cDadGJKhWNmXHhiU34cHR/2jaozq2vLOLmlxeSuVfTTEXKikimjx7n7ksPP3D3ZUC36EWSiqhp3aq8en1v7jy9HR8t+57TH57GjG+3BR1LRIisEKw0swlmNtDMBpjZk8DKaAeTiic+zhg1qDVv3XgSyZXiueypr7jnveXsz9I0U5EgRVIIrgaWA6OBW4EV4W0iRdK5cU3ev7kfV/VpzjMz13H2YzNYvknLYooERdNHJVBTv9nKna8tZsfeg9w+uB0j+rckXtNMRUpckS4Wi5SGAW1T+PjW/gzu0IB/fbSKi8fPYeP2vUHHEokpKgQSuNrJSYy55AQeuKArKzbvZMgj03l9fgbl7WxVpLxSIZAywcw4r3tjPhzdjw7H1uCO1xZz44sL2LHnYNDRRCq8o/UaAsDM2gJ3EupC+t/93f3kKOaSGNWkTlVeHtGLJ6en88Anq5m/fgf3nd+Fge0iWQJDRIqiwIvFZrYYGEuo+dx/5/m5+/zoRsufLhbHjuWbMrntlUV888NurujdjN8PaU+VpPigY4mUS0e7WFzgGQGQ7e5PlHAmkQJ1bFiTd2/qy78/Xs1TM9YyM20bD190PJ0b1ww6mkiFEsk1gvfM7EYzO7aQ6xGIFFvlxHjuPqsDL17Xkz0Hcjj38Zk89vm3ZOUcCjqaSIURydDQ2nw2u7sfdT2CaNHQUOzK3JvFH95ZxnuLN9GmfjX+dHZH+rb5xYJ2IpKPInUfLatUCOST5d9z7wcr2bB9L4M7NODuMzvQtG7VoGOJlGnFKgRmlgjcAPQPb5oCjHP3QNpHqhAIwIHsHJ6asZb/fJFGdo5zXb8WjBrUmuRKkVz2Eok9xS0EE4BE4LnwpsuBHHe/rkRTRkiFQHL7Yed+/vXhKt5c+B0NalTiriHHMaxbI613IJJHcQvBYnfvWtC20qJCIPlZsGEH97y7nMUZmZzQtBZ/OrsjXZvUCjqWSJlR3F5DOWbWKtfBWpLrfgKRsuCEprV568aT+Pf5XdiwfR9Dx8zkztcWs2XX/qCjiZR5kQyo3gl8aWbphJaqbIbaUEsZFBdnXJDahDM6HcN/vkjj6Zlr+XDZ99xySmuu6tOCpAR1VBHJT0SzhsysEtCOUCFY5e4Hoh3sSDQ0JJFK37qbv32wks9XbaFFvWTuPqs9Jx/XIOhYIoEoiTbU3YFOQFfgIjO7oqTCiURLy5RqPHXViTx79YmYwTXPzuOqZ75mzdbdQUcTKVMiuVg8EWgFLOLnawPu7rdEOVu+dEYgRXEw+xDPz17HI599y76sHK7q05xbTm1DjcqJQUcTKRXFnTW0EujgZeTOMxUCKY5tuw9w/8ereWXeRuomJ3Hn6e24oHsT4rQqmlRwxR0aWgYcU7KRRIJRr1ol/nleF94d1ZfmdZP53RtLGTpmJvPWbQ86mkhgIikE9YAVZvaxmb17+CfawUSiqXPjmrw2sjePDO/G1l0HOH/sbEZPWsjmzH1BRxMpdZFMH/1ztEOIBMHMGNqtEYM7NOCJKWsYNy2dT5b/wKhBrbiuX0sqJ2rtA4kNajonErZx+17+9sFKPlr+PY1rV+EPZ7bn9I7HqF2FVAjFukZgZr82s2/NLNPMdprZLjPbWfIxRYLVpE5Vxl7enZeu60lyUgIjX1jApRO+YtX3+usuFVsks4bSgLPdfWXpRDo6nRFIacjOOcRLX2/ggU++Ydf+LC7r1YzbB7elVtWkoKOJFElxZw39UFaKgEhpSYiP44rezZlyx0Au69WMF+asZ+D9U5g4ex3ZWh1NKphICsE8M3vFzC4ODxP92sx+HcnBzewMM1ttZmlmdlc+z19qZkvCP7PMLJCOpiJHUjs5ib8M7cTk0f1of0wN7n5nOWc9NoNZa7YFHU2kxERSCGoAe4HTgLPDP2cV9CIziwfGAEOADsDFZtYhz25rgQHu3gX4KzA+8ugipee4Y2rw0m968sSlJ7BrfzaXPPkVN7wwn43b9wYdTaTYCpw+6u5F7TTaA0hz93QAM5sEDAVW5Dr2rFz7zwEaF/G9RKLOzBjS+VgGHVefJ6el8/iUNXyxagvX92/JyIGtqJqk1dGkfCrwb66ZVQauBToClQ9vd/drCnhpI2BjrscZQM+j7H8t8GFBeUSCVjkxnptPacN53Rvzzw9X8egXabw2P4O7hhzHOV0barqplDuRDA1NJNRi4nRgKqFv7bsieF1+/xrynaJkZoMIFYLfHeH5EWY2z8zmbd26NYK3Fom+hrWq8OjFx/PayN7USU5i9KRFnPfELOav3xF0NJFCiaQQtHb3u4E97v4ccCbQOYLXZQBNcj1uDGzKu5OZdQEmAEPd/cf8DuTu49091d1TU1JSInhrkdJzYvM6vHtTX/51Xmc27tjHeU/MYtRLC3T9QMqNSApBVvjPn8ysE1ATaB7B6+YCbcyshZklAcOB/9OjyMyaAm8Cl7v7NxGnFilj4uOMi05sypQ7BnLLKW34fOUPnPLAVP4xeSWZ+7IKPoBIgCIpBOPNrDZwN6EP8hXAfQW9yN2zgZuAj4GVwKvuvtzMRprZyPBufwTqAo+b2SIz051iUq4lV0rg9sFtmXLHIM7p1pDx09MZ+O8veW7WOrJ0/4GUUeo1JBJFy77L5G8frGR2+o+0TEnm90Pac2r7+rqgLKWuuAvTVALOIzQc9N9ZRu7+lxLMGDEVAilv3J3PV27h7x+uJH3rHnq1rMMfzuxAp0Y1g44mMaS4LSbeITT/PxvYk+tHRCJgZpzaoQEf39qfvwztyDc/7Obs/8zgt68u5vvM/UHHE4nojGCZu3cqpTwF0hmBlHc792cx5ss0npmxjrg4GNGvJdcPaEVyJd2QJtFT3DOCWWYWyXRREYlAjcqJ/H5Iez7/7QBObd+AR79IY+D9U5j09QZyDpWva3ZSMRzxjMDMlhK6ASwBaAOkAwcI3Sjm4f5ApU5nBFLRLNiwg3vfX8GCDT9x3DHV+X9ntqdfG90vIyWrSBeLzazZ0Q7q7utLIFuhqRBIReTuTF76Pf/8aCUbt+9jYLsU/vdX7WnboHrQ0aSCKGohqApkuXtW+HE74FfAend/M1phC6JCIBXZgewcnpu1jse+SGPPgWyG92jKbae2JaV6paCjSTlX1GsEHxG+g9jMWgOzgZbAKDP7R0mHFBGolBDPiP6tmHrnIK7o3ZxX525k0P1TGPNlGvuzcoKOJxXUUa8RuHvn8O9/Beq4+6hwu4j5h58rbTojkFiSvnU3//hwFZ+u+IGGNStz5xntGNq1EXFxuiFNCqeoZwS5K8TJwKcA7n4Q0L3yIqWgZUo1nrwilZd/04s61ZK47ZXFDHt8Jl+v3R50NKlAjlYIlpjZ/WZ2G9Aa+ATAzGqVSjIR+a/erery7qi+PHBBV7bsPMCF42Zz/cR5rN2mezul+I5WCH4DbCN0neA0dz/cU7cDcH+Uc4lIHnFxxnndG/PlHQP57eC2TP92G6c9NJW/vLeCn/YeDDqelGOFajpnZie4+4Io5imQrhGIhGzZtZ+HPv2GV+ZupHrlRG4+uTVX9G5OUkIk94lKrCnuncW5TSiBPCJSAupXr8w/ft2FyaP70aVxTe79YCWDH5rKh0s3U966CkuwClsINFVBpIw57pgaTLy2J89efSKVEuK44cUFXDhuNos2/hR0NCknClsI7olKChEptoHt6jP5ln78/dzOrN22h2FjZjJ60kIydmjJTDm6iK4RmFkjoBn/dz2CaVHMdUS6RiBSsN0Hshk7ZQ1PTk/HgSt6NWPkwFbUq6Y7lGNVcRem+RdwEaElKg/f2ujufk6JpoyQCoFI5Db9tI8HPvmGtxZmUCkhniv7NOf6/i2pnZwUdDQpZcUtBKuBLu5+IBrhCkuFQKTw0rfu5pHPv+XdxZuomhjPNX1bcF3fltSsmhh0NCklxZ01lA7ob4tIOdYypRqPDD+ej2/tz4B2KTz2RRp97/uCRz77ll37s4KOJwGL5IzgDaAr8Dmh9QgAcPdbohstfzojECm+FZt28tBn3/Dpih+oVTWREf1bcmXv5lolrQIr7tDQlfltd/fnSiBboakQiJScpRmZPPjpar5cvZW6yUmMHNCKy3o1o0pSfNDRpIQVqxCED1AFaOruq0s6XGGpEIiUvAUbdvDQp98w/dttpFSvxI0DW3Fxj6ZUTlRBqCiKdY3AzM4GFhFanwAz62Zm75ZsRBEJ0glNazPx2p68en1vWqUkc897Kxj47ylMnLOeg9lqNlzRRXKx+M9AD+AnAHdfBLSIYiYRCUiPFnWYNKI3L13Xk0a1q3D328sYdP8UJn29gawcFYSKKpJCkO3umXm2qZGJSAXWp3U9Xh/Zm+eu6UG96pW4682lnPLAVN6Yn0G2CkKFE0khWGZmlwDxZtbGzB4DZkU5l4gEzMwY0DaFt2/sw1NXplK9cgK/fW0xpz08jXcWfUfOIX0frCgiKQQ3Ax0JTR19GdgJ3BrNUCJSdpgZp7RvwPs392XsZd1Jio9j9KRFDHlkGpOXbuaQCkK5V9j1COKBZHffGb1IR6dZQyLBOnTImbxsMw9/9i1pW3bT/tga3HZqGwZ3aICZGhSXVcWdNfSSmdUws2RgObDazO4s6ZAiUj7ExRlndWnIx7f25+GLurE/K4cRE+dzzn9m8uWqLVoLoRyKZGioQ/gMYBgwGWgKXB7VVCJS5sXHGcOOb8Snt/Xn3+d34ad9B7n62bn8+olZzPh2mwpCORJJIUg0s0RCheAdd89Cs4ZEJCwhPo4LUpvw+e0D+fu5nfkhcz+XPfUVF42fw5z0H4OOJxGIpBCMA9YBycA0M2tG6IJxgczsDDNbbWZpZnZXPs8fZ2azzeyAmd1RmOAiUrYkJcRxSc+mfHnnQP4ytCPrtu1h+Pg5XDphDvPXbw86nhxFoS4W//dFZgnunl3APvHAN8BgIAOYC1zs7ity7VOf0II3w4Ad7n5/Qe+ti8Ui5cP+rBxe/GoDT0xJY9vugwxom8Jtg9vSrUmtoKPFpKNdLC6w1aCZVQLOA5rn2f8vBby0B5Dm7unh40wChhJa4AYAd98CbDGzMwvKISLlS+XEeK7t24KLezRh4uz1jJ26hmFjZnJq+/rcempbOjWqGXRECYtkaOgdQh/g2cCeXD8FaQRszPU4I7xNRGJI1aQErh/Qium/O5k7T2/H3HU7OOuxGYycOJ9l3+VtWiBBiKT5eGN3P6MIx85vQnGRLjKb2QhgBEDTpk2LcggRCVi1SgmMGtSay3s34+kZa3lq+lo+Wv49J7Wuy4j+rejfpp7uQwhIJGcEs8yscxGOnQE0yfW4MbCpCMfB3ce7e6q7p6akpBTlECJSRtSonMitp7Zl5u9P5vdDjiNty26ufPprhjwynbcWZqi5XQAiKQR9gfnh2T9LzGypmS2J4HVzgTZm1sLMkoDhgNpXiwgQKgjXD2jF9P85mX+f34VD7tz2ymIG3PclE6ans/vAUeejSAmKZIWyZvltd5S7T90AAA0tSURBVPf1BR7c7FfAw0A88LS7/83MRoZfP9bMjgHmATWAQ8Bufr6BLV+aNSRSMbk7U1ZvZdy0NcxJ3071yglc2rMZV5/UnAY1Kgcdr9wriRXK+gJt3P0ZM0sBqrn72hLOGREVApGKb0nGT4ybls6HSzeH7mDu1ogR/VvSpkH1oKOVW8Vds/hPQCrQzt3bmllD4DV3P6nkoxZMhUAkdmz4cS8TZqTz6ryN7M86xMnH1WdE/5b0bFFHF5YLqbiFYBFwPLDA3Y8Pb1vi7l1KPGkEVAhEYs/2PQeZOHs9z81ex/Y9B+napBbX92/J6R2PIT5OBSESxeo+Chz0ULXw8MGSSzKciEhB6iQnMfrUNsy662TuHdaJn/Ye5MYXF3DyA1OYOHsd+w7mBB2xXIvkjOAOoA2hVhH/AK4FXnL3R6Mf75d0RiAiOYecT5Z/z7hp6Sza+BN1kpO4vFczrujdjLrVKgUdr0wqiYvFg4HTwg8/dvfPSjBfoagQiMhh7s7cdTsYP20Nn63cQuXEOM7v3pjr+rakeT0NXuRWpEJgZrv4+U7gvINw+4E1wP9z989LKmgkVAhEJD9pW3bx5LS1vLXwO7IOHWJIp2MY0b+VmtyFFfuMIJ8DxgOdgBfdvVMx8xWKCoGIHM2Wnft5ZtY6Xpiznl37s+nRog7X92/JoHb1iYvhC8slXghyHfh6dx9X5AMUgQqBiERi94FsJn29gadnrGVT5n5a16/GiH4tGXp8QyolxAcdr9RFrRAEQYVARAojK+cQHyzZzLhp6azcvJP61Stx9UktuKRnU2pWSQw6XqlRIRCRmOfuzEjbxrip6cxI20ZyUjwX92jKNX1b0LBWlaDjRZ0KgYhILsu+y+TJ6em8v2QzBpzdtSEj+rek/bE1go4WNSoEIiL5yNixl6dnrGPS3A3sPZhDvzb1uL5/K05qXbfCtbBQIRAROYrMvVm88NV6npm5jm27D9C7ZV3+OqwTretXCzpaiVEhEBGJwP6sHF6Zu5EHPlnNvqwcRvRvyU2D2lAlqfzPMipuryERkZhQOTGeK/s05/PfDuTsrg0Z8+UaBj80lc9X/hB0tKhSIRARySOleiUevLAbk0b0okpiPNc+N4/fPD+PjB17g44WFSoEIiJH0KtlXT64pR93DTmOGd9uY/CD03hiyhoOZlesdZVVCEREjiIpIY6RA1rx6e396demHv/6aBVnPjqdOek/Bh2txKgQiIhEoHHtqoy/IpWnrkxlX1YOw8fP4fZXFrF114GgoxWbCoGISCGc0r4Bn942gJsGtea9JZs45YEpTJyznpxD5WsGZm4qBCIihVQlKZ47Tm/Hh6P706lRTe5+exm/fnwmSzMyg45WJCoEIiJF1Lp+NV68riePDO/Gdz/t55wxM/jjO8vI3JcVdLRCUSEQESkGM2Not0Z8cccAruzdnBfmrOeUB6by9sLvKC837KoQiIiUgBqVE/nzOR15Z1RfGtWqzK2vLOKSJ78ibcuuoKMVSIVARKQEdW5ckzdvPIl7h3Vi+aZMhjwynfs+WsW+gzlBRzsiFQIRkRIWH2dc1qsZX9wxkHO6NuLxKWs49cGpfLaibLaqUCEQEYmSetUq8cCFXXllRC+qJsVz3fPzuO65steqQoVARCTKerasy+TRoVYVM9O2ceqDU3l8SlqZaVWhQiAiUgoS40OtKj777QAGtE3hvo9W86tHpzN7TfCtKlQIRERKUaNaVRh3eahVxf6sHC5+cg63BdyqQoVARCQAuVtVvL9kEyc/MIWJs9cF0qpChUBEJCCHW1V8dGt/Ojeqyd3vLOfcx2eyJOOnUs0R1UJgZmeY2WozSzOzu/J53szs0fDzS8zshGjmEREpi1ql/NyqYnPmfoaOmcndb5deq4qoFQIziwfGAEOADsDFZtYhz25DgDbhnxHAE9HKIyJSlh1uVfH5b0OtKl78aj2nPDCFtxZmRL1VRTTPCHoAae6e7u4HgUnA0Dz7DAWe95A5QC0zOzaKmUREyrTDrSrevakvjWpX5bZXFnPxk3Oi2qoimoWgEbAx1+OM8LbC7iMiEnM6NarJmzf04W/ndmLFpp0MeWQ6E6anR+W9EqJy1BDLZ1ve85tI9sHMRhAaOqJp06bFTyYiUg7ExxmX9mzG6R2P4R+TV9GsbnJU3ieahSADaJLrcWNgUxH2wd3HA+MBUlNTy0dfVxGREnK4VUW0RHNoaC7QxsxamFkSMBx4N88+7wJXhGcP9QIy3X1zFDOJiEgeUTsjcPdsM7sJ+BiIB5529+VmNjL8/FhgMvArIA3YC1wdrTwiIpK/aA4N4e6TCX3Y5942NtfvDoyKZgYRETk63VksIhLjVAhERGKcCoGISIxTIRARiXEqBCIiMc6i3cyopJnZVmB9EV9eD9hWgnFKSlnNBWU3m3IVjnIVTkXM1czdU/J7otwVguIws3nunhp0jrzKai4ou9mUq3CUq3BiLZeGhkREYpwKgYhIjIu1QjA+6ABHUFZzQdnNplyFo1yFE1O5YuoagYiI/FKsnRGIiEgeMVMIzOwMM1ttZmlmdlfQeQDM7Gkz22Jmy4LOkpuZNTGzL81spZktN7PRQWcCMLPKZva1mS0O57on6Ey5mVm8mS00s/eDznKYma0zs6VmtsjM5gWd5zAzq2Vmr5vZqvDfs95lIFO78P+nwz87zezWoHMBmNlt4b/zy8zsZTOrXKLHj4WhITOLB74BBhNaDGcucLG7rwg4V39gN6F1mzsFmSW38LrRx7r7AjOrDswHhpWB/18GJLv7bjNLBGYAo8PrXQfOzG4HUoEa7n5W0HkgVAiAVHcvU3Pizew5YLq7TwivV1LV3X8KOtdh4c+M74Ce7l7U+5ZKKksjQn/XO7j7PjN7FZjs7s+W1HvEyhlBDyDN3dPd/SAwCRgacCbcfRqwPegcebn7ZndfEP59F7CSMrCWtIfsDj9MDP+UiW8yZtYYOBOYEHSWss7MagD9gacA3P1gWSoCYacAa4IuArkkAFXMLAGoSj4rORZHrBSCRsDGXI8zKAMfbOWBmTUHjge+CjZJSHj4ZRGwBfjU3ctELuBh4H+AQ0EHycOBT8xsfnjt77KgJbAVeCY8lDbBzKKzGG/RDQdeDjoEgLt/B9wPbAA2E1rJ8ZOSfI9YKQSWz7Yy8U2yLDOzasAbwK3uvjPoPADunuPu3Qitb93DzAIfUjOzs4At7j4/6Cz5OMndTwCGAKPCw5FBSwBOAJ5w9+OBPUCZuG4HEB6qOgd4LegsAGZWm9AIRgugIZBsZpeV5HvESiHIAJrketyYEj61qmjCY/BvAC+6+5tB58krPJQwBTgj4CgAJwHnhMfjJwEnm9kLwUYKcfdN4T+3AG8RGiYNWgaQkets7nVChaGsGAIscPcfgg4Sdiqw1t23unsW8CbQpyTfIFYKwVygjZm1CFf74cC7AWcqs8IXZZ8CVrr7g0HnOczMUsysVvj3KoT+gawKNhW4++/dvbG7Nyf0d+sLdy/Rb2xFYWbJ4Yv9hIdeTgMCn6Hm7t8DG82sXXjTKUCgExHyuJgyMiwUtgHoZWZVw/82TyF03a7ERHXN4rLC3bPN7CbgYyAeeNrdlwccCzN7GRgI1DOzDOBP7v5UsKmA0Dfcy4Gl4fF4gP8Nr0EdpGOB58IzOuKAV929zEzVLIMaAG+FPjtIAF5y94+CjfRfNwMvhr+YpQNXB5wHADOrSmh24fVBZznM3b8ys9eBBUA2sJASvsM4JqaPiojIkcXK0JCIiByBCoGISIxTIRARiXEqBCIiMU6FQEQkxqkQSEwzs4dyd5g0s4/NbEKuxw+Y2f+Gp+8V5rhXmdl/SjKrSLSoEEism0X4Lk0ziwPqAR1zPd8H+Nzdzw8gm0ipUCGQWDeTn2/X70jozttdZlbbzCoB7YEdh9eMCH/Tf9PMPjKzb83svsMHMrOrzewbM5tK6Ka8w9ubmdnnZrYk/GfTcPO8dAupZWaHDvcBMrPpZta6lP77RVQIJLaFe/Fkm1lTQgVhNqFOq70JrS2wBDiY52XdgIuAzsBFFlrI51jgHkIFYDDQIdf+/yG05kQX4EXgUXfPIbRGRgegL6E1H/qFi09jd0+Lxn+vSH5UCER+Pis4XAhm53o8K5/9P3f3THffT6hHTjOgJzAl3BjsIPBKrv17Ay+Ff59I6IMfYDqhvvz9gX+Et59IqDeWSKlRIRD5+TpBZ0JDQ3MIfXj3IVQk8jqQ6/ccfu7ZFWm/lsP7TQf6EeoIOhmoRaj31LTIo4sUnwqBSOjD/ixge3i9g+2EPpR7Ezo7iMRXwEAzqxtu4X1BrudmEepKCnApoWUHD7+mD3AofHaxiFCzs+nF+Y8RKSwVAhFYSmi20Jw82zIjXevX3TcDfyZUOD4j1CnysFuAq81sCaGurqPDrzlAaOW8w+87Hagefm+RUqPuoyIiMU5nBCIiMU6FQEQkxqkQiIjEOBUCEZEYp0IgIhLjVAhERGKcCoGISIxTIRARiXH/HxOgIdcT+NYFAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "ces_fig, ces_ax = plt.subplots()\n", "plt.plot(ces_conv)\n", "ces_ax.set_xlabel('Window')\n", "ces_ax.set_ylabel('Jensen-Shannon divergence')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Comparing different clustering methods" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You may want to try different clustering methods, or use different parameters within the methods. `encore.ces_convergence` allows you to pass a list of `clustering_methods` to be applied, much like [normal clustering ensemble similarity methods](clustering_ensemble_similarity.ipynb#Calculating-clustering-similarity-with-multiple-methods).\n", "\n", "
\n", " \n", "**Note**\n", "\n", "To use the other ENCORE methods available, you need to install [scikit-learn](https://scikit-learn.org/stable/).\n", "\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The KMeans clustering algorithm separates samples into $n$ groups of equal variance, with centroids that minimise the inertia. You must choose how many clusters to partition. [(See the scikit-learn user guide for more information.)](https://scikit-learn.org/stable/modules/clustering.html#k-means)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "km1 = clm.KMeans(12, # no. clusters\n", " init = 'k-means++', # default\n", " algorithm=\"auto\") # default\n", "\n", "km2 = clm.KMeans(6, # no. clusters\n", " init = 'k-means++', # default\n", " algorithm=\"auto\") # default\n", "\n", "km3 = clm.KMeans(3, # no. clusters\n", " init = 'k-means++', # default\n", " algorithm=\"auto\") # default" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "When we pass a list of clustering methods to `encore.ces_convergence`, the similarity values get saved in `ces_conv2` in order." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(9, 3)" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ces_conv2 = encore.ces_convergence(u, # universe\n", " 10, # window size\n", " selection='name CA',\n", " clustering_method=[km1, km2, km3]\n", " )\n", "ces_conv2.shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As you can see, the number of clusters partitioned by KMeans has an effect on the resulting rate of convergence." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEJCAYAAACZjSCSAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzdd1hURxfA4d/QBQXsiqBgwQ6oiL3F3jAmpmpi7yaaxCSmfJpmqoklauzG9GasscTeIzZQ7GLFiiIgKH2+PxYrKAu7S5HzPg9P3Lv3zj1rEs7OnZkzSmuNEEKIgssqtwMQQgiRuyQRCCFEASeJQAghCjhJBEIIUcBJIhBCiAJOEoEQQhRwFk0ESqkOSqmjSqkTSqkxGbzfUikVrZQKTvsZa8l4hBBCpGdjqYaVUtbANKAtEA7sUkot1VofeuDULVrrLpaKQwghxKNZLBEAAcAJrfVJAKXUb0A34MFEkCUlSpTQnp6epkcnhBAFyJ49e65qrUtm9J4lE0E54Nw9r8OBBhmc10gpFQJcAEZrrQ8+qlFPT092795tviiFEKIAUEqdedh7lkwEKoNjD9az2AtU0FrHKqU6AYuBKukaUmoQMAigfPny5o5TCCEKNEsOFocDHve8dsfwrf8OrXWM1jo27c8rAFulVIkHG9Jaz9Ja+2ut/UuWzLBnI4QQIpssmQh2AVWUUl5KKTvgeWDpvScopcoopVTanwPS4rlmwZiEEEI8wGKPhrTWyUqpEcBqwBqYp7U+qJQakvb+DKAHMFQplQzcAp7XUg5ViMdSUlIS4eHhxMfH53YojzUHBwfc3d2xtbU1+hqV337v+vv7axksFiL/OXXqFEWKFKF48eKkPQgQZqa15tq1a9y4cQMvL6/73lNK7dFa+2d0nawsFkLkiPj4eEkCFqaUonjx4lnudUkiEELkGEkClpedv+MCkwgi4xL5cNlBbiWm5HYoQgiRpxSYRLDtxFW+336aZ2fu4FK0DFYJURD169ePUqVKUatWrfuOv/nmm1SrVg0fHx+6d+9OVFSU0W326dOHv/76K8uxnD59ml9++SXL11lCgUkEXX3dmPOyPycjYgmcupXgc8b/ixZCPB769OnDqlWr0h1v27YtoaGh7N+/H29vbz777DOLx5KdRJCSYpknGgUmEQC0rl6av4c1wd7Wimdn7mBJ8PncDkkIkYOaN29OsWLF0h1v164dNjaG2fQNGzYkPDw8w+u//PJLateuja+vL2PGpCuojKenJ1evXgVg9+7dtGzZEoBNmzbh5+eHn58fderU4caNG4wZM4YtW7bg5+fHxIkTSUlJ4c0336R+/fr4+Pgwc+ZMADZu3EirVq148cUXqV27NnFxcXTu3BlfX19q1arF77//bvLfiyVLTOQ9WlO1TBGWDG/KkJ/2MPK3YI5dvsEbbatiZSWDWELklA+XHeTQhRiztlnDzZlxXWua3M68efN47rnn0h1fuXIlixcvZufOnTg6OhIZGWl0mxMmTGDatGk0adKE2NhYHBwc+Pzzz5kwYQLLly8HYNasWbi4uLBr1y4SEhJo0qQJ7dq1AyAoKIjQ0FC8vLxYuHAhbm5u/PPPPwBER0eb/JkLTo/gzA6Y1RJiIyjmZMdP/RvwfH0Ppm0IY/BPe4hLSM7tCIUQuWz8+PHY2NjQs2fPdO+tXbuWvn374ujoCJBhz+JhmjRpwuuvv86UKVOIioq60/u417///ssPP/yAn58fDRo04Nq1axw/fhyAgICAO+sCateuzdq1a3n77bfZsmULLi4u2fmo9yk4PQI7R7hyGBYNgp4LsbOx4rOnalO1TBE+Xn6Ip7/bzpze/rgXdcztSIV47Jnjm7u5LViwgOXLl7Nu3boMp2BqrTOdmmljY0NqairAfXP5x4wZQ+fOnVmxYgUNGzZk7dq1Gbb/7bff0r59+/uOb9y4EScnpzuvvb292bNnDytWrOCdd96hXbt2jB1r2p5eBadHUNYXOnwGYeth20TAMN+2bxMv5vcN4HzULbpN3cbu08Z394QQj4dVq1bxxRdfsHTp0jvf+B/Url075s2bx82bNwEyfDTk6enJnj17AFi4cOGd42FhYdSuXZu3334bf39/jhw5QpEiRbhx48adc9q3b893331HUlISAMeOHSMuLi7dPS5cuICjoyO9evVi9OjR7N27N/sfPE3BSQQA/v2g5lOw/hM4s/3O4RbeJVk0rAlFHGx4YfZ//Ln73CMaEULkVy+88AKNGjXi6NGjuLu7M3fuXABGjBjBjRs3aNu2LX5+fgwZMiTdtR06dCAwMBB/f3/8/PyYMGFCunPGjRvHyJEjadasGdbW1neOT5o0iVq1auHr60uhQoXo2LEjPj4+2NjY4Ovry8SJExkwYAA1atSgbt261KpVi8GDB5OcnP6R9YEDBwgICMDPz4/x48fz/vvvm/z3UvBqDcXHwKwWkHQLhmwFp7tVr6NuJjL8l71sO3GNgc28GNOxOtYyiCyEWRw+fJjq1avndhgFQkZ/11Jr6F4OzvDM93AzEv4eBGnP8wBcHe34vm8AvRtVYPaWUwxYsIuY+KTci1UIIXJAwUsEcM94wbo74wW32Vpb8WG3WnzyZC22HL/KU9O3c+Za+ud0QgjxuCiYiQDSxgu6w/rx940X3NarYQV+6B9AxI0Euk3bxvawq7kQpBBCWF7BTQRKQdcp4Foe/uoPcel/0TeuVIIlw5tQorA9L88N4uedD937WQgh8q2CmwjAMF7w7AK4eRUWDb5vvOA2zxJO/D2sMU2rlOC9RaGMWxJKckr684QQIr8q2IkA7o4XnFgL2yZleIqzgy1ze9dnYDMvFuw4Q5/5u4i+KYPIQojHgyQCAP/+aeMFnxhKUWTA2krxXucafNnDh52nrvHk9G2ERcTmcKBCCFNERUXRo0cPqlWrRvXq1dmxI+P/3zNSuHDhbN1z8eLFHDp0KFvX5hRJBPDAeEE/iLv20FOf9ffgl4ENibmVxJPTtrH5WEQOBiqEMMXIkSPp0KEDR44cISQkJEfWNWQnEWS0kMySJBHcZsR4wW31PYuxeHgTyrkWos/8IOZvO0V+W5gnREETExPD5s2b6d+/PwB2dna4urqmO+/y5ct0794dX19ffH192b79/lmFGzdupEuXLndejxgxgu+//x4w1BSqUaMGPj4+jB49mu3bt7N06VLefPNN/Pz8CAsLIywsjA4dOlCvXj2aNWvGkSNHAMNeCa+//jqtWrXi7bffzrB0taUUnKJzxijrC+0/hRWjYftkaPraQ0/1KObIwqGNGfV7MB8uO8Sxyzf4MLAWdjaSW4XI1MoxcOmAedssUxs6fv7Qt0+ePEnJkiXp27cvISEh1KtXj8mTJ99X0A3g1VdfpUWLFixatIiUlBRiY417BBwZGcmiRYs4cuQISimioqJwdXUlMDCQLl260KNHDwBat27NjBkzqFKlCjt37mTYsGGsX78eMNQXWrt2LdbW1nTt2jVd6WpLkd9aD6o/AGo8Ces+fuh4wW1O9jbM7FWPYS0r8WvQOXrN3UlkXGIOBSqEyIrk5GT27t3L0KFD2bdvH05OTnz+efrEsX79eoYOHQqAtbW10WWenZ2dcXBwYMCAAfz9998ZFq+LjY1l+/btPPPMM/j5+TF48GAuXrx45/1nnnnmTo0iY0pXm4v0CB6kFAROgYshhvGCIVvBqfhDT7eyUrzVoRrepYvw1sL9dJu2lTkv16dqmSI5GLQQ+cwjvrlbiru7O+7u7jRo0ACAHj16ZJgIMnNvqWm4W27axsaGoKAg1q1bx2+//cbUqVPvfNO/LTU1FVdXV4KDgzNs+97eSUalq6tVq5bleI0hPYKMOLik1SO6CouHPHK84LYn65Tj90ENiU9K5anp21h3+LLl4xRCGK1MmTJ4eHhw9OhRANatW0eNGjXSnde6dWu+++47wLBHcEzM/TupVahQgUOHDpGQkEB0dDTr1q0DDN/2o6Oj6dSpE5MmTbrzy/7ectPOzs54eXnx559/AoY9CEJCQjKMN6PS1ZYiieBh3PwM4wXH/4XtU4y6pE75oiwd0QSvkk4M+GE3MzeFySCyEHnIt99+S8+ePfHx8SE4OJh333033TmTJ09mw4YN1K5dm3r16nHw4MH73vfw8ODZZ5/Fx8eHnj17UqdOHQBu3LhBly5d8PHxoUWLFkycaKhj9vzzz/PVV19Rp04dwsLC+Pnnn5k7dy6+vr7UrFmTJUuWZBhrRqWrLaXglaHOCq3hz95weDn0XQHlGxp12a3EFEb/GcI/By7ydF13Pn2qFvY21plfKMRjTMpQ5xwpQ21OSkHgt+Dqken6gnsVsrNm6ot1GNWmCgv3hvPCrP+IuJFg4WCFECJ7JBFk5vZ4QVyE0eMFYNgGc1Qbb6a9WJdDF2PoNnUrBy9EWzZWIYTIBkkExnCrc3e8YMe3Wbq0s09Z/hrSmFQNPb7bwarQi5lfJIQQOUgSgbHqD4Aa3WDth3B2Z5YurVXOhaUjmlC1TBGG/LSXb9cdl0FkIUSeIYnAWPeNF/Q1bHWZBaWcHfhtUEOe9HPj6zXHePW3YOKTUiwUrBBCGM+iiUAp1UEpdVQpdUIpNeYR59VXSqUopXpYMh6T3TtesMj48YI7l9taM/E5P97qUJXl+y/w7MwdXI6Jt0ysQghhJIslAqWUNTAN6AjUAF5QSqVbvZF23hfAakvFYlZudaDdeDi+OsvjBWAYRB7WsjIze9XjxJVYAqduJeRclAUCFULcKz4+noCAgDvz98eNG5el66UMNaCUcsr8rPsEACe01ie11onAb0C3DM57BVgIXMli+7knYCBUD8zWeMFt7WqWYeHQxthYWfHszB2sPSQrkYWwJHt7e9avX09ISAjBwcGsWrWK//77z+L3fSzKUCulGiulDgGH0177KqWmG9F2OeDcPa/D047d23Y5oDsww+iI8wKloNvUu+sLsjhecFv1ss4sGdGEamWKMPinPSzaF27mQIUQtyml7nyrT0pKIikpCaVUuvOkDHXGJgLtgaUAWusQpVRzI65L/zcMD06VmQS8rbVOyehfyJ2GlBoEDAIoX768EbfOAbfHC+a2M4wXvPAbWGX9SVuJwvb8PLAhAxfs5rXfQ4i5lUzvxp5mD1eIvOSLoC84Emne2jnVilXj7YC3H3lOSkoK9erV48SJEwwfPvxOAbp7SRnqh9Ban3vgkDHTXcIBj3teuwMXHjjHH/hNKXUa6AFMV0o9mcH9Z2mt/bXW/iVLljQm5Jxx33jB1Gw3U9jehvl969O2RmnGLT3IFJleKoRFWFtbExwcTHh4OEFBQYSGhqY7R8pQZ+ycUqoxoJVSdsCrpD0mysQuoIpSygs4DzwPvHjvCVprr9t/Vkp9DyzXWi82Mva8IWAgnN4Maz8w1CLyCMhWMw621nzXsy5v/bWfb9YcI/pWEu91qo6V1cN7SkLkV5l9c7c0V1dXWrZsyapVq6hVq1aWri2oZaiHAMMxPN8PB/zSXj+S1joZGIFhNtBh4A+t9UGl1BCl1JDsh5zHKAWBU8HFHf7M+vqCe9lYWzHhGV/6NPZk7tZTvLVwP8kpWZuiKoTIWEREBFFRhhl6t27deugvVilDnQGt9VWtdU+tdWmtdSmtdS+ttVHV17TWK7TW3lrrSlrr8WnHZmit0w0Oa637aK3/yvpHyAMKuRrGC2Ivw+Khhqql2WRlpRjXtQYjW1fhrz3hDP9lLwnJsvBMCFNdvHiRVq1a4ePjQ/369Wnbtu19g763SRnqjE5QagEwUmsdlfa6KPC11rqfxaJ6hBwtQ51VO2fCyreg3SfQ+BWTm5u39RQfLT9Ek8rFmfWSP072sqGcyL+kDHXOsUQZap/bSQBAa30dqGNSlI+rgEFQvathvOBckMnN9WvqxYRnfPnvZCQ95+wk6qbshyyEMD9jEoFVWi8AAKVUMWSv44zdHi9wLmfyeMFtPeq5M71nXQ5diJGSFEIIizAmEXwNbFdKfayU+hjYDnxp2bDysfvGC4aZNF5wW/uaZfi+b33Cr9+ix4ztnL120/Q4hcgFMi3a8rLzd2zMYPEPGOb4X8ZQBuIprfWPWb5TQVKurmGc4NhK2DHNLE02rlyCXwY25EZ8Mj1mbOfoJcutMhTCEhwcHLh27ZokAwvSWnPt2rUsLz4zas/itMJwpbnnkZDW+mxWgzSHPD1YfC+t4fdecGwV9F0FHvXN0uyxyzd4ae5O4pNSmd+3PnXLF838IiHygKSkJMLDw+/MuxeW4eDggLu7O7a2tvcdf9RgsTGzhl4BxmHoEaRgKB2htdY+Zok6i/JNIgC4FQUzmxmSwuDN4FjMLM2ei7xJr7k7ibiRwKyX/GlapYRZ2hVCPL5MnTU0Eqiqta6ptfbRWtfOrSSQ79weL7hxyWzjBQAexRz5c3AjyhdzpN/3u2T7SyGESYxJBOcA2XU9u8rVM/t4Adzd8axmOWeG/byXP3Y/WA5KCCGMY8w00JPARqXUP0DC7YNa628sFtXjpsFgOL0F1o4z1CNyz7B3lmWujnb8PKABg3/cw1t/7SfmVhIDmlU0S9tCiILDmB7BWWANYAcUuedHGEsp6DYNnN3Mtr7gNkc7G+b09qdT7TJ88s9hvv73qMzKEEJkSaY9Aq31h2DYoUxrHWf5kB5ThVyhx/cwrz0sGQ7P/2JIEGZgb2PNty/UpYj9Ab5df4KYW0mM61pTKpcKIYxizA5ljbK5Q5l4kHs9aPcxHF0B/5n3r9DaSvH507UZ2MyLBTvO8PofwSRJ5VIhhBGMGSOYRPZ2KMtTzseeZ/elPDDttERZqNIYtn8GNklQ1AsnWydaerTExsq0yh1KKd7tVB1XRzu+Wn2U2IRkpr5YFwdbazMFL4R4HBn1m0drfe6BrSTzXV3k0KuhvL/t/dwO464SReHQ3Dsvm7s356vmX+Fom35Xo6xQSjG8VWWcC9kydkkovecFMae3P0UcbDO/WAhRIFlyh7I8pVm5Zqx8amVuh3HXpVD4szd4NmNLnaf5fNfn9FnVh6mtp1LKsZTJzb/UsALODja88UcIL87eyfd961O8sL0ZAhdCPG6MWVlcApgMtMGwqvhfDPsTGLU5jbnlq5XFmdkxHVa/A+0/Y7NHLUZvGo2LvQvTWk/Du6i3WW6x/shlhv60F/eihfixfwPcXAuZpV0hRP5iUomJvOaxSgRaw289DfWIyjfkcGlvRkTt4qZO4ZtWE2nk1sgst9l58hoDFuzGuZAtPw1ogFcJp8wvEkI8VkytNTQlg8PRwG6tdcZ7rFnQY5UIAOKjYfMEOL0VLu3nktIMK1OSU7Z2jHWqSvfK3cGjARSvZNJ009Dz0fSeF4RSsKBfADXdXMz4IYQQeZ2piWAWUA34M+3Q08BBwAM4qbUeZcZYM/XYJYJ7Jd6EC3uJPb2ZN84uZbuOY2BUNK9cj0Y5ljAkhPINDP8s6we2WSs1GxYRy0tzdnIjIZl5fepT39M8RfCEEHmfqYlgPdBOa52c9toGwzhBW+CA1rqGmeN9pMc6EdwjKTWJ8Ts+YeGJv+nkUo2PdXHszgVBZJjhBGs7cKsDHgHg0dCQHAqXzLTd81G3eGnOTi5E32JGr3q0rGr6wLQQIu8zNREcBQK01tFpr12AnVrrakqpfVrrHN2/uKAkAjBsMjE3dC6T906mbqm6THliCi5JiRAeBGf/g3M74cI+SEnby7hYRUNSuN1rKFEVrNKvGbwam0DveUEcu3yDb571o6uvWw5/MiFETjM1EfQH3gc2Ypg11Bz4FPgV+EBr/aZZo81EQUoEt608tZL3tr5HucLlmN5mOh5FPO6+mZwAF4Lh3H9wdqchOdy8anjPwTWtx5DWayhXD+wM6xRi4pPo//0udp+5zvgna/Nig/K58MmEEDkl24lAGVaRuQPJQACGRBCktb5giUCNURATAcCey3sYuWEk1sqaKU9Mwbekb8Ynag2RJ+/2GM7thIgjhvesbKBM7Tu9hltl/Bm69CIbj0bwdodqDG1ZKec+kBAiR5naI9ijta5nkciyoaAmAoDT0acZtm4YV25e4bNmn9G2QlvjLrwZCeG7DEnh7E44vweSbwGgXTzYnVqVJZEeVKzzBH27d0JZm1bqQgiR95iaCKYB32utd1kiuKwqyIkAIDI+klfXv8r+iP284f8GL9d4GZXVaaUpSXBpf9qjpP/QZ3eiYi8BEG/liL1nA5RHA/BsAp7NzFYlVQiRe0xNBIeAqsBpIA7ZszjXxSfH8+7Wd1lzZg3PVX2OMQFjTCtYpzU66gz//LOY60e20MrpFOUSTqLQ0HA4tB8vyUCIfO5RicCY3x4dzRyPMJGDjQMTWkxg0p5JzD84n4txF00rWKcUqqgnXXqNYuamrjRdeYT2lR2ZVvofbP6bBtY20OZDSQZCPKYy3Y9Aa30Gw+KxJ9L+fNOY64RlWSkrXvd/nfcbvM/W81vps6oPV25eMbndwS0q8flTtVkTdpPnz3Xnlm8f2DYZNow3PWghRJ5kzMY044C3gXfSDtkCP1kyKGG856o9x7dPfMvpmNP0XNGTY9ePmdzm8wHlmfpiXfZfiKHZgU6c83oGNn8FG78wQ8RCiLzGmG/23YFADOMDpE0dlT2L85Dm7s1Z0GEBqamp9F7Zmx0XdpjcZqfaZVn+SlPKFHWk+eFuBLl0gI2fwpavzRCxECIvMSYRJGrDiLIGw97Flg1JZEf14tX5ufPPlC1clmFrh7Ho+CKT2/QuXYRFw5rwSuuq9Ix4iVVWzWHdR7D9WzNELITIK4xJBH8opWYCrkqpgcBaYLZlwxLZUcapDD90+IGAsgGM3T6WKXunYGqZcVtrK15v681fQ5vyteMolqc0hH/fJ3HrNDNFLYTIbcYMFk8A/gIWYphGOlZrbdRXQqVUB6XUUaXUCaXUmAze76aU2q+UClZK7VZKNc3qBxD3K2xXmKmtp/J0laeZfWA2Y7aMIfF2LSIT+Hq4smxkS4Lrf8XKlPrYrX2Xs6szqlAuhMhvjFlH8Brwp9Y6PEsNK2UNHMNQpTQc2AW8oLU+dM85hYE4rbVWSvkAf2itqz2qXVlHYJx7C9bVK12Pya0m42Jvnj0Idhy7QMpvL9M0dRerK75LqxffxM5GJpIJkZc9ah2BMf/3OgOrlVJblFLDlVKljbxvAHBCa31Sa50I/AZ0u/cErXWsvpuJnEgbhxCmU0oxoPYAvmz+Jfsj9tNrRS/O3ThnlrYbebvh8/piDhduQNuwz5j6zYccuRRjlraFEDnPmEdDH2qtawLDATdgk1JqrRFtlwPu/c0TnnbsPkqp7kqpI8A/QL+MGlJKDUp7dLQ7IiLCiFuL2zp6dWR2u9lcT7hOrxW9CIkIMUu7zoULU33kUq6Xacyom5OZM/VzZmwKIyVVcrkQ+U1W+vNXgEvANcCY3UwyWoaa7reE1npR2uOgJ4GPM2pIaz1La+2vtfYvWTLzzVfE/eqVrsdPHX/CydaJ/qv7s+bMGvM0bOtA8f5/keLRhC9tviN09Xyem7mDM9fizNO+ECJHGLOgbKhSaiOwDigBDDSyzlA4hhXJt7kDDy1frbXeDFRSSpUwom2RRZ4unvzU6SeqFavGGxvfYMHBBSbPKALAzhHbl/5AlW/AFPvpuF9eS8fJW/jpvzPmaV8IYXHG9AgqAKO01jW11uPuHezNxC6gilLKSyllBzwPLL33BKVU5bQ9D1BK1QXsMPQ4hAUUcyjGnHZzaFOhDRN2T2D8zvEkpyab3rCdE6rnn1i5+zPRajIDSh7h/cWh9Jm/i0vR8aa3L4SwqIcmAqWUc9ofvwTOKqWK3fuTWcNpexyPAFYDhzHMCDqolBqilBqSdtrTQKhSKhiYBjyn5WukRd0uWNe3Zl9+P/o7IzeM5GbSTdMbti8CPf9ElfXlteufMLdxJDtPXaP9pM0sCT4vvQMh8rCHTh9VSi3XWndRSp3C8Gz/3mf+WmtdMScCfJBMHzWf34/8zqdBn1K1aFWmtp5KKUczbGR/Kwp+CIQrR7jUeT5D/3Nl39koOvuU5ZNutSjqZGf6PYQQWWbSfgR5jSQC89ocvpnRm0bjYu/CtNbT8C7qbXqjNyNhQSBcO07K878z45w7k9Yew9XRji+ers0T1YydgSyEMJdsrSNQStV91I/lwhU56XbBupTUFLMVrMOxGLy8BIpVxPr3FxjudZklw5tS3MmOft/vZszC/cQmmGFsQghhFo96NLQh7Y8OgD8QguHxkA+wU2udK+UgpEdgGZfiLjFs3TBORZ1ibKOxdK/S3fRGY6/A950h+jy8tIgEN38mrjnOrM1huLkWYsIzvjSsWNz0+wghMpWtHoHWupXWuhVwBqibNo+/HlAHOGGZUEVuuV2wrn6Z+mYrWEfhUtB7GRQpAz89jf2lfYzpWI0/BjfC2krxwuz/+GT5IeKTUszzIYQQ2WLM9NFqWusDt19orUMBP8uFJHJLYbvCTGszzbwF64qUMSQDp+Lw41NwYR/+nsVY8WozejYoz5ytp+j67VYOhEeb50MIIbLMmERwWCk1RynVUinVQik1G8N0UPEYsrWyZVyjcYysO5IVp1YwasMoElISTGvUpRz0Xg6FXOCHJ+Hifpzsbfjkydos6BdATHwS3advY/La4ySlpJrngwghjGZMIugLHARGAqOAQ2nHxGPqdsG6sY3GsuX8Fl5Z9wq3km+Z1qirh6FnYFcYfugGlw3rElt4l+TfUS3o4lOWiWuP8fR32zlx5YYZPoUQwlgyfVQ80qLjixi3fRz+ZfyZ+sRUHG0dTWvwWphhADk1Gfr8AyWr3nnrn/0XeX/xAW4mpvBWh2r0beyJlVVGJauEEFllahlqUYB1r9KdT5t9yp7Lexi6diixibGmNVi8kuExkbKCBV3h6vE7b3X2Kcvq15rTtHIJPl5+iBfn/Me5SDOsehZCPJIkApGpLhW78GXzLwmJCGHwmsHEJJq490CJyobHRKkphmRwLezOW6WKODCntz9fPu3DgfBoOk7ewh+7z0mJCiEsSBKBMEp7z/Z83fJrDkUeYuC/A4lOMHGWT8mqhmSQnGBYhXz99J23lFI8W9+DVaOaU9PNmbf+2s/AH3Zz5YYUsCAScwUAACAASURBVBPCEowpQ+2tlJqtlPpXKbX+9k9OBCfyltblWzO51WROXD9Bv9X9iIyPNK3B0jUMK5ATYw09g6j7d1DzKObIrwMb8n7n6mw+fpX2Ezez8sBF0+4phEjHmD2LQ4AZwB7gzsofrfUey4aWMRkszn3bz2/n1Q2v4l7YnTnt51CikIlbSFzYBwu6GUpT9F0Bzm7pTjl++Qav/xHCgfPRPOnnxtiuNSkmBeyEMJpJRefSLq5nkciyQRJB3hB0MYgR60dQ2rE0c9rNobSTiYXkwncb1hgULmVIBkXKpDslKSWVaRtO8O36EzjZWTOqjTcvNaqArbU84RQiM6bOGlqmlBqmlCqblf0IxOMtoGwAM9rMIOJWBH1X9+VirImPbNz9oddCuHHJ8Jgo9kq6U2ytrRjVxpuVI5vh4+7KR8sP0WnyFrYcl32shTCFMT2CUxkclv0IBAAhESEMXTMUZ3tn5rSbg3sRd9MaPL0Nfu4BRT0N00ydMi5Kp7VmzaHLfPLPYc5G3qRN9dL8r0t1KhR3Mu3+QjymZD8CYVEHrx1k0L+DKGRTiLnt51LBuYJpDZ7cBL88C8WrQO+lhrGDh0hITmHu1lNMXX+C5BRNv6ZejHiiMoXtbUyLQYjHjEmPhpRStkqpV5VSf6X9jFBK2Zo/TJFf1Sxek3nt55GYkkjfVX05GXXStAYrtoDnf4Grx+DHJw27nj2EvY01w1pWZsPolnTxLcuMTWG0mrCRv/aEk5qav77kCJFbjBkj+A6oB0xP+6mXdkyIO6oWq8q89vNI1an0Xd2X49ePZ37Ro1RuDc/9BFcOw09PQfyj1y2Udnbgm2f9WDSsMW6uhRj9Zwjdv9vOvrPXTYtDiALAqOmjWmvfzI7lFHk0lLedij7FgNUDSExNZHa72VQrVs20Bo+sgD9egnL1DIPJ9kUyvSQ1VbNo33k+X3WEiBsJPFW3HGM6VKOUs4NpsQiRj5k6ayhFKVXpnsYqcs96AiHu5eXixfwO83GwcaD/6v4cvHrQtAardYIe8wzTS39+FhLjMr3EykrxdD13NoxuydCWlVgecpFWEzYyfeMJEpLlP10hHmRMj6A1MB84iWGrygpAX631hkdeaCHSI8gfzseep//q/kQnRPNdm+/wK2XiXkahC2HhAKjQBF78A+yMr4J6+moc41ccZs2hy1Qo7sh7narTtkZplJLKpqLgMHnWkFLKHqiKIREc0VqbuFNJ9kkiyD8uxV2i/+r+XL11lWmtp+FfJsP/Bo23/w/4exBUbAnPLgAHlyxdvuV4BB8tO8TxK7E0q1KCsV1qUKV05o+ahHgcmCMRNAY8gTtz8rTWP5grwKyQRJC/XLl5hQH/DuBS3CWmPDGFhmUbmtZg8C+weBgUKgot3gL/fmBjb/TlSSmp/PTfGSauOUZcYgovNazAa228cXGUiXDi8WZqiYkfgUpAMHfHBrTW+lWzRmkkSQT5z9VbVxn470DO3TjH5FaTaVKuiWkNXtgHa8bBqU3gWgFaj4WaT4GV8aUmIuMS+frfo/wadBaXQra80a4qLwSUx1o2whGPKVMTwWGghs4jK88kEeRP1+OvM2jNIMKiwpjYciItPFqY1qDWELYO1nwAlw9AWV9o8yFUapWlZg5diOHDZQfZeSqSamWKMK5rTRpVyng1sxD5mamzhkKB9BXAhMiCog5FmdNuDt5FvRm1cRTrzqwzrUGloHIbGLwZus+Em5GGxWc/doeL+41upoabM78Nasj0nnW5EZ/MC7P/Y9jPewi/LjujiYLDmB7BBsAPCALuDBJrrQMtG1rGpEeQv91IvMHQtUMJvRrK580+p4NXB/M0nBQPu+bAlglw6zrUfhaeeB+KGl/uIj4phVmbTzJ94wm0hsHNKzKkZSUc7aRchcj/TH00lGEfXmu9yQyxZZkkgvwvLimOYWuHERwRzCdNPqFrpa7ma/xWFGydCDtngE6F+gOh+ehH1it60IWoW3y+8ghLQy5Q1sWBdzpVp6tPWZluKvI1KTon8pybSTd5df2rBF0K4sPGH9K9Snfz3iD6PGz81DDLyK4INB0FDYeCbSGjm9h1OpIPlx0k9HwM/hWK8kFgTWqVy9qUVSHyClN7BE8BXwClMKwjUBhmDTmbO1BjSCJ4fMQnxzNqwyi2XdjG/xr+j2erPmv+m1w+BOs+hGOroIgbtHoHfF8Ea+Me96Skav7ac46vVh/lWlwiz/l7MLp9VUoUNn7KqhB5gamJ4ATQVWt92BLBZZUkgsdLQkoCb2x8g03hmxgTMIae1Xta5kant8GasXB+N5SsBm0+AO8OhkFnI8TEJ/HtuuPM33aaQrbWjGxThZcbeWJnI7ujifzB1FlDl7ObBJRSHZRSR5VSJ5RSYzJ4v6dSan/az3alVK4UshO5x97anoktJ9KmfBs+D/qc+aHzLXMjzyYwYC08+wOkJMGvz8P8TnBul1GXOzvY8l7nGqx+rTn1PIvyyT+H6TB5MxuOpt9JTYj8xpgewWQM00cXc/+sob8zuc4aOAa0BcKBXcALWutD95zTGDistb6ulOoIfKC1bvCodqVH8HhKSk3i3S3vsur0Kkb4jWCw72DL3SwlCfYugI1fQNwVqN4VWo+DElWMbmLDkSt8vPwQJ6/G8US1UrzfuToVSxa2XMxCmMjUR0MZfUXTWut+mVzXCMMv9vZpr99Ju/Czh5xfFAjVWpd7VLuSCB5fyanJjN02lmUnlzHYZzDD/YZbdqZOQizsmAbbp0DSLajXG1qMgSKljbo8MTmVBdtPM2XdceKTU+jT2JMRrapIuQqRJ+XKrCGlVA+gg9Z6QNrrl4AGWusRDzl/NFDt9vkPvDcIGARQvnz5emfOnLFIzCL3paSm8OGOD1l0YhH9avVjVN1Rlp+2GXsFNn0Je+aDtR00GgFNXjVq7wOAiBsJTFh9lD/2nKOwvQ39m3rRr6kXzg6SEETeYWqPwAHoD9QE7uzsYUSP4Bmg/QOJIEBr/UoG57bCsPtZU631tUe1Kz2Cx1+qTuXTnZ/y+9Hf6VW9F2/Vfytn5vBfC4P1H8PBReBYAlq8DfX6gI2dUZcfvhjD5LXHWXXwEs4ONgxqXpE+Tbxk/2SRJ5g6WPwjhjGC9sAmwB24YcR14YDHPa/dgQsZBOcDzAG6ZZYERMFgpax4r8F79Krei58O/8T4neNJ1amWv3HxSvDM9zBwPZSqDivfhGkBhr0QjOg5Vy/rzIyX6rH8laYEeBVnwr/HaPbFer7bGEZcQrLl4xcim4zpEezTWtdRSu3XWvukbVy/Wmv9RCbX2WAYLG4NnMcwWPyi1vrgPeeUB9YDL2uttxsTsPQICg6tNRP3TmR+6HyervI0YxuNxUrl0HRNreH4Glg7Dq4cArc60PYj8GpudBMh56KYtPYYG45GUNzJjiEtKtGrYQUK2VlbMHAhMmZqjyAp7Z9RSqlagAuGvQkeSWudDIwAVgOHgT+01geVUkOUUkPSThsLFAemK6WClVLyG17coZTitbqvMchnEAuPL+R/2/5HSmoObTWpFHi3gyFb4cnvIDYCFnSFn3rApVCjmvD1cGV+3wD+HtaYGm7OjF9xmGZfbmDe1lPEJ8mWmSLvMKZHMABYCPhg2LKyMDBWaz3D8uGlJz2CgmlGyAymBU+jo2dHxjcbj61VDg/EJt2CoFmw5WuIjwHfF6DVu+Dqkfm1aXadjmTimmNsD7tGaWd7hreqzHP1PbC3kR6CsDypNSQeC/NC5zFxj2Hx2VC/oXi5eOV8Qrh1HbZ8AztnGl43GARNX89SUbsdYdeYuOYYQacjKeviwPBWlXnW30NWKQuLMnXWkD3wNOm3qvzIjDEaTRJBwfbjoR/5cteXANha2VLZtTLVilWjarGqVCtWDe+i3hSxy4F9iKPOwYZPIeRXcHA2JIMGg40uaqe1ZnvYNb7+9yh7z0ZRzrUQrzxRmafruWNrLQlBmJ+piWAVEA3s4e5WlWitvzZnkMaSRCDOxpwl9GooRyKPcCTyCEevHyUyPvLO++6F3e9LDtWKVaO0Y2nLTEG9FAprP4ATa8DZHRoNA7+eUMjVqMu11mw+fpVv1hwj5FwU5Ys58mrrKjzp54aNJARhRqYmglCtdS2LRJYNkgjEg7TWRNyKMCSFyKN3ksOZmLsLD13sXahW9G5yqFqsqnkfLZ3abOghnN0Bto7g8xwEDILSNYz+DBuOXuGbNccIPR+DVwknXm1dmUDfcrKPsjALUxPBLOBbrfUBSwSXVZIIhLHikuI4fv343Z5D5FGORx0nIcVQMiujR0tVi1alsJ0JNYMuhhgGlQ/8Bcnx4NnMkBCqdjKq9LXWmjWHLjNx7XEOX4yhUkknRrbxpkvtslhJQhAmyFYiUEodADSGcYEqwEkMRedu70fgY5lwH00SgTBFcmoyZ2LO3Nd7OBJ5hOsJ1++c41HE405SuJ0ksvxo6WYk7P0Bds2F6LOGx0b1+0HdPuBUPNPLU1M1qw9eYuLaYxy7HIt36cKMauNNh5plJCGIbMluInjkZq9a61wp+COJQJibMY+WXO1dDb2Gex4vebp4Zv5oKTXFsCnOzplwahNY20Otpw2zjdzqZBpbaqrmnwMXmbT2GGERcVQrU4TX2nrTroaFxjzEYyu7icARSNJaJ6W9rgp0As5kVoLakiQRiJyS2aMlOys7KhetfKf3UL14dWqXqI2N1UMeAV05ArtmQ/CvkBQH7vUhYDDU6JZpPaOUVM2ykAtMXnecU1fjqOnmzOttvXmiWilJCMIo2U0Em4H+WuvjSqnKQBDwM1ADCNJav2OpgB9FEoHITZk9WirmUIxOXp0IrBRItWLVMv4lHR9tSAZBsyAyDJxKgX9fqNcXnMs++v4pqSwOvsCUdcc5G3kTX3cXXmvrTQvvkpIQxCNle4xAa1077c8fA8W01sOVUnbAntvv5TRJBCKvuf1oKSQihBUnV7ApfBNJqUlUdq1MYKVAOnl1orRTBnscpKZC2HpDQjj+L1hZQ/VAw3oEjwaP3EYzKSWVv/eGM2XdCc5H3aJueVdea+tN08olJCGIDGU3Eey/PSCslNoGfKW1Xpz2OkRrnSvbSkoiEHlddEI0q0+vZmnYUkIiQlAoGpZtSNdKXWldvjWOto7pL7oWBrvnwd4fISEayvgYZhvV7vHIRWqJyan8ueccU9ef4GJ0PPU9i/JaW28aVyphwU8o8qPsJoKfgEsYKoeOAby01jeVUq7AJkkEQmTuTMwZlp9czrKwZZyPPU8hm0K0rdCWrpW6Ur90faytHqgzlBgH+3+HoNmGqqeFikLdl6H+AHAt/9D7JCSn8Puuc0zbcILLMQk0qlic19p6E+BlfOkL8XjLbiIoBIwEygLztNYhaccbA5W01j9aKN5HkkQg8qNUncq+K/tYFraM1adXE5sUS2nH0nSu2JnASoFUcq10/wVaw+mthsdGR/4BtGEtQsBA8Grx0MdG8Ukp/LLzLNM3hnE1NoGmlUvwWltv6lUoavkPKfI0sxWdU0rV1VrvNVtk2SCJQOR38cnxbAzfyLKwZWw7v40UnUKN4jUIrBRIB88OFC/0wDqD6HDDeoS9C+DmNShR1ZAQfF8A+4wXv91KTOGn/84wY1MY1+IS6Vy7LGO71qC0s0OG54vHnzkTwV6tdV2zRZYNkgjE4+TarWusPLWSpWFLORx5GGtlTdNyTelaqSstPVpib21/9+SkeDj4t2FNwsVgsHc21DWqPwBKVM6w/biEZGZvOcn0jWHYWVvxRjtvXmpYQeoYFUDmTAT7tNaZr4KxIEkE4nF14voJlp1cxvKTy7ly8wpFbIvQzrMdgZUCqVOqzt3ZQFpD+G4ImgkHF0NqElRuYxhcrtwWrNL/kj99NY7/LQlly/Gr1CrnzPgna+PrYVxhPPF4MGciePL2zKHcIolAPO5SUlMIuhTE8pPLWXNmDbeSb+Fe2J2ulbrStWJXPJzv2QznxmXY871hxlHsJSjqZXhslEEFVK0Nq5Q/WnaIiNgEejWowOj2VXEplMN7OohcYXIiUEqVAypw/34Em80WYRZIIhAFyc2km6w7u46lYUvZeXEnGo1fST+6VupKe8/2uNi7GE5MToQjy2DnLDj33yMroMbEJ/HNv8f4YcdpijnZ878u1Qn0dZP1B485U6uPfgE8Bxzi7n4EWmsdaNYojSSJQBRUl+IuseLUCpaFLeNE1AlsrWxp6dGSrhW70rRcU2yt077ZZ1QBtd3H6WobHQiP5r3FB9gfHk2TysX5uFstKpY0ofKqyNNMTQRHAR+tdYIlgssqSQSioNNacyTyCEvDlrLi1Aoi4yMpal+UDl4dCKwUSM3iNQ3f7m9XQP3vO8NsozYfQMNh940hpKRqftl5hi9XHSUhOZUhLSsxrGUlHGxlH+XHjamJYCXwjNY61hLBZZUkAiHuSkpNYseFHSwLW8b6s+tJTE3Ey8WLwEqBdPbqTNnCZQ0JYekrcGS5YTD5ye+gcMn72rlyI57x/xxmSfAFPIs78lG3WjT3LvmQu4r8yNREsBDwBdZh2I8AAK31q+YM0liSCITIWExiDGtOr2Fp2FL2XtmLQlG/TH26VupKZ89O2O77AVa9Cw4u8NRMqPREuja2Hr/K/5aEcupqHF18yjK2Sw1KydqDx4KpiaB3Rse11gvMEFuWSSIQInPhN8LvlLY4e+MsLd1b8k3Lb7C9egz+6gcRR6DJSHjif2B9/6yh+KQUZmwKY/rGMOytrRjdviq9GlaQLTPzOXPMGioElNdaHzV3cFkliUAI42mt+fXIr3wW9BltK7Tly+ZfYpOcCKvfhT3zoVw9eHouFPNKd+2pq3GMTVt7ULucC+O718LHXdYe5FePSgSZLi9USnUFgoFVaa/9lFJLzRuiEMISlFK8WP1F3qr/FmvOrOHdLe+SYmMPXSfBMwvg2gmY0cwww+gBXiWc+KFfAFNeqMOlmHi6TdvGuCWhxMQn5cInEZZkzDrzD4AAIApAax0MpP/6IITIs16q8RKj6o5i5emVjN0+llSdCjWfhCFboXRNWNgfFg+DhPvnhCilCPR1Y90bLXi5YQV++O8Mrb/exNKQC2RlMarI24xJBMla6+gHjsl/AULkM/1r92eY3zCWhi3lox0fGZKBa3no8w80fwuCf4FZLQzrEB7g7GDLh91qsWR4E8o4O/Dqr/t4eV4Qp67G5cInEeZmTCIIVUq9CFgrpaoopb4Ftls4LiGEBQzxGcLA2gNZeHwhnwd9bvhWb20DT7wHvZdB4k2Y0wZ2TDfUNHqAj7sri4c34cPAmgSfjaL9pM1MWnuM+KSUDO4m8gtjEsErQE0MU0d/BWKAUZYMSghhGUopXqnzCr1r9ObXI7/y9e6v7z7i8WoGQ7cZCtitfgd+eQ7irqZrw9pK0buxJ+veaEH7mmWYtPY4HSdvYevx9OeK/CGrReesASetdYzlQno0mTUkhOm01nwW9Bm/HvmVgbUH8kqdV+6vbrprDqx+z7BD2lMzoWLLh7a1+VgE/1sSyplrNwn0deP9LtUpVUTWHuQ1ps4a+kUp5ayUcgIOAkeVUm+aO0ghRM5RSjEmYAxPV3ma2QdmM2P/jHvfNFQwHbjesPjshydh7YeQkvFsoebeJVk9qjkjW1dhVeglWn+9iR92nCYlVYYS8wtjHg3VSOsBPAmsAMoDL1k0KiGExVkpK8Y2GktgpUCmB09n7oG5959QphYM2mDYM3nrNzC/I1w/nWFbDrbWvNbWm1WjmuHr7srYJQfpPn0bB8IfnGci8iJjEoGtUsoWQyJYorVOwshZQ0qpDkqpo0qpE0qpMRm8X00ptUMplaCUGp210IUQprJSVnzU+CM6enVk0t5J/HDwh/tPsHOCwCnQYz5EHDOsOQhd+ND2KpYszI/9A5j8vB8XouLpNm0rHyw9KGsP8jhjEsFM4DTgBGxWSlXAMGD8SGnjCdOAjkAN4AWlVI0HTosEXgUmZCFmIYQZWVtZ82nTT2lboS1f7f6K3478lv6kWk/BkC1QspqhRMWS4ZCY8dRRpRTd/Mqx7o0W9GxQgQU7TtPm600s3y9rD/KqTBOB1nqK1rqc1rqTNjgDtDKi7QDghNb6pNY6EfgN6PZA21e01rsA+bogRC6ysbLhi2Zf0NK9JeN3jufv43+nP6loBei7ApqNhn0/w8wWcHH/Q9t0KWTLx0/WYvGwJpRytmfEL4a1B6dl7UGeY8xgsb1S6kWl1LtKqbFKqbHAu0a0XQ44d8/r8LRjQog8yNbalq9bfk2Tck34YPsHLAtblv4ka1to/T/ovRQSY2FOa/hvRoZrDm7z9XBlyfCmfNC1BvvORtFu0mamrDtOQrKsPcgrjHk0tATDN/lkIO6en8xkVKowW/1CpdQgpdRupdTuiIiI7DQhhDCCnbUdk1pOIqBsAO9ve59Vp1ZlfKJXcxiyzVDKetXb8OsLEHftoe1aWyn6NPFi3RstaFujNN+sOUbHSVv4Legs0TflgUBuM6YMdajWulaWG1aqEfCB1rp92ut3ALTWn2Vw7gdArNY607ECWUcghOXdTLrJ0LVDCYkI4esWX9O6QuuMT9TasC3mv++DY3F4apYhSWRi07EIPl5+iBNXYrG1VrTwLklXXzfa1iiNo51NpteLrDN1P4JZwLda6wNZvKkNcAxoDZwHdgEvaq0PZnDuB0giECJPiUuKY9CaQRy6dohJLSfRwqPFw0++uN8wiHztBDR7A1q+Yyhd8Qhaaw6cj2Zp8AWW77/IpZh4Ctla06ZGabr5utHcuyR2NsY8tBDGMDURHAIqA6cwlJlQGDav9zHixp2ASYA1ME9rPV4pNQRDAzOUUmWA3YAzkArEcnfdQoYkEQiRc2ISYxj470COXz/O1Cem0rhc44efnBgHK9+CfT+BewA8PccwwGyE1FRN0OlIloZcYOWBi1y/mYRLIVs61ipDoK8bDSoWl41xTGRqIsjw32Ta7KEcJ4lAiJwVnRBN/9X9OR1zmumtpxNQNuDRFxz4C5a/BigInAw1u2fpfkkpqWw9fpWlIRf49+Al4hJTKFXEns4+ZQn0dcPPw/VuOQxhNHPsUNYUqKK1nq+UKgkU1lqfMnOcRpFEIETOi4yPpP/q/pyPPc+MNjOoW7ruoy+4fhoWDoDwXYaVyR2+ADvHLN/3VmIK649cYWnIeTYciSAxJRWPYoUI9HUj0LccVcsUyd4HKoBM7RGMA/yBqlprb6WUG/Cn1rqJ+UPNnCQCIXLH1VtX6buqLxG3IpjVdhY+JTN5OpySBBs+ha0ToYQ39JhnKFuRTTHxSawOvcTSkAtsO3GVVA1VSxch0M+Nrj5ulC+e9URTkJiaCIKBOsBerXWdtGP7jRkjsARJBELknstxl+m7ui9R8VHMaT+HGsUfLBaQgZMb4e/BcOs6tB8P9QcYCtuZ4GpsAisOXGRJ8AX2nLkOgJ+HK4G+bnTxKUspZ6l++iBTE0GQ1jpAKbVXa103rQrpDkkEQhRMF2Mv0mdVH+KS45jbbi5Vi1XN/KK4q7B4KBz/F6p2gm7TwLGYWeIJv36TZSEXWRpygcMXY7BS0LBicQJ93ehYqywujrZmuU9+Z2oiGA1UAdoCnwH9gV+01lPMHagxJBEIkfvO3ThHn1V9SE5NZl77eVRyrZT5RVrDf9/BmrHgVBKeng2eTc0a14krN1gafIGlIRc4fe2mrFG4hzkGi9sC7dJertZarzVjfFkiiUCIvOF09Gn6ru4LwPz28/F08TTuwgvBsLA/XAsDjwbg7g/l6kI5f8MeymaYEaS1JvR8DEuCz6dboxDo60Zz7xLY21ibfJ/8JFuJQCl1g7slIR78NxMPhAHvaa3XmStQY0giECLvCIsKo9/qfthY2fB9h+/xKOJh3IUJsYZB5FOb4WIIpCQYjjuVhHL1DEmhXF3DT6GiJsWY0RoFZwcbOtYqS6CfGw0LyBoFk3sEGTRoDdQCfs5O+QlTSCIQIm85GnmU/v/2x8nGifkd5uNW2C1rDaQkweVQOL8HwvcY/nn16N33i1dOSwz1wL0elK4NNnbZijWjNQoli9jTubYhKdR5jNcomD0R3NPwYK31zGw3kA2SCITIew5dO8SA1QNwsXfh+w7fU9qptGkNxkfD+b2GpHB+D4Tvhrgrhves7aCMT9ojpXqGn2IVs/xI6b41CkcjSEw2rFHo6uNGn8aej93MI4slgtwgiUCIvGl/xH4GrRlEyUIlmd9hPiUKlTBf41pDdHhaYtht6DlcDIakm4b3CxW9mxRu9x6cihvd/L1rFLaHXaOIgw2fdq9Np9plzfcZcpkkAiFEjth7eS9D1g6hXOFyzG0/l2IO5pkimqGUZIg4fLfHcH6v4bVONbxf1POeR0r+UKY22BbKtNmwiFhe/z2YkPBonqpbjg8Ca+LskP+noEoiEELkmKCLQQxbNwxPZ0/mtp+Li71Lzt084YZhVtLtnsP5vRBz3vCelQ2UrnU3MZTzN4w/WKWvcJqUksq3608wbcMJyjg78M2zvjSoaHwPIy+SRCCEyFHbz29nxPoReBf1Zna72RSxy8WaQDEX05JCWs/hwj7D7moA9i5Qrs79j5SK3B3f2Hv2Oq/9HszZyJsMal6R19t659tpp5IIhBA5btO5TYzaOIoaxWswq+0snGydcjskg9QUuHrsnkdKu+HyIdBpW2e6eEDd3tB4BNgWIi4hmU/+OcyvQWepXtaZSc/55ctid5IIhBC5Yt2Zdbyx6Q18S/ryXZvvcLTNo4XhEm8a1jOc3wMnN8CJteBSHtp9DDW6gVKsPXSZMX/vJyY+mbfaV6VfEy+s8tH6A0kEQohcs+rUKt7e8jb1S9dnauupONjkg2mZpzbDqncM6xsqNIUOn0FZH67GJjBm4X7WHr5Ck8rFmfCML2VdMh+AzgselQhkHzghhEV18OrAJ00+IehSEKM2jiIxJTG3Q8qc/NbWlAAADRRJREFUV3MYtAk6fwNXDsHM5rBsJCXUDWa/7M/nT9Vm39ko2k/czJLg87kdrckkEQghLK5rpa580PgDtp3fxhub3iApJSm3Q8qctQ3U7w+v7oWGQw1bcE75f3v3HhxVecZx/PuQECDBkACBhKxILEaFpggil5iltpYOXupttFDb2vpH661qWxUtVqvjtHZa22LvtWCrLWq94IzTouKNsiGSIjGaAopICyYQQMDEJEBI8vSP92xYYqAJ7OYcPM9nJpO9nN19kkn2t+c973mfScjK3zJnUj5LbojyiRGDufGxam549HUaWo6Bn+kQLAiMMX3ikpMu4fapt7PsvWXcGruVto42v0vqmUG5bmjomgo37fT5efC7UsbsruCJq6Zz08xiltRsZdb9y6nY8L7f1R4RCwJjTJ+Zc8oc5p4xlxc2vcC88nnHThgA5J0MX3kKLn/cnbS26FLSH5vN9RNg8bWlDMpI4/IFldzz97Xs3d/ud7W9YgeLjTF9bmHNQuZXzWfIgCGUjiolWhjlzMIzU3smcjK1tcK//gD//Ilb5mLKVewpvZl7X9nCw69uonjkYObPnsi4Udl+V9rJZg0ZYwJnee1ylv53KeV15ezcuxNBKBleQlmkjBmFMzh12Kn0k4APWjTtgJfvgaqHXce1z97BsqxZ3LJ4DR+0tHLT50/mG9ETA7HMtQWBMSawOrSDdbvWEauNEauLUbOjBkUZOnAoZYVlRCNRSkeVkp0RnE/XH7H1DXj2NthcASNLaPzMPcxdlc1za+qZUjSUn39xApFcf8+hsCAwxhwzdu/dzYotK4jVxlixZQUN+xpIkzQm5E0gGokSLYxSnFscvL4BqrDmadeKs+E9dNxFPFtwLXNf/AAB7rpgPJdMKvStbgsCY8wxqb2jnZr3a4jVxYjVxli3ax0AIzJHEC2MMiMyg2kF04J1xnJrC1T8ynVgQ2mYeDXXbfo05Zv3cG5JPj+8qITcrCNrrHM0LAiMMR8L21u2s6JuBbG6GBVbKmje30x6v3Qmj5xMtDBKNBJlTPaYYOwtNNTCi3dBzRPocaN4KXIt17x5IrmZA7jvsgnMKM7r03IsCIwxHzv7O/ZTvb2689jChg82ABAZHOkcQjoj/wz/l7TYvBKevRW2VtMy4nRuab6cf+ws4GvTT+C2c05lUEbfrGZqQWCM+djb0rSF8rpyYrUxKusr2dO2hwFpA5iSP6UzGCLHRfwprqMD3ngEXrwbmrdTPew8vlF3Ltl5EebPnkhJJPU9GywIjDGhsq99H6vrV7tjC3UxNjVuAqBoSFHnENLpI06nf1ofdx7b2wix++DV39LWL4M/6MX8Zs9Mrv3ceK45a2xKp5laEBhjQm1T46bOvYVV9ato7WglMz2TaQXTiEailBWWkZ+V33cF7XwXlt4Bb/+DHf1HcXvzbHZGZvKL2RMZPSw1B74tCIwxxtOyv4VV9auI1cVYXrucrc1bASjOLe7cWyjOLe6brmrvvuyWu97xFisp4d6OK/jyF87hssmRpB/wtiAwxphuqCobGzZ2HnCu2lZFm7r1j7L6Z5GfmU9+lvsamTWy83r8clKmrba3wWsP0vHyD2FfI39pO5vqT1zD9y8rY9jgAUf//B4LAmOM6YGm1iYq6yvZ3LiZbS3bqG+u7/zauXfnR7bPzsjuDIrE0IhfH5k1koy0Hp4z0LILfeVH6GsP0tgxkD+mzWHypd/lM+OSc4DbtyAQkVnA/UAasEBVf9zlfvHuPxdoAb6uqlWHe04LAmOMH1rbWw8Kh65BUd9ST8O+ho88bujAod0HhXfb8Mzh9O+XcNB621qan7mFrLpy1ncUUj72ZubMuYLMjPSjqt+XIBCRNGA9MBOoBVYBX1LVtQnbnAtcjwuCqcD9qjr1cM9rQWCMCaqW/S0HB0RLPduaD77evL/5oMf0k34MHzj8wJCTFxB5u2rJrvwrJ+2p520mMfTinzL+kxOPuLbDBcHRRczhTQE2qOpGr4jHgAuBtQnbXAg8rC6NVopIjogUqOrWFNZljDEpkdk/k6IhRRQNKTrkNk2tTZ2h0HWP4p3d7xCrjbG3fa/bOC8NKCRN6xmx8nLKXp/EnV9dlPS6UxkEhcB7CddrcZ/6/982hcBBQSAi3wS+CTB69OikF2qMMX1lcMZgxmaMZWzu2G7vV1Ua9jUcFBSb31/P+rVLycsZk5KaUhkE3c196joO1ZNtUNUHgAfADQ0dfWnGGBNMIkLOwBxyBuZwytBTDtxRdmfKXjOVXR9qgeMTrkeALUewjTHGmBRKZRCsAk4SkSIRyQDmAM902eYZ4ApxpgENdnzAGGP6VsqGhlS1TUS+BTyPmz76oKquEZGrvft/DyzBzRjagJs+emWq6jHGGNO9VB4jQFWX4N7sE2/7fcJlBa5LZQ3GGGMOL+CdoY0xxqSaBYExxoScBYExxoScBYExxoTcMbf6qIjsADYd4cOHA+8nsZxkCWpdENzarK7esbp65+NY1wmqmtfdHcdcEBwNEXntUIsu+SmodUFwa7O6esfq6p2w1WVDQ8YYE3IWBMYYE3JhC4IH/C7gEIJaFwS3Nqurd6yu3glVXaE6RmCMMeajwrZHYIwxpovQBIGIzBKRt0Vkg4jc5nc9ACLyoIhsF5F/+11LIhE5XkReEZF1IrJGRG70uyYAERkoIv8SkTe8uu72u6ZEIpImIq+LyN/9riVORP4rIjUiUi0igenx6nUjfFJE3vL+zqYHoKaTvd9T/KtRRL7td10AIvId72/+3yLyqIgMTOrzh2FoqCf9k32qawbQhGvX+Uk/a0kkIgVAgapWichxwGrgogD8vgTIUtUmEekPlAM3qupKP+uKE5HvApOBbFU93+96wAUBMFlVAzUnXkQeAmKqusBbpj5TVT/wu6447z2jDpiqqkd63lKyainE/a2PU9U9IvI4sERV/5ys1wjLHkFn/2RVbQXi/ZN9parLgV1+19GVqm5V1Srv8ofAOlwLUV+p0+Rd7e99BeKTjIhEgPOABX7XEnQikg3MABYCqGprkELAczbwrt8hkCAdGCQi6UAmSW7gFZYgOFRvZPN/iMgYYCJQ6W8ljjf8Ug1sB15Q1UDUBcwH5gIdfhfShQJLRWS11/s7CE4EdgB/8obSFohIlt9FdTEHeNTvIgBUtQ64D9iM6+feoKpLk/kaYQmCHvVGNgcTkcHAU8C3VbXR73oAVLVdVU/DtTWdIiK+D6mJyPnAdlVd7Xct3ThTVScB5wDXecORfksHJgG/U9WJQDMQiON2AN5Q1QXAE37XAiAiubgRjCJgFJAlIl9J5muEJQisN3IveWPwTwGLVHWx3/V05Q0lLANm+VwKwJnABd54/GPAZ0Xkr/6W5KjqFu/7duBp3DCp32qB2oS9uSdxwRAU5wBVqrrN70I8nwP+o6o7VHU/sBgoTeYLhCUIetI/2Xi8g7ILgXWq+nO/64kTkTwRyfEuD8L9g7zlb1Wgqt9T1YiqjsH9bb2sqkn9xHYkRCTLO9iPN/TyecD3GWqqWg+8JyInezedDfg6EaGLLxGQYSHPZmCaiGR6/5tn447bJU1KW1UGxaH6J/tcFiLyKHAWMFxEaoEfqOpCf6sC3CfcrwI13ng8wDyv9aifCoCHvBkd/YDHVTUwUzUDaCTwtHvvIB14RFWf87ekTtcDi7wPZhsJSL9yEcnEzS68yu9a4lS1UkSeBKqANuB1knyGcSimjxpjjDm0sAwNGWOMOQQLAmOMCTkLAmOMCTkLAmOMCTkLAmOMCTkLAhNqIvKLxBUmReR5EVmQcP1nIjLPm77Xm+f9uoj8Opm1GpMqFgQm7CrwztIUkX7AcGB8wv2lwEuqeqkPtRnTJywITNit4MDp+uNxZ95+KCK5IjIAOBXYHe8Z4X3SXywiz4nIOyLyk/gTiciVIrJeRP6JOykvfvsJIvKSiLzpfR/tLZ63UZwcEemIrwMkIjERGdtHP78xFgQm3Ly1eNpEZDQuEF7FrbQ6Hddb4E2gtcvDTgNmAyXAbHGNfAqAu3EBMBMYl7D9r3E9Jz4FLAJ+qartuB4Z44AyXM+HqBc+EVXdkIqf15juWBAYc2CvIB4EryZcr+hm+5dUtUFV9+LWyDkBmAos8xYGawX+lrD9dOAR7/JfcG/8ADHcuvwzgHu928/ArY1lTJ+xIDDmwHGCEtzQ0Ercm3cpLiS62pdwuZ0Da3b1dL2W+HYxIIpbEXQJkINbe2p5z0s35uhZEBjj3uzPB3Z5/Q524d6Up+P2DnqiEjhLRIZ5S3hflnBfBW5VUoAv49oOxh9TCnR4exfVuMXOYkfzwxjTWxYExkANbrbQyi63NfS016+qbgXuwgXHi7iVIuNuAK4UkTdxq7re6D1mH65zXvx1Y8Bx3msb02ds9VFjjAk52yMwxpiQsyAwxpiQsyAwxpiQsyAwxpiQsyAwxpiQsyAwxpiQsyAwxpiQsyAwxpiQ+x+LkuIGH1xn/AAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "labels = ['12 clusters', '6 clusters', '3 clusters']\n", "\n", "ces_fig2, ces_ax2 = plt.subplots()\n", "for data, label in zip(ces_conv2.T, labels):\n", " plt.plot(data, label=label)\n", "ces_ax2.set_xlabel('Window')\n", "ces_ax2.set_ylabel('Jensen-Shannon divergence')\n", "plt.legend()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Using default arguments with dimension reduction ensemble similarity" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "See [dimension_reduction_ensemble_similarity.ipynb](dimension_reduction_ensemble_similarity.ipynb#Calculating-dimension-reduction-similarity-with-default-settings) for an introduction on comparing trajectories via dimensionality reduction." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "dres_conv = encore.dres_convergence(u, # universe\n", " 10, # window size\n", " selection='name CA') # default" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Much like `ces_convergence`, the output is an array of similarity values." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0.53383413],\n", " [0.41463847],\n", " [0.30425993],\n", " [0.24385849],\n", " [0.18008969],\n", " [0.12072474],\n", " [0.05691098],\n", " [0.02318104],\n", " [0. ]])" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dres_conv" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, 'Jensen-Shannon divergence')" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3yV9fn/8deVRQgjbAgjbJlhGVARFzgAB462ghOkIrXuandtbX9t/bbirBtBxbpq3UJRUcA6gATZG5QhK7IhrJDr98c5YIohnIyTO8l5Px+P80juc+5zzhuFc537M83dERGR2BUXdAAREQmWCoGISIxTIRARiXEqBCIiMU6FQEQkxiUEHaC4GjRo4K1atQo6hohIpZKdnf2tuzcs7LFKVwhatWpFVlZW0DFERCoVM1t9rMfUNCQiEuNUCEREYpwKgYhIjFMhEBGJcSoEIiIxToVARCTGqRCIiMS4mCkEObv2c887CzmQlx90FBGRCiVmCsGsr7cy/tOv+eXr89AeDCIi34mZQjA4I43bzz6B12d/wyMfrQg6johIhVHplpgojVsGtGP11j3c/8Ey0uulcHHPZkFHEhEJXEwVAjPj3ku7sX77Xn7+2jzSUpM5qU39oGOJiAQqZpqGDktKiOPJqzJpUa86oyZkszJnd9CRREQCFXOFACA1JZHxw/uQEGdc9+wstuzeH3QkEZHAxGQhAEivn8LT12ayccc+Rk3IZt/BQ0FHEhEJRMwWAoBe6XV54PIeZK/exp3/mkt+voaVikjsielCAKFhpb8a1JF3523gvveXBh1HRKTcxdSooWMZdXobvt6Sy2NTV5JeL4WhfdKDjiQiUm6iekVgZgPNbKmZrTCzXxby+JlmtsPM5oRvd0czTxE5+dOQLpx+QkN+8+YCPlmeE0QMEZFARK0QmFk88CgwCOgMDDOzzoWc+om79wjf/hitPMeTEB/Ho1f0pH2jmtz4wmyWbtwVVBQRkXIVzSuCPsAKd1/l7geAl4EhUXy/UquVnMi44b2pnhTPdc/OYvPOfUFHEhGJumgWgmbA2gLH68L3He0UM5trZpPMrEthL2Rmo8wsy8yycnKi22zTtE51xg3vzbbcA/z4+SxyD+RF9f1ERIIWzUJghdx39PjM2UBLd+8OPAK8WdgLuftT7p7p7pkNGzYs45jf17VZKo8M68mCb3Zw68tzOKRhpSJShUWzEKwDWhQ4bg6sL3iCu+90993h3ycCiWbWIIqZIjagU2PuvqAzHyzaxJ/fWxx0HBGRqInm8NFZQHszaw18AwwFrih4gpk1ATa5u5tZH0KFaUsUMxXL8FNbs3prLuM+/YqW9VO4tm+roCOJiJS5qBUCd88zs5uAyUA8MM7dF5rZ6PDjTwA/AH5iZnnAXmCoV7BdY357fmfWbt3LPe8spHnd6gzo1DjoSCIiZcoq2OfucWVmZnpWVla5vmfugTwuf/ILVubs5tUbTqFrs9RyfX8RkdIys2x3zyzssZhfYiISKUkJPHNtJnVTkrju2Vms37436EgiImVGhSBCjWonM254b/YeOMR1z85i176DQUcSESkTKgTF0KFJLR69shfLN+/mphe/JO9QftCRRERKTYWgmE4/oSF/vrgr05blcPfbC6lsfSwiIkfT6qMlMLRPOqu35vL41JW0qp/CqNPbBh1JRKTEVAhK6K5zO7Bmay5/mbiEFnVTGJSRFnQkEZESUdNQCcXFGWN+2J1e6XW47ZU5zF6zLehIIiIlokJQCsmJ8Tx9TSaNaydz/XNZrN2aG3QkEZFiUyEopfo1qzF+RG/y8p3h42eyI1fDSkWkclEhKANtG9bkyatPZM3WXEa/kM2BPA0rFZHKQ4WgjJzcpj5/+0E3Pl+1hV++Pk/DSkWk0tCooTJ0Sc/mrN6Sy4MfLqdV/RrcMqB90JFERI4r4kJgZjXcfU80w1QFtw5oz5qtudz/wTLS66Vwcc/CNmUTEak4jts0ZGZ9zWwRsDh83N3MHot6skrKzLj30m6c3KYeP39tHjNWVZjtFUREChVJH8EDwHmEN4xx97nA6dEMVdklJcTx5FWZNK9XnVETslmZszvoSCIixxRRZ7G7rz3qrkNRyFKlpKYk8uzwPiTEGdc9O4stu/cHHUlEpFCRFIK1ZtYXcDNLMrM7CTcTSdHS66fw9LWZbNyxj1ETstl3UPVTRCqeSArBaOCnQDNCG9L3CB9LBHql1+X+H/Uge/U27vzXXPLzNaxURCqW444acvdvgSvLIUuVdX63NNZu68i9k5aQXi+Fnw/sGHQkEZEjIhk19JyZ1SlwXNfMxkU3VtVzw+ltGNYnncemruTlmWuCjiMickQk8wi6ufv2wwfuvs3MekYxU5VkZvxxSBfWbcvlN28uoFnd6pzWvmHQsUREIuojiDOzuocPzKwempFcIonxcTx2ZS/aN6rJjS/MZunGXUFHEhGJqBCMAT4zsz+Z2Z+Az4C/RTdW1VUrOZFxw3tTPSme656dxead+4KOJCIx7riFwN2fB34AbAI2A5e6+4RoB6vKmtapzjPX9mbrngOMfC6L3AN5QUcSkRgW6eqjS4DXgbeA3WaWHr1IsSGjeSqPDOvJwvU7uPXlORzSsFIRCUgko4ZuJnQ18AHwLvBe+KeU0tmdG3P3BZ35YNEmfvvmfC1dLSKBiKTT91agg7tr9bQoGH5qazbv2s9jU1eSWj2JXw7SHAMRKV+RFIK1wI5oB4lld53XgR17D/LEtJWkVk/kJ2e2DTqSiMSQSArBKmCqmb0HHFk5zd3vP94TzWwg8BAQD4x193uPcV5v4Avgcnd/LZLgVUlojkFXdu7L4//+s4TU6olccZK6YUSkfERSCNaEb0nhW0TMLB54FDiH0BpFs8zsbXdfVMh5/wdMjvS1q6L4OGPMD7uza99BfvPmfGolJ3Bh96ZBxxKRGBDJWkP3QIl2KOsDrHD3VeHnvwwMARYddd7NwL+B3sV47SopKSGOx688kWvGzeCOV+dQKzmBMzs0CjqWiFRxkYwaOqWEO5Q1I9S/cNi68H0FX7sZcAnwxHEyjDKzLDPLysnJieCtK6/qSfGMvbY37RvVYvQL2WR9vTXoSCJSxUUyj+BBSrZDmRVy39HjIx8EfuHuRS7U7+5PuXumu2c2bFj11+dJrZ7Ic9f1IS21OiOencWi9TuDjiQiVVg0dyhbB7QocNwcWH/UOZnAy2b2NaHZy4+Z2cWRZKrqGtaqxoSRfahZLYFrxs3kq2+L0yonIhK5aO5QNgtob2atzSwJGAq8XfAEd2/t7q3cvRXwGnCju79ZvD9C1dW8bgoTRp5EvjtXjZ3Bxh1al0hEyl7Udihz9zzgJkKjgRYDr7r7QjMbbWajSx45trRrVJPnRvRhx96DXPXMDLbuORB0JBGpYqyyLWuQmZnpWVlZQccod1+s2sI142bSsUktXrz+ZGpW00rgIhI5M8t298zCHjvup4mZPVzI3TuALHd/q7ThJDInt6nPY1f04oYXsrn+uSzGj+hNcmJ80LFEpAqIpGkomVBz0PLwrRtQDxhpZg9GMZsc5ezOjbnvh934fNUWbn7pS/IO5QcdSUSqgEjaF9oB/cNt/pjZ48D7hGYMz49iNinEJT2bs3NvHr9/eyE///c87vtBd+LiChupKyISmUgKQTOgBt8tPFcDaOruh8xs/7GfJtFybd9W7Nh7kPs/WEbt5ER+f2FnzFQMRKRkIikEfwPmmNlUQpPETgf+YmY1gA+jmE2KcHP/dmzPPci4T7+ibkoSt57dPuhIIlJJFVkILPQ1831gIqG1gwz4tbsfnhh2V3TjybGYGb89vxM79x3kgQ+XUbt6AiNObR10LBGphIosBO7uZvamu59IaJtKqUDi4ox7L81g596D3PPOIlKrJ3Jpr+ZBxxKRSiaSUUNfhPcLkAooIT6Oh4f1pG/b+tz12jw+WLQp6EgiUslEUgjOIlQMVprZPDObb2bzoh1MIpecGM9T12TStVkqP31xNp+t/DboSCJSiURSCAYBbYD+wIXABeGfUoHUrJbAs8N707JeCtc/l8W8dduDjiQilcRxC4G7rya0imj/8O+5kTxPyl/dGklMGHkSdWskce24mazYvCvoSCJSCUSyMc3vgV8AvwrflQi8EM1QUnJNUpN5YeRJxMfFcdXYmazblht0JBGp4CL5Zn8JcBGwByA8dLRWNENJ6bRqUIMJI/uQeyCPq8bOIGeX5v2JyLFFUggOeGiJUofQ3sXRjSRloVNabcaP6M2mnfu5ZtxMduw9GHQkEamgIikEr5rZk0AdM7ue0Gzip6MbS8rCiS3r8cTVJ7Ji8y5GPjuLvQci2VhORGJNJJ3F9xHaPezfQAfgbnd/JNrBpGyccUJDHhrak9lrtjH6hWwO5GnFUhH5X5F0Ft8OLHb3u9z9Tnf/oBxySRkanJHGXy7JYNqyHO54dQ6H8ivXZkQiEl2RLDpXG5hsZluBl4HX3F3TVyuZoX3S2bH3IH+dtITa1RP588VdtWKpiACRNQ3d4+5dCO1T3BSYZmZadbQSuuGMttx4ZltenLGGv01eGnQcEakgirPx7WZgI7AFaBSdOBJtd53Xge17D/L41JWkVk9k9Bltg44kIgGLZM/inwCXAw0JdRpf7+6Loh1MosPM+NOQruzce5B7Jy0htXoiw/qkBx1LRAIUyRVBS+A2d58T7TBSPuLjjPt/1IPd+/P49RvzqZ2cyPnd0oKOJSIBOWYfgZnVDv/6N2CNmdUreCufeBItSQlxPH7liZyYXpfbXvmSactygo4kIgEpqrP4xfDPbCAr/DO7wLFUctWT4nlmeG/aNarF6AnZZK/eGnQkEQnAMQuBu18Q/tna3duEfx6+tSm/iBJNqdUTef66PjRJTWbE+Fks3rAz6EgiUs6KahrqVdStPENKdDWsVY0JI/uQkpTA1c/M5Otv9wQdSUTKUVFNQ2PCt0eBGcBThNYYmgE8HP1oUp6a103hhR/34VB+Plc9M4ONO/YFHUlEyklRTUNnuftZwGqgl7tnhjex7wmsKK+AUn7aNarFc9f1YXvuQa5+Zgbb9hwIOpKIlINIVh/t6O7zDx+4+wKgRyQvbmYDzWypma0ws18W8viQ8D7Ic8wsy8z6RR5doqFb8zo8fU0mq7fmMnz8TDbv0pWBSFUXSSFYbGZjzexMMzvDzJ4GFh/vSWYWT6hZaRDQGRhmZp2POm0K0N3dewDXAWOLF1+i4ZS29Xn0il4s3rCL/vdN4+npq7RqqUgVFkkhGAEsBG4FbgMWhe87nj7ACndf5e4HCC1YN6TgCe6+O7zpDUANwpvfSPDO6dyYybefTu9WdfnzxMUMfGg60zXXQKRKimTRuX3u/oC7XxK+PeDukbQXNAPWFjheF77vf5jZJWa2BHiP0FXB95jZqHDTUVZOjj6MykvrBjUYP6IP44Znkp/vXDNuJtc/n8WaLdoHWaQqieSKoKQKW+P4e9/43f0Nd+8IXAz8qbAXcvenwp3VmQ0bNizjmHI8/TuGrg5+PrADn674lrMfmMaY95dqxzORKiKahWAd0KLAcXNg/bFOdvfpQFszaxDFTFJC1RLiufHMdnz0szMZ1LUJj3y0ggFjpvLuvPV817onIpVRNAvBLKC9mbU2syRgKPB2wRPMrJ2Fd0cJT1JLIrTMtVRQTVKTeWhoT/41+hTqpCRx04tfMuzpL1iyUTOSRSqrSJahPgG4i9AqpEfOd/f+RT3P3fPM7CZgMhAPjHP3hWY2Ovz4E8BlwDVmdhDYC1zu+npZKfRuVY93bu7HizPXMOb9pZz/8H+5+uSW3H72CaSmJAYdT0SKwY73uWtmc4EnCC02d6RR2N2zoxutcJmZmZ6VpTXvKpJtew5w/wfL+OeM1dRJSeKu8zrwo8wWxMdpK0yRisLMst09s9DHIigE2eEZxRWCCkHFtXD9Du55exEzv95KRrNU/nBRF05sWTfoWCJC0YUgkj6Cd8zsRjNL034EUpQuTVN55YaTeWhoD3J27eeyxz/jjlfmsHmnZieLVGSRXBF8VcjdHtRS1LoiqBz27M/j0Y9XMPaTr0iMN24Z0J4Rp7YmKSGa4xNE5FhK1TRU0agQVC5ff7uHP727iClLNtOmQQ3uvrAzZ3ZoFHQskZhTqqYhM0s0s1vM7LXw7SYz07AQiUirBjV4Znhvxg/vjQPDx8/ix89lsXqL9jwQqSgiaRoaCyQCz4Xvuho45O4/jnK2QumKoPLan3eI8Z9+zSNTlnMw3xl1WhtuPKstKUnHHcUsIqVU2lFDc929+/HuKy8qBJXfpp37uHfSEt748hvSUpP59eBOXNAtjfDcQhGJgtKOGjpkZm0LvFgbCswnECmuxrWTeeDyHrw2+hTqpiRx80tfMvSpL7RfskhAIrkiGACMB1YRWkiuJTDC3T+Ofrzv0xVB1XIo33l51hrum7yUHXsPhmYnn3MCdVKSgo4mUqWUetSQmVUDOhAqBEvcfX/ZRoycCkHVtD03NDv5hS9Wk1o9kTvP68DQ3umanSxSRsqiEPQFWvG/aw09X1YBi0OFoGpbtH4nf3hnITO/2krXZrW556IunNhS8xdFSqu0w0cnAPcB/YDe4VuhLyZSWp2b1uaVUSfz8LCefLvrAJc9/jm3vzKHTZqdLBI1kYzbywQ6a1VQKS9mxkXdm3J2p0Y89vFKnpq+ivcXbuTmAe25TrOTRcpcJP+iFgBNoh1E5GgpSQnceV4HPrjjdE5pW597Jy1h4IPT+Xjp5qCjiVQpkRSCBsAiM5tsZm8fvkU7mMhhLevXYOy1vRk/ojcAI8bP4oYJWWzYsTfgZCJVQyTDR88o7H53nxaVRMehzuLYdiAvn7H/XcXDU5YTb8Yd53bg2lNakhCv5iKRomjROaly1m7N5XdvLWDq0hy6NK3NXy7JoHuLOkHHEqmwSjtq6FIzW25mO8xsp5ntMjNNAZVAtaiXwvjhvXnsyl7k7NrPxY99yt1vLWDnvoNBRxOpdCK5nv4bcJG7p7p7bXev5e61ox1M5HjMjMEZaXz4szO45uSWTPhiNWePmca789ZT2a50RYIUSSHY5O6Lo55EpIRqJydyz5CuvPXTU2lUuxo3vfglw8fPYs2W3KCjiVQKkXQWP0Ro+OibwJGlJdz99ehGK5z6CKQoeYfyef7z1Yx5fyl5+c4tA9pz/WltNPdAYl5pVx+tDeQC5wIXhm8XlF08kbKTEB/Hdf1a8+HPzuCsDo34++SlnP/wJ8z6emvQ0UQqLI0akiptyuJN3P3WQr7ZvpfLM1vwy0EdqVtDK5tK7CnqiuC4S0yYWTIwEugCJB++392vK7OEIlEyoFNjTmlbn4emLGfsJ1/xweJN/HpwJy7r1Uwb4YiERdI0NIFQH8F5wDSgObArmqFEylJKUgK/GtSJ927pR6v6Kdz5r7kMe/oLVmzeHXQ0kQohkkLQzt1/B+xx9+eA84GM6MYSKXsdm9TmtdF9+cslGSxav5NBD03n/veXsu+gNtyT2BZJITg8Q2e7mXUFUgntTSBS6cTFGVeclM6Un53J+RlpPPzRCgY+OJ1PlucEHU0kMJEUgqfMrC7wO+BtYBGhSWYilVbDWtV4cGhPXhh5EmbG1c/M5JaXvmTzLu17ILEnqqOGzGwg8BAQD4x193uPevxK4Bfhw93AT9x9blGvqVFDUtb2HTzEY1NX8sTUlVRLjOMXAztyRZ904rRNplQhpVp0Lrxf8WV8f6vKPx7nefHAMuAcYB0wCxjm7osKnNMXWOzu28xsEPAHdz+pqNdVIZBoWZmzm9++sYDPV22hZ3od/nJJBp3StJqKVA2lnVD2FjAEyAP2FLgdTx9ghbuvcvcDwMvh1znC3T9z923hwy8IjUgSCUTbhjV58fqTuP9H3Vm9JZcLHvkvf5m4mNwDeUFHE4mqSLaqbO7uA0vw2s2AtQWO1wFFfdsfCUwq7AEzGwWMAkhPTy9BFJHImBmX9mpO/46NuHfSEp6avor35m3gnou6cHbnxkHHE4mKSK4IPjOzkgwXLayBtdB2KDM7i1Ah+EVhj7v7U+6e6e6ZDRs2LEEUkeKpk5LEvZd141+jT6FGtXh+/HwWN0zIYv127YomVc8xC4GZzTezeUA/YLaZLTWzeQXuP551QIsCx82B9YW8TzdgLDDE3bcUL75IdPVuVY93bz6NXwzsyLRlOZxz/zTGfrKKvEP5QUcTKTPH7Cw2s5ZFPdHdVxf5wmYJhDqLBwDfEOosvsLdFxY4Jx34CLjG3T+LJLA6iyUoa7fmcvdbC/g4vCvany/JoId2RZNKoqSdxTnAendfHf7QTwYuBU48XhEAcPc84CZgMrAYeNXdF5rZaDMbHT7tbqA+8JiZzTEzfcJLhdWiXgrjwruifbt7P5doVzSpIoq6IpgOjHT35WbWDpgJ/BPoDMx091+VX8zv6IpAKoJd+w4y5v1lPPf51zSsWY27L+zM+RlpWshOKqySXhHUdffl4d+vBV5y95uBQWg/AolxtZIT+cNFXb63K9pX30YyslqkYimqEBS8VOgPfAAQnhOgnjIRoFvzOrx546n8/sLOZK/exrkPTOOvExezS81FUokUVQjmmdl9ZnY70A54H8DM1DsmUkBCfBwjTm3NR3eewcU9mvHk9FWcdd9UXp21lvz8yrXxk8SmogrB9cC3hJaWONfdD+8E3hm4L8q5RCqdRrWS+fsPu/PWT08lvV4KP//3PIY8+ilZ2iZTKrhiLTpnZr3cfXYU8xyXOoulMnB33p67nr9OXMLGnfu4qHtTfjmoI03rVA86msSo0q41VNDYMsgjUuWZGUN6NOOjO8/glv7tmLxwI/3HTOXBD5ex94A2wpGKpbiFQGPjRIohJSmBO87twJSfncGATo158MPlDBgzlXfmrieaS8CLFEdxC8E9UUkhUsU1r5vCo1f04pVRJ1MnJYmbX/qSHz35OQu+2RF0NJHI+gjMrBnQkv/dj2B6FHMdk/oIpLI7lO+8mrWW+yYvZWvuAS7PbMGd53WgQc1qQUeTKqyoPoLjLkNtZv8HXE5oi8rDjZsOBFIIRCq7+DhjWJ90Bmek8ciU5Tz72de8N28Dtwxoz7V9W5GUUNwLdZHSiWSHsqVAN3ffXz6RiqYrAqlqVubs5s/vLeajJZtp3aAGv7ugE2d1aKTlKqRMlXbU0CogsWwjichhbRvWZNzw3owf0RszuO7ZLIaPn8WKzbuCjiYxIpIdynKBOWY2BThyVeDut0QtlUgMOqtDI/q1a8Dzn6/mwQ+XMfDBT7jmlFbcOqA9qSn6LibRE0kheDt8E5EoS4yPY2S/1lzcoyljPljG+M++4s053/Czc09gaO904uPUXCRlL9JRQ9WBdHdfGv1IRVMfgcSShet38Md3FjHjq610bFKLuy/sTN+2DYKOJZVQqfoIzOxCYA7wn/BxDzPTFYJIOejSNJWXR53MY1f2Yte+PK54egY/eSGbtVtzj/9kkQhF0ln8B6APsB3A3ecAraOYSUQKMDMGZ6Qx5WdncOe5JzB1aQ4D7p/GfZOXsmd/XtDxpAqIpBDkufvR0x81N16knCUnxnNT//Z8fOeZnJ+Rxj8+XkH/MVN548t1Wu5aSiWSQrDAzK4A4s2svZk9AkS00byIlL0mqck8cHkP/v2TvjSpncztr8zlsic+Y87a7UFHk0oqkkJwM9CF0NDRl4CdwG3RDCUix3diy7q8ceOp3PfD7qzbtpeLH/2UO16dw6ad+4KOJpVMcfcjiAdquPvO6EUqmkYNiXzf7v15PPrxCp755CsS4o2fntWOkf1ak5wYH3Q0qSBKO2roRTOrbWY1gIXAUjO7q6xDikjJ1ayWwC8GduSDO06nX7sG/H3yUs55YBr/WbBRy13LcUXSNNQ5fAVwMTARSAeujmoqESmRlvVr8NQ1mbww8iSqJ8Yz+oVsrhw7g+WbtFyFHFskhSDRzBIJFYK33P0gGjUkUqH1a9+Aibecxp+GdGHRhp0MfvgT/j55CfsOanc0+b5ICsGTwNdADWC6mbUk1GEsIhVYQnwcV5/Siil3nMFF3Zvx6McrOfeB6UxblhN0NKlgitVZfORJZgnuHshMFnUWi5TM5yu38Js357MqZw8XdEvj7gs606h2ctCxpJwU1VkcyX4E1YDLgFb87w5lfyzDjBFTIRApuf15h3hy2ir+8fEKqsXH8fOBHbjipJZazC4GlHY/greAIUAesKfATUQqmWoJ8dwyoD2Tbzud7i3q8Lu3FnLp45+xcL32To5lkVwRLHD3riV6cbOBwENAPDDW3e896vGOwHigF/Abd7/veK+pKwKRsuHuvD13PX96dxHbcg8yom8rbj/nBGpUi2R1eqlsSntF8JmZZZTgTeOBR4FBQGdgmJl1Puq0rcAtwHELgIiULTNjSI9mTLnjTC7v3YKx//2Kc+6fxvsLNwYdTcpZJIWgH5BtZkvNbJ6ZzTezeRE8rw+wwt1XufsB4GVCTUxHuPtmd58FHCx2chEpE6kpifzlkgz+/ZO+1K6eyKgJ2Vz/fBbrt+8NOpqUk0iuAQeV8LWbAWsLHK8DTirJC5nZKGAUQHp6egnjiEhRTmxZl3du7se4/37Fgx8u5+z7p3HHOScwvG8rEuIj+c4oldVx/++6+2qgBdA//HtuJM8DChuGUKKJaO7+lLtnuntmw4YNS/ISIhKBxPg4bjijLe/ffjont6nP/3tvMRf941OtbFrFRbLW0O+BXwC/Ct+VCLwQwWuvI1RADmsOrC9uQBEpfy3qpfDMtZk8cVUvtu45wCWPfcrv3lzAzn1qxa2KIvlmfwlwEeEho+6+HqgVwfNmAe3NrLWZJQFDAW1xKVJJmBkDu6bx4c/OYHjfVvxzxmoGjJnGO3PXayG7KiaSQnDAQ//XHSC8CulxhWce3wRMBhYDr7r7QjMbbWajw6/VxMzWAXcAvzWzdWZWuyR/EBGJjprVEvj9hV1466f9aFI7mZtf+pJrx89izRbtm1xVRDKP4E6gPXAO8FdgJPCiuz8c/Xjfp3kEIsE5lO9M+Pxr7nt/GQcP5XPLgPZcf1obkhLUmVzRlWqJifALnAOcGz6c7O4flmG+YlEhEAnexh37+OO7C5k4f5jWvWAAAAyYSURBVCPtG9Xkz5dk0Kd1vaBjSRFKVAjMbBffjfI5egTQPmAlodnAU8oqaCRUCEQqjo+WbOJ3by7km+17+VFmc341qBN1ayQFHUsKUVQhOOY8Anc/ZodweNZwV+Cf4Z8iEoP6d2zMyXfU56Epy3nmk6/4cPFmfj24E5f1aoaZFrKrLErUsOfuh9x9LvBIGecRkUomJSmBXw3qxLu39KN1gxrc+a+5DHv6C1Zs3h10NIlQqXp43P3JsgoiIpVbxya1+dcNp/DXSzNYtH4ngx6azv3vL9WuaJWAuvpFpMzExRnD+qTz0Z1nckG3pjz80QoGPjid/y7/NuhoUgQVAhEpcw1qVuOBy3vwwsiTMDOuemYGt778JTm79gcdTQqhQiAiUdOvfQMm3Xoatwxoz6T5GxkwZir/nLGa/HzNTK5IVAhEJKqSE+O545wTmHTbaXRuWpvfvLGAHzzxGUs27gw6moSVaPP6IGkegUjl5e68Pvsb/jxxMTv2HuTsTo0YnJFG/46NqJWcGHS8Kq1E8whERMqamXHZic3p37ERj368gnfnbWDywk0kxcdx+gkNGNQ1jbM7Nya1uopCedIVgYgEJj/f+XLtdibO38Ck+RtYv2MfifHGqe0aMLhrGud0bqyZymWk1GsNVSQqBCJVk7szd90OJs3fwMQFG1i7dS/xcUbftvUZ1DWNc7s0pkHNakHHrLRUCESkUnF3Fq7fycT5G5g4fwNfb8klzuCk1vUZnNGE87o0oVHt5KBjVioqBCJSabk7SzbuYtL8Dbw3fwMrc/ZgBr1b1mNQRhMGdm1CWmr1oGNWeCoEIlJlLN+0i4nzNzJpwQaWbNwFQK/0OgzOSGNg1yY0r5sScMKKSYVARKqklTm7+c+CjUycv4GF60PzEro3T2VQRhqDu6aRXl9F4TAVAhGp8lZv2cOkBRuZNH8Dc9ftAKBL09oMzkhjUNcmtGlYM+CEwVIhEJGYsnZrbuhKYcEGvlyzHYCOTWoxOCONwRlNaNfomNutVFkqBCISs9Zv38t/FoT6FLJWb8Md2jeqGWo+ymhCh8a1YmITHRUCERFg0859TF4Y6lOY+dVW8h3aNKjBoIwmDOqaRpemtatsUVAhEBE5Ss6u/by/aCOT5m/k81VbOJTvpNdLYVBGE05r15DuLVKr1PpHKgQiIkXYuucAHyzayMT5G/l0xbfk5TtxBic0rsWJLeseuaXXS6m0VwwqBCIiEdq17yBz1m4ne/U2sldvY86a7ezanwdAg5pJ9Ez/rjBkNEslOTE+4MSR0eqjIiIRqpWcyGntG3Ja+4YAHMp3lm/exezVoeIwe802Pli0CYDEeKNz01ROLFAcmqRWvqUvdEUgIlJMW3bv58s128leE7pqmLt2O/vz8gFomppMrwLNSZ3SapMYH/weYLoiEBEpQ/VrVuPszo05u3NjAA4eymfxhp1HmpNmr97Gu/M2AJCcGEe35nVChSG9Lr1a1qVeBVtaW1cEIiJRsGHH3iPNSdlrtrHwmx3khfdqbt2gBr0KNCe1b1STuLjodkIH1llsZgOBh4B4YKy733vU4xZ+fDCQCwx399lFvaYKgYhURvsOHmL+Nzv+56phy54DANSqlkCP9DpHCkOPFnXKfOhqIE1DZhYPPAqcA6wDZpnZ2+6+qMBpg4D24dtJwOPhnyIiVUpyYjy9W9Wjd6t6QGh57TVbc48UhuzV23hoynLcwQw6NK4V6msIXzm0rB+9oavR7CPoA6xw91UAZvYyMAQoWAiGAM976LLkCzOrY2Zp7r4hirlERAJnZrSsX4OW9Wtwaa/mQGjo6ty1O440J70zdz0vzlgDQL0aSfzkjLZcf3qbMs8SzULQDFhb4Hgd3/+2X9g5zYD/KQRmNgoYBZCenl7mQUVEKoJayYn0a9+Afu0bAKE9nVfk7D5yxdA4SkNTo1kICruGObpDIpJzcPengKcg1EdQ+mgiIhVfXJxxQuNanNC4FsP6RO9LcDQHt64DWhQ4bg6sL8E5IiISRdEsBLOA9mbW2sySgKHA20ed8zZwjYWcDOxQ/4CISPmKWtOQu+eZ2U3AZELDR8e5+0IzGx1+/AlgIqGhoysIDR8dEa08IiJSuKjOLHb3iYQ+7Ave90SB3x34aTQziIhI0YJfAENERAKlQiAiEuNUCEREYpwKgYhIjKt0q4+aWQ6wuoRPbwB8W4ZxykpFzQUVN5tyFY9yFU9VzNXS3RsW9kClKwSlYWZZx1p9L0gVNRdU3GzKVTzKVTyxlktNQyIiMU6FQEQkxsVaIXgq6ADHUFFzQcXNplzFo1zFE1O5YqqPQEREvi/WrghEROQoKgQiIjEuZgqBmQ00s6VmtsLMfhl0HgAzG2dmm81sQdBZCjKzFmb2sZktNrOFZnZr0JkAzCzZzGaa2dxwrnuCzlSQmcWb2Zdm9m7QWQ4zs6/NbL6ZzTGzrKDzHBbelvY1M1sS/nt2SgXI1CH83+nwbaeZ3RZ0LgAzuz38d36Bmb1kZmW6VVlM9BGYWTywDDiH0GY4s4Bh7r6oyCdGP9fpwG5C+zZ3DTJLQWaWBqS5+2wzqwVkAxdXgP9eBtRw991mlgj8F7jV3b8IMtdhZnYHkAnUdvcLgs4DoUIAZLp7hZocZWbPAZ+4+9jwfiUp7r496FyHhT8zvgFOcveSTmAtqyzNCP1d7+zue83sVWCiuz9bVu8RK1cEfYAV7r7K3Q8ALwNDAs6Eu08Htgad42juvsHdZ4d/3wUsJrSXdKA8ZHf4MDF8qxDfZMysOXA+MDboLBWdmdUGTgeeAXD3AxWpCIQNAFYGXQQKSACqm1kCkEIZ7+QYK4WgGbC2wPE6KsAHW2VgZq2AnsCMYJOEhJtf5gCbgQ/cvULkAh4Efg7kBx3kKA68b2bZZjYq6DBhbYAcYHy4KW2smdUIOtRRhgIvBR0CwN2/Ae4D1gAbCO3k+H5ZvkesFAIr5L4K8U2yIjOzmsC/gdvcfWfQeQDc/ZC79yC0v3UfMwu8Sc3MLgA2u3t20FkKcaq79wIGAT8NN0cGLQHoBTzu7j2BPUCF6LcDCDdVXQT8K+gsAGZWl1ALRmugKVDDzK4qy/eIlUKwDmhR4Lg5ZXxpVdWE2+D/DfzT3V8POs/Rwk0JU4GBAUcBOBW4KNwe/zLQ38xeCDZSiLuvD//cDLxBqJk0aOuAdQWu5l4jVBgqikHAbHffFHSQsLOBr9w9x90PAq8DfcvyDWKlEMwC2ptZ63C1Hwq8HXCmCivcKfsMsNjd7w86z2Fm1tDM6oR/r07oH8iSYFOBu//K3Zu7eytCf7c+cvcy/cZWEmZWI9zZT7jp5Vwg8BFq7r4RWGtmHcJ3DQACHYhwlGFUkGahsDXAyWaWEv63OYBQv12ZieqexRWFu+eZ2U3AZCAeGOfuCwOOhZm9BJwJNDCzdcDv3f2ZYFMBoW+4VwPzw+3xAL8O70EdpDTgufCIjjjgVXevMEM1K6DGwBuhzw4SgBfd/T/BRjriZuCf4S9mq4ARAecBwMxSCI0uvCHoLIe5+wwzew2YDeQBX1LGS03ExPBRERE5tlhpGhIRkWNQIRARiXEqBCIiMU6FQEQkxqkQiIjEOBUCiWlm9kDBFSbNbLKZjS1wPMbMfh0evlec1x1uZv8oy6wi0aJCILHuM8KzNM0sDmgAdCnweF9girv/IIBsIuVChUBi3ad8N12/C6GZt7vMrK6ZVQM6AdsO7xkR/qb/upn9x8yWm9nfDr+QmY0ws2VmNo3QpLzD97c0sylmNi/8Mz28eN4qC6ljZvmH1wEys0/MrF05/flFVAgktoXX4skzs3RCBeFzQiutnkJob4F5wIGjntYDuBzIAC630EY+acA9hArAOUDnAuf/g9CeE92AfwIPu/shQntkdAb6Edrz4bRw8Wnu7iui8ecVKYwKgch3VwWHC8HnBY4/K+T8Ke6+w933EVojpyVwEjA1vDDYAeCVAuefArwY/n0CoQ9+gE8Irct/OvDX8P29Ca2NJVJuVAhEvusnyCDUNPQFoQ/vvoSKxNH2F/j9EN+t2RXpei2Hz/sEOI3QiqATgTqE1p6aHnl0kdJTIRAJfdhfAGwN73ewldCH8imErg4iMQM408zqh5fw/mGBxz4jtCopwJWEth08/Jy+QH746mIOocXOPinNH0akuFQIRGA+odFCXxx1345I9/p19w3AHwgVjg8JrRR52C3ACDObR2hV11vDz9lPaOe8w+/7CVAr/N4i5Uarj4qIxDhdEYiIxDgVAhGRGKdCICIS41QIRERinAqBiEiMUyEQEYlxKgQiIjHu/wMaKivcDETsUAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "dres_fig, dres_ax = plt.subplots()\n", "plt.plot(dres_conv)\n", "dres_ax.set_xlabel('Window')\n", "dres_ax.set_ylabel('Jensen-Shannon divergence')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Comparing different dimensionality reduction methods" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Again, you may want to compare the performance of different methods.\n", "\n", "Principal component analysis uses singular value decomposition to project data onto a lower dimensional space. [(See the scikit-learn user guide for more information.)](https://scikit-learn.org/stable/modules/decomposition.html#pca)\n", "\n", "The method provided by MDAnalysis.encore accepts any of the keyword arguments of [sklearn.decomposition.PCA](https://scikit-learn.org/stable/modules/generated/sklearn.decomposition.PCA.html) *except* `n_components`. Instead, use `dimension` to specify how many components to keep." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "pc1 = drm.PrincipalComponentAnalysis(dimension=1,\n", " svd_solver='auto')\n", "pc2 = drm.PrincipalComponentAnalysis(dimension=2,\n", " svd_solver='auto')\n", "pc3 = drm.PrincipalComponentAnalysis(dimension=3,\n", " svd_solver='auto')" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "(9, 3)" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dres_conv2 = encore.dres_convergence(u, # universe\n", " 10, # window size\n", " selection='name CA',\n", " dimensionality_reduction_method=[pc1, pc2, pc3]\n", " )\n", "dres_conv2.shape" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Again, the size of the subspace you choose to include in your similarity comparison, affects the apparent rate of convergence over the trajectory." ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzdeXhNV/fA8e/KzU1CBiESQ0JiJoYYYop5KlpKVYuaS1W1qBqqv85zq9UaWlO1hlaLavu21aJmiiLmeQohhESMCZFp//44oaHBJXfIlf15nvNw7z337MX71srZe5+1RCmFpmmalne5ODoATdM0zbF0ItA0TcvjdCLQNE3L43Qi0DRNy+N0ItA0TcvjXB0dwL0qXLiwCgkJcXQYmqZpTmXLli1nlVL+2X3mdIkgJCSEyMhIR4ehaZrmVEQk+naf6akhTdO0PE4nAk3TtDzOpolARNqIyAEROSwio7P5fKSIbM88dotIuogUsmVMmqZp2s1stkYgIibgS6AVEANsFpHflFJ7r5+jlPoE+CTz/PbAMKXUOVvFpGmadiepqanExMSQnJzs6FDum4eHB0FBQZjNZou/Y8vF4jrAYaVUFICIzAU6AHtvc3434AcbxqNpmnZHMTExeHt7ExISgog4Opx7ppQiISGBmJgYSpUqZfH3bDk1FAicyPI6JvO9/xCR/EAb4KfbfD5ARCJFJDI+Pt7qgWqapgEkJyfj5+fnlEkAQETw8/O75zsaWyaC7P4mb1fqtD2w7nbTQkqpaUqpcKVUuL9/tttgNU3TrMJZk8B19xO/LaeGYoASWV4HAaduc25XbDwtFHUhikXHFhHoFUigVyAlvEvgn88fk4vJlsNqmqblerZMBJuBciJSCjiJ8Y/9U7eeJCIFgCZADxvGwsHzB5m6Yyoqy02Jq4vrjcQQ5BVEoHfm772DCPIKwsfNx+l/OtA0zbk8/fTTLFy4kICAAHbv3g1Anz59WL16NT4+Ply9epV69erx4YcfEhiY7Wz7PbNZIlBKpYnIC8ASwAR8o5TaIyIDMz+fknnqY8BfSqkkW8UC0KZUG1qUbEFsUiwxl2OISYzhZOJJYi4bv+5J2MPFaxdv+o632ZtA78wkkZkgAr0CbyQMd5O7LUPWNC0P6tOnDy+88AK9evW66f1PPvmEzp07o5Ri3LhxNGvWjN27d+Pm5pbjMW1aYkIp9Sfw5y3vTbnl9Uxgpi3jACAxDvOBRZSs2I6SgRHZn5KSeCM5xCTG3EgSRy4eYe3JtVxLv3bT+QH5Am4khxtJIvP3AfkDcBH9vJ6mafemcePGHDt27LafiwjDhg3jl19+YdGiRXTo0CHHYzpdraH7dugv+H0ILBwGIQ0gtANUbA/eRW6c4uXmRYVCFahQqMJ/vp6hMki4mnAjQcQkxnDy8klOJp5k85nNLIxaeNO0k9nFfFNiuDVZFHAvYJc/tqZp9+ft3/ew99Qlq14ztLgPb7avbJVr1axZk/379+tEcE+qd4ciVWDfb7D3V/hjOPwxAkrWM5JCpfZQIOi2X3cRF/zz++Of358aATX+83lKegqxSbGcvHzSSBaZiSImMYbdCbv/O+3k5k2QV5CxHuEdRNXCVanuXx3//HpXlKZpd2fNfvN5JxGIQPHqxtH8dYjb929SWDzaOALDIfRRqPQoFLL8YQwAN5MbwT7BBPsEZ/v55ZTLnEw8+W+iyLyrOHzhMKtOrCI1IxWAQK9AqgdUp4Z/DaoHVKesb1m9s0nTHMBaP7nbyrZt22jRooVVrpV3EkFWIlAk1Diajoazh2Hfr0ZSWPqGcRStZiSF0I5QuFyOh/R286ZioYpULFTxP5+lpqey79w+tsdtZ3v8djbGbuSPqD8A8DR7UrVwVWoE1KC6f3Wq+lfF2807x/FomuaclFJMnDiR2NhY2rRpY5VrijVvL+whPDxc2bQfwfljsO93IynEbDbe869kTB+FPgoBoUYisSGlFKeSTrEtbhvb47azI34HB88fJENlIAjlCpajun91qgcYR5BXkN7mqmlWsG/fPipVquTQGLp168aqVas4e/YsRYoU4e2332bt2rU3to9euXLlxvbRoKDsp7Oz+3OIyBalVHh25+tEcCcXTxpJYd9vEL0eUFCozL9JoVh1myeF6xJTEtl1dhfb47ezPW47O+N3kpiaCICfh58xnRRQgzD/MEL9QnEz5XxLmablNbkhEViDTgS2cvkM7F9oJIWja0Glg29JYz0htIOxvuBiv+2i6RnpHLl4xJhOitvOtrhtxCTGAODm4kblwpWp7l+dsIAwqvtXxy+fn91i0zRnpROBk3BYIsgqKQEO/GlMH0WtgoxU8C5u7DwK7WDsRHLAAu/Zq2fZEbfDmFKK387ehL03FqFLepe8MZVU3b86ZXzL6OccNO0WOhE4iVyRCLK6egEOLjGSwuFlkH4NPAOgUjvjbiGkEZgcsyZ/Lf0aexP23rhr2B6/nXPJRl0/b7M31QKq3VhrqFa4GvnN+R0Sp6blFnk1EeTNXUPWlM8XwroYx7VE48G1vb/CjrkQ+Q3kKwgVH4FKHaB0U3C139y9u8mdGgE1bjz3oJTixOUTN+4YtsdtZ9LJSSgULuJChYIVCPMPM3YoBVSnmGcxvQitaXmAviOwlZQrcGQ57P0NDi6Ga5fAvQBUaGNMH5VpDuZ8jo6SSymX2Bm/88Zdw86zO7madhWAgPwBNA1qSquQVoQXCcfVRf/coD3Y9B2BZl1u+Y01g0rtIe2asZaw9zc48AfsnAdmTyj/kJEUyrYCdy+HhOnj5kPDwIY0DGwIQFpGGofOH2Jb3DYiz0Tye9TvzD84H193X5qXbE7Lki2pV6weZpPlbfA0Tcvd9B2BvaWnwrG1RlLYvxCS4sHkBoUrQEClzCPU+LVACbvuRMrO1bSrrDu5jqXRS1kds5qk1CS8zd40LdGUlsEtiSgegYerh0Nj1DRryQ13BCdOnKBXr16cPn0aFxcXBgwYwNChQ++pFLW+I8jtTGZjWqhMc3hkLBzfYKwrnNlrPKuwa/6/57p5gX/Fm5NDQCh4Bdjt+YV8rvloGdySlsEtSUlPYcOpDSyNXsrKEyv5Pep38rnmo0lQE1oGt6RRYCO94KxpOeTq6srYsWOpWbMmly9fplatWrRq1QqwXSlqnQgcycUEIQ2N47rkixC3H+L2GvWQ4vYaW1W3ffvvOfkKZUkM15NERWNh2obcTG40KdGEJiWakJqRyubYzSw9vpQVx1ew+Nhi3E3uNAxsSMvgljQJaqJLYWjafShWrBjFihUDwNvbm0qVKnHy5MmbzrF2KWqdCHIbjwJQsq5xZJUYf3NyiNtn7ExKufzvOd7F/zu95F/RWK+wMrOLmYjACCICI3it7mtsjdvK0uilLI9ezvLjyzG7mKlfvD4tS7akecnmuuy25nwWjYbTu6x7zaJVoe1HFp9+7Ngxtm3bRt26dZkzZ85/PrdWKWqdCJyFlz94NYHSTf59Tym4GHNzcojbC5v+Np5nAECgYMjNdxBFKhulMqy0ldXkYqJ20drULlqb0XVGszN+J0ujl7IsehlrYtbw9oa3qVO0Di2DjaRQOF9hq4yraQ+yxMREHn/8ccaNG4ePj0+251hrjVcnAmcmAr4ljKP8Q/++n5EO547+9w7i4GKjNAaAi9moqnrrHYRvSI4WqF3E5cYTzCPCR7A3YS9Lo5eyNHop7/7zLu9vfJ+aATWNdYeSLSniWeTuF9U0R7iHn9ytLTU1lccff5zu3bvTqVOn255nrVLUOhE8iFxMULiscYQ++u/7adfg7KGbk0NMJOz+6d9zXPMZ6w0BocadQ5XON3VxuxciQuXClalcuDJDaw7l4PmDLDu+jGXRy/ho00d8tOkjwvzDaBXcipbBLQn0sk4jbk1zZkop+vXrR6VKlXjppZdue441S1Hr7aOa8UR0/IH/3kEknja2toZ1hfqDwb+81YaMuhjFsmgjKew7tw+AUL9QIymUbElIgRCrjaVplsoN20f//vtvGjVqRNWqVXHJvDv/4IMPmD9/vsWlqHNVrSERaQOMB0zAdKXUf+61RKQpMA4wA2eVUk1uPScrnQjs6Oxh+OdL2DbHWHOo8DA0GGoU1bOiE5dP3EgKO8/uBKBcwXK0KtmKVsGtKONbRpe60OwiNyQCa8g1iUBETMBBoBUQA2wGuiml9mY5xxdYD7RRSh0XkQClVNydrqsTgQMkxsOmabD5K7h6HoLqQIMhRmKwcpXV00mnWRa9jKXRS9kWtw2FIsQnhFbBRlKoWKiiTgqazehEYGUiUh94SynVOvP1KwBKqQ+znDMIKK6Ues3S6+pE4EApScbdwYYv4EK0sfMo4gUI62aTuklnr55lefRylh5fSuTpSNJVOkFeQbQKacXTlZ/G18PX6mNqeVteTQS2rF8QCJzI8jom872sygMFRWSViGwRkV7ZXUhEBohIpIhExsfH2yhc7a7cPKHuABi8FTrPAA8fWDgMPq8Cq8fAlXNWHa5wvsJ0qdiF6Q9NZ+WTK3k74m1CCoTw7Z5veWLhE+yI32HV8TQtr7JlIsju/v3W2w9XoBbwCNAaeF1E/rMiqZSappQKV0qF+/v7Wz9S7d6YXKFKJ3hmJfReCIE1YeX78Hll+HOk0ffZygp6FKRTuU5MbjmZ7x75DpOY6LO4D9/t/c5qe6k1La+yZSKIAUpkeR0EnMrmnMVKqSSl1FlgDRBmw5g0axKBUo2g+4/w3AYI7QiRM2BCDfixL5zaZpNhK/tVZl67eTQKbMTHmz9m+OrhXM76hLWmaffElolgM1BOREqJiBvQFfjtlnN+BRqJiKuI5AfqAvtsGJNmK0VC4bHJ8OJOqP+C0a1tWlOY2Q4OLTWegraiAu4FGN9sPCPCR7Di+Aq6LuzK/nP7rTqGpuUVNksESqk04AVgCcY/7vOVUntEZKCIDMw8Zx+wGNgJbMLYYrrbVjFpduBTHB56F4bthlbvQsIRmNMZJkfA9u8hLcVqQ4kIvSv3ZkabGSSnJ9P9j+78dPAnPVWkObXk5GTq1KlDWFgYlStX5s033wSgT58+lCpVirCwMMqXL0+vXr3+U4zuflm8a0hEPJVSSVYZNQf0riEnk5ZiPLm8foLxoJp3cag3EGr1MQrsWcm55HOMXjOaDbEbaF+6Pa/Ve02XxNbuWW7YNaSUIikpCS8vL1JTU2nYsCHjx49nypQptGvX7qYy1JMnT862DLXVdw2JSISI7CVzykZEwkRk0n3/KbW8xdUNqneD59ZD95+MshdL3zB2Gv31Oly6ddno/hTyKMTklpMZVH0QC6MW0v3P7kRdiLLKtTXNnkQELy+jY2Fqaiqpqan/eXbmehnqokWLsmjRohyPaUmtoc8xdvT8BqCU2iEijXM8spa3iEC5lsZxahusm2A8j/DPZKj6BEQMNtYZcsDkYuK5sOeo7l+d0WtH0/WPrrxV/y0eLv2wlf4QWl7y8aaPrb7uVLFQRV6u8/Jdz0tPT6dWrVocPnyY559/nrp16zJ58uT/nGetMtQWrREopU7c8lZ6jkbV8rbiNeCJGTBkG4Q/DXv/B5Prw3ed4eiaHC8s1y9en/nt5lOpUCVeXvsy7/3zHtdulOXWtNzPZDKxfft2YmJi2LRpE7t3Z790as8y1CdEJAJQmbt/hqB39mjWUDAEHh4DTUfD5q9h01SY1R6KVTdKWFTqYDyzcB+KeBZheuvpTNw2kRm7Z7Azfidjm46lhHeJu39Z08Cin9xtzdfXl6ZNm7J48eJsP7dWGWpL7ggGAs9jPBUcA1TPfK1p1pG/EDQZCS/ugnbjICURFjwNE2vCxmlGaYv7YHYx81Ktl5jQbAIxiTF0+b0LK46vsHLwmmZd8fHxXLhwAYCrV6+ybNkyKlaseNM5SikmTJhgtTLUd00ESqmzSqnuSqkiSqkApVQPpVRCjkfWtFuZ80F4X3h+M3SZA95FYdFI44nlFe8bxe/uQ7OSzZjfbj4lfUoydOVQPt38KakZqVYOXtOsIzY2lmbNmlGtWjVq165Nq1ataNeuHQAjR468sX108+bNrFy5MseN68GC7aMiMgsYqpS6kPm6IDBWKfV0jke/D3r7aB5z/B9jYfnAn+DqbhS4ixgMfmXu+VIp6SmM2TyGeQfmUSOgBmMaj6GoZ1EbBK05q9ywfdQabFF0rtr1JACglDoP1MhRlJpmqZL1oNv38MJmqNbFeCjti3DYteCeL+VmcuO1eq8xpvEYDpw7wJO/P8n6k+ttELSmORdLEoFL5l0AACJSCN3iUrO3wuXg0QnGE8sl6sKvL0Dszvu6VNtSbZnbbi5++fwYuGwgX27/kvQMvRFOy7ssSQRjgfUi8q6IvIvRSGaMbcPStNvwCoAnZxsLzHO7Q9L9LVeVKlCK7x/5nvZl2jNlxxQGLhtIwlW99KVZb0umo9xP/JYsFs8GOgNngDigk1Lq23seSdOsxSsAunwLiWfgx96QnnZfl8nnmo/3G77POxHvsC1uG0/8/gRbzmyxcrCaM/Hw8CAhIcFpk4FSioSEBDw8PO7pexbVGspsO1mELFNCSqnj9xqkNejFYu2G7T/A/wZCvUHQ5sO7n38HB84dYPjq4cRcjmFIzSH0qdwHF7FlcV4tN0pNTSUmJobk5GRHh3LfPDw8CAoKwmw23/T+nRaL7zrXLyKDgTcx7gjSMRrOKKBajiPWtJyo3g1id8A/k6BoNeP1fapQqAJzH5nLG+vf4PMtn7PtzDbea/geBdytVxhPy/3MZjOlSpVydBh2Z8mPPEOBCkqpykqpakqpqkopp0sCcZeSefWXXSSn6kXBB8pD70JII/h9KJzcmqNLebl5MbbJWEbXGc3fp/6my8Iu7D6rq6JrDz5LEsEJ4KKtA7G1rcfP8/2m4zw/Zyup6RmODkezFpMZnpgJXkVgXg9IjMvR5USE7pW6M6vNLDJUBr0W9eKH/T847ZyxplnCkkQQBawSkVdE5KXrh60Ds7Y2VYrxXscqLN8fx7B520nP0P9hPzA8C0PX7+DKOZjfG9Jz/tRwNf9qzG83n3rF6vHBxg94ec3LJKU6vB2HptmEJYngOLAUcAO8sxxOp3vdYP7v4Yos3BnL//28iwydDB4cxcLg0YlwfD0sfsUql/T18OWLFl8wtOZQlkQvoevCrhw8f9Aq19a03OSui8VKqbch93Qoy6kBjcuQeC2dCcsPkd/dxBvtQv/T9EFzUtWegNM7YP1EKFYNavbK8SVdxIX+VfsT5h/GqDWj6P5Hd16t9yody3a0QsCaljtY0qGs/oPWoWxYy3I83aAUM9Yd4/Ol+ie8B0qLt6B0M/hjOJzYbLXL1i5amx/b/0iYfxivr3udN9a9QXKa824x1LSsLJkaGofRoSwBjA5lgEUdykSkjYgcEJHDIjI6m8+bishFEdmeebxxL8HfLxHh9XaV6Fq7BBNWHGbq6iP2GFazB5MrdP4GfIobi8eXT1vt0oXzFWZqq6kMqDaAXw7/Qvc/u3Ps4jGrXV/THMVmHcoyH0L7EmgLhALdRCS7XoRrlVLVM493LInHGkSE9x+rSvuw4ny4aD/f/RNtr6E1W8tfCLp+D9cuwbyekGa97mQmFxODawxmcsvJxF2Jo+sfXVlybInVrq9pjmDR9tGsHcpEZASWdSirAxxWSkUppVKAuUDOGmtamclF+OzJMFpWCuD1X3fzy7YYR4ekWUuRytBxEsRsgkWjrH75hoEN+bH9j5T1LcuI1SP4cOOHpFpht5KmOYItO5QFYjyDcF1M5nu3qi8iO0RkkYhUzu5CIjJARCJFJDI+/v6ak9yO2eTCF0/VpH5pP0b8uJPFu603laA5WOXHoOFLsGUmRH5j9csX9SzKjNYz6Bnak+/3f0/vxb05k3TG6uNomq3ZskNZdltxbt2vuRUIVkqFAROB/90mhmlKqXClVLi/v78FQ98bD7OJr3qFUy2oAEN+2Maag9ZNNpoDNX8NyraCP0dB9AarX95sMjOq9ijGNR1H1MUoei/uzYnLt86kalruZsmuoQnZHO+KyN2meWKArJ3Cg4BTWU9QSl1SSiVm/v5PwCwihe/xz2AVnu6uzOxTh7IBXgz4NpLNx845IgzN2lxM8Ph08C0J83vBxZM2GaZFcAu+fuhrElMT6b2oN4fPH7bJOJpmC5ZMDXlgTAcdyjyqAYWAfiIy7g7f2wyUE5FSIuIGdAV+y3qCiBSVzE38IlInMx6HFYUvkN/M7H51KO6bj6dnbGZXjNNX1tAA8vkai8epV4ydRKm22fZZuXBlZraeCUCfJX10nSLNaViSCMoCzZVSE5VSE4GWQCXgMeCh231JKZUGvAAswVhcnq+U2iMiA0VkYOZpnYHdIrIDmAB0VQ4u6lLYy505/etSIL+ZXt9s5OCZy44MR7OWgIrw2FQ4tRX+eAls9H+zsgXLMqvtLLzMXvRb0o/Np633LIOm2YolzesPAHWUUhczXxcANiqlKorINqWUXfsX26sfQXRCEk9MMeaUfxxYn2A/T5uPqdnByg9g9cfQ9hOoO8Bmw5xJOsOzS58lJjGGz5p+RuMgix690TSbyWnz+jHAdhGZISIzgW3ApyLiCSyzXpi5S7CfJ3P61yU1PYOnvtrIqQtXHR2SZg1NRkP5trB4NBxda7NhingWYUabGZTxLcPQFUNZfHSxzcbStJy6YyLInL//C4jA2NHzP6ChUmq6UipJKTXSDjE6TLki3nzbry6XrqbSY/pGziZa78EkzUFcXKDTNPArY7S5vGC7HT4FPQry9UNfExZg1ClacHCBzcbStJy4YyLInK//n1IqVin1q1Lqf0qpU3f6zoOmSmABZvStTezFZHp+vYmLV/RDQ07Pw8dYPE5PhblPQcoVmw3l5ebF5JaTaRDYgLc3vM2sPbNsNpam3S9Lpob+EZHaNo8kFwsPKcS0XrU4EpdI7xmbSLx2f83StVykcDljW+npXUZ3MxvuUcjnmo8JzSbwUPBDfBr5KRO3TdSNbrRcxZJE0AwjGRwRkZ0isktEdto6sNymUTl/Jj5Vg10nL/LMrEjd8vJBUL41NH8Vds03+h7bkNlkZkzjMXQq14lpO6fx0aaPyFC6U56WO9y1HwFG0TgNaF25KGOfCGPY/O0MmrOVqT1rYTZZVLdPy60ajYDYHfDXaxAQCmWa2Wwok4uJt+q/hZfZi9l7Z5OUmsRbEW/h6mLJf4aaZjuWlJiIxnhCuHnm769Y8r0HVccagbzXsQordMvLB4MIdJwMhSvAgr5w/piNhxNGhI9gUPVB/HrkV0atGUVKeopNx9S0u7GkxMSbwMvA9f5/ZuA7WwaV22VtefnKzzt1y0tn5+4NXeeAyoC53SHFto34RITnwp7j5dovszR6KYNXDOZKqu0WrDXtbiz5yf4x4FEgCSBz15BT9iy2pgGNyzCkRTnmR8bwzsK9evHP2fmVMRraxO2FX5+36eLxdT1Ce/BOxDv8E/sPA5cN5FLKJZuPqWnZsSQRpGRuI1Vg9C62bUjO43rLy5nrj/GZbnnp/Mq2hBZvwp5fYN2dymhZz2PlHmNM4zHsOruLfkv6kXDVYaW2tDzMkkQwX0SmAr4i8gzG08Rf2TYs55C15eXEFYeZolteOr8GQ6FyJ1j2Nhyyz4PzrUNaM7H5RI5dPEafxX04naR7Ymj2Zcli8afAAuAnoALwRmbxOY2bW15+tGg/3+qWl85NBDp8YXQ4++lpSLBPcm8Y2JApraZw9upZei/qzfFLx+0yrqaBZYvFw4B9SqmRSqkRSqmldojLqdzU8vJ/u/l5q2556dTcPI3FY3ExFo+v2acCba0itfi69ddcTbtKr0W9OHheTzdq9mHJ1JAPsERE1orI8yJSxNZBOaPrLS8jyvgxcoFueen0CobAEzPh7AH4ZSBk2Ofhr1C/UGa2mYlJTPRd3Jed8Xnu2U3NASyZGnpbKVUZo09xcWC1iDywVUdz4nrLy7CgAgz+YSurdctL51a6KTz0HuxfCGvH2m9Y39LMajsLHzcf+v/Vn02xm+w2tpY33cuDYXHAaYwOYgG2Ccf5ebq7MqNvHcoFePPst5FsOqpbXjq1eoOgWhdY+T4csF8p6SDvIGa3nU2gVyDPLXuOVSdW2W1sLe+xZI3gORFZBSwHCgPPKKWq2TowZ1Ygn9HyMtA3H0/P3MzOmAuODkm7XyLQfjwUqwY/PwNnD9ltaP/8/sxoPYPyBcvz4soX+SPqD7uNreUtltwRBAMvKqUqK6XeVErttXVQD4LCXu58178uvvnN9PpmEwdO65aXTsucD7rMAZMb/NANku3Xy9rXw5fpradTI6AGr6x9hfkH5tttbC3vuG0iEBGfzN+OAY6LSKGsh33Cc27FCuRjTv+6uJlc6PH1Ro6dtW3pAs2GfEvAk7PgXBT8/KzdFo8BPM2eTG45mUZBjXj3n3f5etfXdhtbyxvudEfwfeavW4DIzF+3ZHmtWeB6y8u09Ay6T9ctL51aSENo8xEcXASrP7Lr0B6uHoxrNo62IW0Zt3Uc47eO12VNNKu5bSJQSrXL/LWUUqp05q/Xj9KWXFxE2ojIARE5LCKj73BebRFJF5HO9/5HyP1ubXkZf1m3vHRadZ6B6j1g9cew73e7Dm12MfNhow95vNzjTN81nQ82fqB7GmhWcaepoZp3Ou52YRExAV9i9DMIBbqJSOhtzvsYWHL/f4zc7+aWlxt1y0tnJQKPjIXAWsbzBXH77Tq8ycXEm/XfpE/lPsw9MJfX/n6NtAzdMU/LmTtNDY3NPL4ENgLTMGoMbQQmWHDtOsBhpVSUUioFmAt0yOa8wRjlK+LuIW6ndL3lZVR8km556czMHtDlOzDnh7nd4Kp9d4WJCC/VeonBNQbze9TvDF81nGvp+i5Tu393mhpqppRqBkQDNZVS4UqpWkAN4LAF1w4ETmR5HZP53g0iEohR5nrKnS4kIgNEJFJEIuPjnfshrawtL/vP2qxbXjorn+LQ5Vu4cAJ+6g8Z9v3fUUQYUG0Ao+uMZsWJFTy//Hnd00C7b5ZsH62olNp1/YVSajdQ3YLvSTbv3bq6NQ54WSl1x/+KlFLTMhNRuL+/vwVD527XW15uPHqOQXO2kpKm53mdUs4PWLIAACAASURBVMl68PAYOLwUVrznkBC6V+rOew3eY/PpzQxYOoCL1+y3tVV7cFiSCPaJyHQRaSoiTUTkK2CfBd+LwWhxeV0QcOqWc8KBuSJyDOgMTBKRjhZc2+nd1PJy/nadDJxV+NNQqw/8/Rns/tkhIXQo24GxTcayJ2EP/Zb04+zVsw6JQ3NeliSCvsAeYCjwIrA387272QyUE5FSIuIGdAV+y3pC5g6kEKVUCEap60FKqf/dQ/xO7XrLyz92xvLoF3+z+6T+ac4ptR0DJeoanc1ObXdICC2DW/Jl8y+JvhRN38V9iU2MdUgcmnOypOhcslLqc6XUY5nH50qpZAu+lwa8gLEbaB8wXym1R0QGisjAnIf+YBjQuAzTe4VzLimFDl+u49MlB7iWptcNnIqrOzw5G/L7wXedIP6AQ8KICIxg2kPTSLiaQK/FvTh28ZhD4tCcjzjbQynh4eEqMvLBe57t4pVU3lm4l5+2xlChiDefPFGNakG+jg5LuxcJR+CbNuBigr6LoFAph4SxL2EfA5cZP2tNazWNCoUqOCQOLXcRkS1KqfDsPruX6qOaDRXIb2bsk2F80yecC1dTeGzSej5Zsl/fHTgTvzLQ61dIS4bZHeDSrUti9lHJrxIz2szA7GKm75K+bI9zzHSV5jx0Ishlmlcswl/DmtCpRiBfrjxCuwl/s+OErl7qNIqEQo+f4Mo5IxkkOWbhtnSB0sxuO5uC7gUZsHQAkacfvLtozXosKUNdXkS+EpG/RGTF9cMeweVVBfKZ+eSJMGb0rc3l5DQem7SOjxfv188cOIvAWvDUPOMZg28fs/sDZ9cV9yrOrLazKOpZlCErhnDovP1KaGvO5a5rBCKyA+OBry3AjX+JlFJbbBta9h7UNYLbuZScynsL9zI/MoayAV580rkaNUoWdHRYmiUOLYMfukJgTej5i9EL2QFOJZ6ix589EBHmPDyHop5FHRKH5lg5XSNIU0pNVkptUkptuX5YOUbtNnw8zIzpHMbMvrVJupbG45PX8+GiffruwBmUawmdv4aYzUYfg9S7brazieJexZnccjJXUq/w3LLn9ENn2n9Ykgh+F5FBIlJM9yNwnKYVAlgyrDFPhpdg6uooHpmwlq3Hzzs6LO1uQjtAh0lwdDUs6Avpjik2WKFQBcY1G8exS8cYunKork2k3cSSqaGj2bytLC1FbW15bWooO2sOxjP6p52cvpRM/0alealVeTzMJkeHpd3Jpq/gzxFQpTN0mmZsMXWARUcXMWrNKFoFt+KTxp9gclAcmv3laGrolj4E99SPQLONxuX9WTKsMV1ql2TamigeHr+WLdHnHB2Wdid1noGWb8HuBbBwGDjo+Z22pdoyInwES6OXMmbzGN3cRgMs2zVkFpEhIrIg83hBRMz2CE67PW8PMx92qsp3/epyLS2DzlM28N7CvVxN0WsHuVbDYdBoOGydBX+95rBk0Ltyb3qG9uT7/d8zY88Mh8Sg5S6uFpwzGTADkzJf98x8r7+tgtIs17BcYZYMa8yHf+5j+t9HWb4/jk86VyM8RC/j5ErNX4eUJNjwBbh5QbNXHBLGiPARnL1yls+3fI5/Pn/al2nvkDi03MGSRFBbKRWW5fWKzC2lWi7h5e7K+49V5eGqxRi1YCdPTN1A34hSjGxdgXxueg44VxGB1h/CtUSj77G7F0QMtnsYLuLCew3fIyE5gTfWvYFfPj8iikfYPQ4td7Bk11C6iJS5/kJESpPleQIt92hQ1rg76FE3mG/WHaXt+DVsOqrXDnIdFxd4dAKEdjSmiCIdMz3jZnJjXLNxlPYtzbCVw9ibsNchcWiOZ0kiGAmsFJFVIrIaWAEMt21Y2v3ycnfl3Y5V+P6ZuqQrRZdpG3jrtz1cSdFtMXMVFxN0+grKPWQsHu/80SFheLt5M7nlZAq4F2DQskGcuHzi7l/SHjgWVR8VEXegAkbXsf1KKYdtQtbbRy2XdC2NMYv3M2tDNMF++RnzeDXqlvZzdFhaVqlXYc4TEL3eaH1Z8RGHhBF1IYqei3pS0KMg37b9loIe+un1B401qo/WAqoAYUAXEellreA02/F0d+XtDlX44Zl6KAVdpv2j7w5yG3M+6PYDFK8BP/aBI44p41XatzRftPiC00mneWH5C7r/cR5jyfbRb4FPgYZA7cwj26yi5U71y/ix+MVG9IkIYeb6Y7QZt5YNRxIcHZZ2nbs3dP8RCpeHud3h+D8OCaNGQA0+bvwxuxN2M2rNKNIy9A8MeYUlTxbvA0JVLnnyRE8N5czGqARG/bST6IQr9KofzMttKuLpbsnmMc3mEuNgRlvj196/Q/HqDglj3v55vLfxPR4v9zhv1n8TEXFIHJp15XRqaDegyxU+IOqW9mPR0Eb0bRDCt/9E03rcGtYf0c3OcwWvAKOxjUcBo3x13H6HhNGlYheeqfoMPx36iSk7pjgkBs2+LEkEhYG9IrJERH67ftg6MM128ru58mb7ysx/tj6uLsJTX23ktf/tIvGangpwuAJBRjIwmY3GNueyK/Vle4NrDKZDmQ5M2jGJBQcXOCQGzX4smRpqkt37SqnVNonoLvTUkHVdTUnn078O8M26oxQvkI8xnavRoGxhR4elndkLMx821g/6LoYCgXYPITUjlSErhrD+1HrGNxtP0xJN7R6DZj05LTq3OrvDwoHbiMgBETksIqOz+byDiOwUke0iEikiDS25rmY9+dxMvN4ulB+frY+bqwvdp2/k/37ZxeVkx5RL1jIVCYUeP8OV8/BtR0iMt3sIZhczY5uMpVKhSoxcPZId8bqgwIPKkl1DnUTkkIhcFJFLInJZRC5Z8D0T8CXQFggFuolI6C2nLQfClFLVgaeB6ff+R9CsITykEIuGNuKZRqX4YdNx2oxby9pD9v/HR8sisCZ0n2+0vPzOMS0v85vz82WLL/HP788Ly1/g6EXHTFVptmXJGsEY4FGlVAGllI9Sylsp5WPB9+oAh5VSUUqpFGAu0CHrCUqpxCy7kTyBXLEzKa/yMJt49ZFQFgysj7vZhZ5fb6Ln1xvZcCRBlyt2lOAI6PqdsXA8p7NRo8jO/PL5MbXlVFzEheeWPcfZq3pzwYPGkkRwRim17z6uHQhkfV49JvO9m4jIYyKyH/gD467gP0RkQObUUWR8vP4p1dZqBRfizyGNGNWmAvtiL9Htq3/oNHk9S/eeISNDJwS7K9sSOn8DJ7fCXMe0vCzhU4JJLSZxLvkcg5YNIjHF/glJsx1LFovHY2wf/R9wo7SEUurnu3zvCaC1Uqp/5uueQB2lVLalFkWkMfCGUqrlna6rF4vtKzk1nR8jTzB1TRQx569SvogXzzUtQ7tqxTGbLH0wXbOKHXPhl2ehfFujHIXJ/m1B1sasZfCKwdQuWptJLSZhdkAM2v3J6XMEPsAV4CGgfebRzoLvxQAlsrwOAk7d7mSl1BqgjIjoLSu5iIfZRM/6Iawa0ZRxXaojCMPm7aDZp6uYveEYyam6EK3dhHWFR8bCwUVGQsiw/999o6BGvBXxFv/E/sPr618nQ2XYPQbN+u76SKlSqu99XnszUE5ESgEnga7AU1lPEJGywBGllBKRmoAboGsf5EKuJhc61gjk0bDirNgfx6RVh3nj1z1MWH6Ivg1K0bN+MD4e+qdDm6vd31gnWPYmuHlC+wlGjwM76li2I3FX4pi4bSIB+QN4qdZLdh1fs767JgIR8QD6AZUBj+vvK6Wync/P8nmaiLwALAFMwDdKqT0iMjDz8ynA40AvEUkFrgJdckspCy17Li5Cy9AitKgUwMaj55i06gifLDnAlFVH6FE/mKcblMLf293RYT7YGr4IKYmw5hOjy1nrD+yeDJ6p+gxxV+KYsXsGRfIXoXul7nYdX7MuS9YIfgT2Y/w0/w7QHdinlBpq+/D+S68R5D67T15k8qoj/Lk7FrPJhSfDg3i2cRlKFMrv6NAeXErB4tGwcQo0eRma/Z/dQ0jPSOelVS+x8sRKPmnyCa1DWts9Bs1yd1ojsCQRbFNK1RCRnUqpapmN65copZrbIti70Ykg9zp6Nompq4/w09YYMhS0r1aM55qWpUJRb0eH9mDKyIDfB8O276DVu9BgiN1DSE5LZsDSAew+u5upraZSu2htu8egWSani8XXHzG9ICJVgAJAiJVi0x4gpQp78tHj1Vg7qjl9I0L4a+8ZWo9bQ/9Zm9kSfd7R4T14XFyMNYLKj8HS1yHyG7uH4OHqwcTmEynhXYKhK4Zy6Pwhu8eg5ZwldwT9gZ+AasAMwAtjm6dDyhLqOwLncT4phVkbjjFz/TEuXEmlTqlCDGpahibl/XVpY2tKS4F5PeDQX/DYVAjrYvcQYhNj6fFnDxCY8/AcinrqgsW5TY6mhnIbnQicz5WUNH7YdIKv1kRx+lIylYv78FzTMrStUgyTi04IVpG15eWTs6BSe7uHcODcAfos7kNRz6LMbDOTAu4F7B6Ddns5XSNwx9jdE0KWXUZKqXesGKPFdCJwXilpGfxv20mmrD5C1NkkShX25NnGpXmsZiDuriZHh+f8rl2G2R3h9E7oNhfKtrB7CJtiN/HssmcJ8w9jaqupuJv0DrLcIqdrBL9i1AhKA5KyHJp2T9xcXXiydgmWvtSESd1r4uluYvTPu2g8ZiVfrYkiSfdDyBl3b+ix4N+Wl9Eb7B5CnWJ1+KDhB2w5s4VX1r5CugMeetPunSV3BLuVUlXsFM9d6TuCB4dSir8Pn2XSyiNsiEqgQD4zvSNC6BsRQkFPN0eH57xuann5GxSvYfcQZu+ZzSeRn9CtYjdeqfOKXhPKBXJ6R7BeRKpaOSZNQ0RoVM6fHwbU4+dBEdQpVYgJyw8R8dEK3vl9L6cuXHV0iM7pRstLX/i2E8TdT83InOlVuRe9Q3vzw/4f+Ga3/XczaffmtncEIrILoyy0K1AOiMIoOieAUkpVs1eQWek7ggfbwTOXmbLqCL/uOIWLQMfqgQxsWoYy/l6ODs35nIuCb9oav++zEAqXs+vwGSqD0WtHs+joIj5o+AHty9h/AVv7130tFotI8J0uqpSKtkJs90wngrzhxLkrTF8bxdzNJ0hJz6BN5aIMalqWqkF6J8o9idsHM9uByoCn5kGJOnYdPiU9hUHLBrHlzBa+bPElEYERdh1f+9f9JoL8QKpSKjXzdQXgYSD6biWobUkngrzlbOI1Zqw7yuwN0VxOTqNRucI817QM9Uv76XlnSyUcge8eh8uxRl+Dio/YdfjLKZfps7gPMZdjmNFmBqF+tzYq1OzhfhPBGqCfUupQZpXQTcAcjLaTm5RSr9gq4DvRiSBvupScypx/jvP130c5m3iNOqUKMbxVeeqW9nN0aM4hMR5+6AKntkHbMVDnGbsOH3cljp5/9iQ5PZnvHv6OEt4l7v4lzaruNxHsUkpVzfz9u0AhpdTzIuIGbLn+mb3pRJC3JaemM2/zCb5YeZj4y9doVK4wwx+qQPUSvo4OLfdLSYIFT8PBxdDgRWjxplGmwk6iLkbRa1EvfN19md12NoU8CtltbO3+dw1lzRDNgaUAmf2HdTcKzSE8zCZ6R4SwZmQzXn24EntOXaLjl+voP2sze05ddHR4uZubJ3SZA7X6wrpxRnObtBS7DV+6QGm+aP4Fp5NO88LyF7iSesVuY2t3dqdEsFNEPhWRYUBZ4C8AEdE/emkOl8/NxDONS7NmVDNGPFSejUfP8ciEv3l+zlYOx112dHi5l8kV2n0OzV+HXfNhzuOQbL8EWj2gOmMaj2FPwh5GrRlFWoZ+iDA3uFMieAY4i1Fa4iGl1PX0HQp8auO4NM0iXu6uvNC8HH+Pas7g5mVZdSCOhz5fw0vzthOdoB+Az5YINB4BHacYtYlmPAyXbttF1uqal2zOq3VfZXXMaoasGMLlFJ24He2eis6JSE2l1FYbxnNXeo1Au5OExGtMXRPFrPXHSMtQPBkexAvNyxHom8/RoeVOR1bAvF7g4QPdF0AR++3ombd/Hh9t+ogg7yDGNx9P6QKl7TZ2XmS16qMislUpVdNqkd0HnQg0S8RdSubLlYf5ftNxBOGpuiUZ1LQMAT4ed/9yXhO706hcmnoVus6BUo3sNnTk6UiGrx5OSnoKHzf+mMZBje02dl5jzUSwTSll/8IlWehEoN2LmPNX+GLFYX7cEoPZJPSuH8KzTcpQSNcyutmF4/BdZzh/FDpOhqqd7TZ0bGIsQ1cOZf+5/QyuMZj+VfvrZ0RsIKe1hrJ6+x4HbiMiB0TksIiMzubz7iKyM/NYLyJh9xiPpt1RUMH8fPR4NZa/1ISHqxRj2tooGn28gs/+OsDFq6l3v0Be4VsS+i2BwHD4qR+sm2D0RbaDYl7FmNV2Fm1KtWHCtgmMWD1C7yiyM4vuCEQkEAjm5n4Ea+7yHRNwEGgFxACbgW5Kqb1ZzokA9imlzotIW+AtpVTdO11X3xFoOXHozGXGLTvEH7ti8fFw5dkmZegTEYKnu+vdv5wXpCYb20r3/g/qDoTWH4CLfXpFKKWYuWcm47aOo6xvWcY3G0+Qd5Bdxs4LctqY5mOgC7AXuF5cXCmlHr3L9+pj/MPeOvP1K5lf/PA25xcEdiulAu90XZ0INGvYc+oiny89yLJ9cRTydGNQ0zL0qBeMh1k3yCEjA/56Df750uh01ukrMNtvsX3dyXWMXDMSk5j4tMmn1C12x58NNQvldGqoI1BBKfWwUqp95nHHJJApEDiR5XVM5nu30w9YZMF1NS3HKhcvwPTetfllUASVi/vw3h/7aDxmJbM3HONaWh5vpuLiAm0+MO4G9i00up5dOWe34RsENuCHR37Az8OPZ5c+y5x9c3C2lrrOxpJEEAWY7+Pa2a32ZPu/pog0w0gEL9/m8wEiEikikfHx8fcRiqZlr0bJgnzbry5zB9QjxM+TN37dQ/NPVzNv83FS0/P4A/T1n4cnZhj1ib5+CM4fs9vQwT7BzHlkDo2DGvPRpo94fd3rXEu/Zrfx8xpLpoZ+AsKA5Rj9CABQSg25y/csmhoSkWrAL0BbpdTBuwWsp4Y0W1FKsfbQWcb+dYAdMRcJ8cvPiy3L0z6sOCaXPLyLJXo9/NAVTO7Q/UcoXt1uQ2eoDCbvmMyUHVOoWrgq45qNIyB/gN3Gf5DkdI2gd3bvK6Vm3eV7rhiLxS2AkxiLxU8ppfZkOacksALopZRaf8dAMulEoNmaUopl++IY+9cB9p++TLkAL15qVZ7WlYviklcTQtx+mNPZmCJ6cjaUa2nX4ZdFL+P//v4/PM2efN70c6oH2C8ZPShy/ByBiOQDSiqlDtzjwA8D4wAT8I1S6n0RGQiglJoiItOBx4HrTW7SbhfodToRaPaSkaH4c3csny89yJH4JCoX92H4Q+VpViEgb+5zvxRrPHgWtxcenQA1eth1+EPnDzFkxRDOXDnDa/Veo1O5TnYd39nl9I6gPUZtITelVCkRqQ68Y+GCsdXpRKDZW3qG4tftJxm37BDHz12hRklfRjxUgYgyebA5TvIlmN8LolZC0/+DJqOM2kV2cvHaRUauHsmG2A10q9iNkbVHYna5nyXMvCeniWALRhnqVdefKs7aq8DedCLQHCU1PYMFW2KYsPwQsReTqVe6EMMfqkDtkDxWVz89FX4bDDt+gJq94JHPjaqmdpKWkca4LeOYtXcW4UXCGdt0rO5tYIGcbh9NU0rdWqdW7+XS8hyzyYVudUqyckRT3mofyuG4JJ6YsoHe32xiZ8wFR4dnPyazUYai8UjYOhvmdoNriXYb3tXFlRG1R/BBww/YdXYXXRd2Zf+5/XYb/0FkSSLYLSJPASYRKSciEwGLFnY17UHkYTbRp0Ep1o5qxittK7Iz5gKPfrGOZ2ZHsvfUJUeHZx8i0Pw1o7fB4WUwqx0kxtk1hPZl2jOr7SwyVAY9/+zJ4qOL7Tr+g8SSqaH8wKvAQxjPBiwB3lVKJds+vP/SU0NabnM5OZUZ647x1dooLien8XDVogxtUZ4KRb0dHZp9HFgEP/YFrwDo8TMULmvX4c9ePcvwVcPZGreVp6s8zZAaQzDZqSyGM7Fm9VET4KmUctiPPToRaLnVxSupfP13FN+sO0ZSShrtqhVnaIuylA3IAwkhZgt8/ySoDHhqHpSoY9fhU9NT+WjTR8w/OJ+GgQ35uPHH+Lj52DWG3C6ni8XfAwMx6gxtAQoAnymlPrF2oJbQiUDL7c4npfDV2ihmrj/G1dR0OoQVZ0iLcpT293J0aLaVcMR41uDSKXj8a6jUzu4hzD8wnw83fkigdyATmk2gtK9udnNdThPBdqVUdRHpDtTCKAOxRSlVzfqh3p1OBJqzSEi8xrS1UcxeH821tHQeqxHEkBZlCfbzdHRotpN01rgzOLUN2o6BOs/YPYStZ7YybNUwrqVf46NGH9G0RFO7x5Ab5XTXkFlEzBjF535VSqWidw1p2l35ebnzSttKrBnVjKcblGLhzlM0H7ualxfs5MS5B7Tevmdh6L0QyrWGP0fA0jeNaqZ2VLNITea1m0ewTzBDVgxh6o6pumjdXViSCKYCxwBPYI2IBAN5ZGuEpuWcv7c7r7ULZe2oZvSsF8wv20/S7NNVvPLzLk5euOro8KzPLT90+Q7Cn4Z144z+Bmkpdg2hqGdRZrWZxSOlH+GL7V8wfPVw3ezmDu5psfjGl0RclVJpNojnrvTUkObsTl9MZtKqw8zddAKFomvtkgxqVoZiBexX898ulIK/P4Pl70CpxkZy8Chg5xAUs/fO5rMtn1HGtwzjm42nhHcJu8aQW+R0jcAdox5QCDd3KHvHijFaTCcC7UFx8sJVvlx5mB8jTyAiPFWnJIOaliHAx8PRoVnXjrnw6/NQuIJRvbTAHXtP2cT6U+sZuXokIsKnTT6lXrF6do/B0XKaCBYDFzF2DN3o2KGUGmvNIC2lE4H2oDlx7oqRELbE4Ooi9KgXzMAmZfD3dnd0aNZzZCXM6wkePtB9ARQJtXsIJy6dYMjKIRy9eJTh4cPpUalHnqoVldNEsFspVcUmkd0HnQi0B1V0QhITVxzm560xuLm60Lt+CAMal8bP6wFJCLE7jeqlqVeh63fGdJGdJaUm8erfr7L8+HIeLfMob9R/A3fTA/L3exc5TQTTgIlKqV22CO5e6USgPeii4hOZuOIwv24/iYfZRO+IEAY0Kk1BTzdHh5ZzF04YzxqcizLqFVXtbPcQMlQG03ZO48vtX1LFrwqfN/ucop5F7R6HveU0EewFygJHMTqUCUbzev0cgabZ0OG4RCYsP8TvO0/h6eZK3wYh9G9YmgL5nbzs8tXzMLc7RK8zdhY1eRm87f8P8YrjK3hl7Svkc83HuGbjHvhmNzlNBMHZva+Uis7ufVvTiUDLaw6eucz4ZYf4Y1cs3u6u9GtUiqcblsLHw4kTQmoyLHsTNk8HFzPUew4aDIV8vnYN48iFIwxZMYRTSad4te6rdC5v/zsUe7FGh7KGQDml1AwR8Qe8lFJHrRynRXQi0PKqfbGXGL/sEIv3nMbHw5VnGpWmT4MQvJ05IZyLghXvw+4F4OELjYYbTyOb7beV9uK1i7y85mXWnVpHlwpdeLn2y5hNTvx3ehs5vSN4EwgHKiilyotIceBHpVQD64d6dzoRaHnd7pMXGbfsEMv2ncE3v5kBjUvTu34Inu72aw5jdbE7YfnbRklrn0BoOhrCnrJbw5v0jHTGbxvPjN0zqBlQk8+afoZfPj+7jG0vOa41BNQAtmbpULZTrxFommPtjLnA50sPsvJAPIU83Xi2cWl61g8mv5sTJ4Sja2HZW3AyEgqXh+avQ6X2dmuH+WfUn7y5/k0K5yvMlFZTCPbJdmbcKeW01lCKMrKFyrzYA1wxS9OcR7UgX2b0rcPPgyKoXNyHDxftp/GYlUxfG0VyavrdL5AblWoE/ZcZTyEDzO8J01saCcIOHi79MF+3/poraVfo8WcPdsTvsMu4jmZJIpgvIlMBXxF5BlgOTLfk4iLSRkQOiMhhERmdzecVRWSDiFwTkRH3FrqmaQA1Sxbk2351WTCwPhWKevPeH/toNGYlM9cddc6EIGLcBTy3AR79Ai7HGh3Qvu0Esbb/h7mafzW+bfstPm4+9FvSjxXHV9h8TEezdLG4FUaHMoAlSqllFnzHBBwEWgExwGagm1Jqb5ZzAoBgjMqm55VSn97tunpqSNPu7J+oBD5bepBNR89R1MeD55uV4YnwEniYnbRrV+pV2PQVrB0LyRegyuNGm8xCtu01cC75HIOXD2bX2V28UvcVulXsZtPxbO2+1ghE5DL/lpu+dYIuGTgCvKqUWn6b79cH3lJKtc58/QqAUurDbM59C0jUiUDTrEMpxYYjRkKIjD5PIU83nqpTkh71gilawElrGV29AOsnwIZJkJEKtfpA41HgXcR2Q6Zd5eU1L7PyxEr6VunLizVfxEUsmUjJfazWqjLLBU1AFWDO7cpPiEhnoI1Sqn/m655AXaXUC9mc+xZ3SAQiMgAYAFCyZMla0dEOeYRB05yOUooNUQnMWHeMZfvOYBKhbdVi9IkIoWZJX+estXP5NKz+GLbMAld3qDcIGgyxWWXT9Ix0Ptz0IfMOzKNtSFvea/gebibne8rb6okgy4WfVUpNvc1nTwCtb0kEdZRSg7M59y30HYGm2dTxhCvM3nCMeZEnuJycRlhQAfo0COGRqsVxc3XCn3ITjsCK92DPz5CvkPEMQu3+YLb+HY9Sihl7ZvD5ls8JLxLO+Objna4nss0SwV0G1VNDmpYLJV1L4+etMcxYf4yo+CT8vd3pXrckT9UtSYC3E04bndpm9Dw4sgJ8gqDZK1Ctq02eQfgj6g9eW/cawd7BTG45mWJexaw+hq04KhG4YiwWtwBOYiwWP6WU2pPNuW+hE4Gm2VVGhmLt4bPMWHeUVQfiMZuE9tWK06dBCNWC7FvqwSqiVhvPIJzaavQ+aPEGVHzE6s8gqyQDhgAAD8RJREFUbIrdxIsrX8TD1YNJLSdRsVBFq17fVhySCDIHfhgYB5iAb5RS74vIQACl1BQRKQpEAj5ABpAIhCqlbtsKUycCTbO+qPhEZm+I5sfIEySlpFMruCB9IkJoU6UoZpMTTRspBft+M+4QEg5DUG1o+RaENLTqMIfOH2LQ8kFcTrnMZ00/I6J4hFWvbwsOSwS2oBOBptnOpeRUFkTGMGvDMaITrlDUx4Oe9YPpVqckhZypDHZ6Gmz/DlZ9ZDyHULYVtHwTila12hBnks4waPkgoi5E8VbEW3Qo28Fq17YFnQg0Tbsn6RmKVQfimLHuGH8fPoubqwsdqxenT0QpQos70SJp6lXYONXonZx8Cao+Ac3+DwqVssrlE1MSGbZqGP/E/sP/t3fv0VHWdx7H398kEMgFkSQghEuIIhbrBRYhAUFWwENbT7VVT2FXXDg9SwVEuyv1WLfrqj21ni3alcqGQ0FUili1aGnlgCteQEi4yqWA5RIghIvkQiAQQm7f/eP3TBhCBJLM5Bmc7+ucHGaeeWbmO5zk+czze57n+3vk1keYdPOkiD0Ty4LAGNNsu78q57U1+1m86RBnqmsZ3LsTE4dmMOpbXYi7UoaNzhyH1S9D3myoq4GBE2H4zyCpc4tfurq2mmdyn2HJ3iXc1+c+fpH1C+JiIq/fkwWBMabFTlRU88cNBby+5gCHys6Q3rE9D2X3YuxtPa+cyXJOHoHPXoBNCyCuHWRPhSHT3FzKLaCqvLL5FeZsncOw9GHMuGMGCW0SQlR0aFgQGGNCpqa2jo92HuO1NfvIyy+lfZtYfjAgnQlDMri+S7Lf5V2e4j3w8S9hx/vuGoTh02Hgj1t8DcI7u97hV3m/om+nvswaOYvU9qkhKrjlLAiMMWGx4/BJXl+zn/c2H6Kqpo7br0tlwpAM7ryhMzExkTlWfp5Dm9wpp/s+g6t6wIifwy1jIab5fZlWFq5k+mfT6dSuEzmjcuh9VWiOR7SUBYExJqxKT1exaF0BC3IPcPRkJb1SEngoO4MHBna/MqbU3PuJC4QjmyHtBnf8oN890MyZyrYXb2fKiinUai2/u/N39O/cP7T1NoMFgTGmVVTX1rF8+1Hmr97PxgPHSWwby/3/0J1/GZJBZlqS3+VdXF0d7PyzmzqzZLebKW3Qv7rmdu2vbvLLHSw/yJSPpnD41GFeGP4Co3uNDn3NTWBBYIxpdVsLy3ht9X7+svUw1bXKiL5pTBiSwfA+aZE9bFRXB7uXQ97/wr6V0CYBbv0nGDwZUq9r0kuVVZYx7eNpbCnawhO3PcGD/R4MU9GXZkFgjPHNsfJKFq09yB/WHqCo/CyZaYlMGJLBfQO6R/48y0e3QV4ObHsHaqvg+jGQNRl633HZrSsqayp5ctWTrChYwfh+45k+cLovrawtCIwxvquqqWPptiPMX72PLYUnSI6P446+aWRlppCVmcK1aYkRezEW5V/Bhnmwfh5UFEOXb7tAuOkB1wr7EmrravnNht+wcOdC7up1F88Pe5742Es/L5QsCIwxEUNV+eJgGW+uLeDz3cUcPVkJQGpSPFmZncjKTCH72hQyUyMwGKorYdvbbnKcop2QmOZaXw/8MSSlXfSpqsobO95gxoYZDOg8gJl3zuSq+PDModAYCwJjTERSVQpKK8jdW0Jefgl5+aX1wZCWHO/tLbhwiKhgUIX8T9yw0e4PITYebn7ATZLT5caLPnXZvmU89flTdE/uTs6oHNKT0lulZAsCY8wVQVU5UFLhhUIJufklfHXyLHB+MGRnptA7UoKhaBeszYHNi6DmjDt+kD3VNbqLafxYwIajG3j0k0eJj41n1shZ9EvpF/YyLQiMMVekQDDkBoJhbwnHyl0wdK4PBhcOvgdDRSlsnA/rfu86nqb0gayH4ZZx0DbxgtX3lu1l8keTKTtbxksjXuL29NC2ym7IgsAY842gquwP3mOIxGCoqXKtK3JnuQvU2nV0Te4GTYIO3c5b9VjFMaaumMru47t5Ovtpftjnh2Ery4LAGPONpKrsKz5NXn5pfTgEgqFLh+BgSCEjJaF1g0EVCvIgbxZ8+QFIDNz4A3ccIX1A/Wqnq0/z+KePs/rwaibfMpnJt0wOS50WBMaYqNAwGHLzSyiKhGAo3Qfr5riup1Xl0CMLsqfADXdDTCzVddU8l/sc7+95n3uvu5ens5+mTUxoW3NYEBhjopKqkl98uv6MpLygYLimQ7v6M5KyMlPo1RrBUHkSvlgAa2dDWQF07AmDH4b+49H4ZHK25JCzJYeh3Yby4ogXSWxz4bGF5rIgMMYYLgyG3L0lFJ+6MBh6pyaSmhxPWnI8yfFxoQ+I2hr4+wfueoSDedA2GQaMh8E/YXHxJp7LfY7rr76eWSNnkZZw8esTLpcFgTHGNEJV2Vt0uv74Ql5+aX0wBLSNiyEtKd4FQ1Jb0pLjSU1yP+duu+VJzQmNQxvd9Qjb3wOtgxu+x6o+w3l851w6xndk9qjZZHbMbPFn9S0IRGQM8DIQC8xV1RcaPC7e498FKoAJqrrpYq9pQWCMCZfAMYbDZZUUnaqkuLyK4lNnKSo/S9GpsxSfqqKo/Cylp89S18imMz4u5oKgSEtq64WIC5PA8sS2seeHxolDsP73sGE+VJaxI/0mpiTUUBUTw8x/nMnAaxrdhl82X4JARGKBXcBooBBYD4xT1R1B63wXmIYLgsHAy6o6+GKva0FgjPFbbZ1yvMKFQiAoioOCInhZyekqGtvMtmtzYWikJsXTtV0NN5cu49q9b1BSUcDkrt0ojIvl+az/ZEzf+5td88WCIJyt/wYBe1Q13yviLeAeYEfQOvcAb6hLozwR6SgiXVX1SBjrMsaYFomNkfoN+KXU1imlp6uCwuLC4DhYWsGmA8cprQiERl+EXzIiZgtTCz7gzW7F/CzvWTau+zP/MX5ByD9POIMgHTgYdL8Q963/UuukA+cFgYhMAiYB9OzZM+SFGmNMuMTGiBsiSr50aNTU1lFav6dRRXF5f/adGsego5uJO/FbrkkLz7SX4QyCxo6YNNxBupx1UNU5wBxwQ0MtL80YYyJPXGwMnZPb0Tm5XYNHrgXuC9v7hnN2hEKgR9D97sDhZqxjjDEmjMIZBOuBPiLSW0TaAmOBJQ3WWQI8JE4WcMKODxhjTOsK29CQqtaIyCPActzpo6+q6nYRedh7fDawFHfG0B7c6aMTw1WPMcaYxoV1wlBVXYrb2Acvmx10W4Gp4azBGGPMxbX+DMrGGGMiigWBMcZEOQsCY4yJchYExhgT5a647qMiUgQcaObTU4HiEJYTKpFaF0RubVZX01hdTfNNrKuXqjba0/qKC4KWEJENX9d0yU+RWhdEbm1WV9NYXU0TbXXZ0JAxxkQ5CwJjjIly0RYEc/wu4GtEal0QubVZXU1jdTVNVNUVVccIjDHGXCja9giMMcY0YEFgjDFRLmqCQETGiMjfRWSPiDzpdz0AIvKqiBwTkb/5XUswEekhIp+IyE4R2S4ij/ldE4CItBORdSKyxavrWb9rCiYisSLyhYj81e9aAkRkv4hsE5HNIhIxk31709K+KyJfer9n2RFQU1/v/ynwc1JEfup3XQAi8m/e7/zfRGSRiDScuaZlrx8NxwhEJBbYBYzGTYazHhinqjsu+sTw1zUcOIWbt/nbftYSTES6Al1VdZOIJAMbgXsj4P9LgERVPSUibYDPgcdUNc/PugJE5N+BgUAHVb3b73rABQEwUFUj6uIoEXkdWKWqc735ShJUtczvugK8bcYhYLCqNvcC1lDVko77Xe+nqmdE5G1gqaq+Fqr3iJY9gkHAHlXNV9Uq4C3gHp9rQlVXAqV+19GQqh5R1U3e7XJgJ24uaV+pc8q728b7iYhvMiLSHfgeMNfvWiKdiHQAhgPzAFS1KpJCwDMS2Ot3CASJA9qLSByQQIhncoyWIEgHDgbdLyQCNmxXAhHJAPoDa/2txPGGXzYDx4D/U9WIqAv4H+AJoM7vQhpQ4EMR2Sgik/wuxpMJFAHzvaG0uSKS6HdRDYwFFvldBICqHgJmAAXAEdxMjh+G8j2iJQikkWUR8U0ykolIEvAn4KeqetLvegBUtVZVb8XNbz1IRHwfUhORu4FjqrrR71oaMVRVBwDfAaZ6w5F+iwMGADmq2h84DUTEcTsAb6jq+8A7ftcCICJX40YwegPdgEQReTCU7xEtQVAI9Ai6350Q71p903hj8H8CFqrqYr/racgbSvgUGONzKQBDge974/FvAXeKyB/8LclR1cPev8eA93DDpH4rBAqD9ubexQVDpPgOsElVv/K7EM8oYJ+qFqlqNbAYGBLKN4iWIFgP9BGR3l7ajwWW+FxTxPIOys4DdqrqS37XEyAiaSLS0bvdHvcH8qW/VYGq/lxVu6tqBu5362NVDek3tuYQkUTvYD/e0MtdgO9nqKnqUeCgiPT1Fo0EfD0RoYFxRMiwkKcAyBKRBO9vcyTuuF3IhHXO4kihqjUi8giwHIgFXlXV7T6XhYgsAkYAqSJSCPyXqs7ztyrAfcMdD2zzxuMBnvLmoPZTV+B174yOGOBtVY2YUzUjUBfgPbftIA54U1WX+VtSvWnAQu+LWT4w0ed6ABCRBNzZhT/xu5YAVV0rIu8Cm4Aa4AtC3GoiKk4fNcYY8/WiZWjIGGPM17AgMMaYKGdBYIwxUc6CwBhjopwFgTHGRDkLAhPVROS3wR0mRWS5iMwNuv+iiDzlnb7XlNedICKvhLJWY8LFgsBEuzV4V2mKSAyQCtwY9PgQYIWq3u9Dbca0CgsCE+1Wc+5y/RtxV96Wi8jVIhIPfAs4Hpgzwvumv1hElonIbhH578ALichEEdklIp/hLsoLLO8lIitEZKv3b0+veV6+OB1FpC7QB0hEVonIda30+Y2xIDDRzevFUyMiPXGBkIvrtJqNm1tgK1DV4Gm3Aj8CbgJ+JG4in67As7gAGA30C1r/FdycEzcDC4GZqlqLmyOjH3A7bs6HYV74dFfVPeH4vMY0xoLAmHN7BYEgyA26v6aR9Veo6glVrcT1yOkFDAY+9RqDVQF/DFo/G3jTu70At+EHWIXryz8c+LW3/DZcbyxjWo0FgTHnjhPchBsaysNtvIfgQqKhs0G3aznXs+ty+7UE1lsFDMN1BF0KdMT1nlp5+aUb03IWBMa4jf3dQKk330EpbqOcjds7uBxrgREikuK18H4g6LE1uK6kAP+Mm3Yw8JwhQJ23d7EZ1+xsVUs+jDFNZUFgDGzDnS2U12DZicud61dVjwDP4ILjI1ynyIBHgYkishXX1fUx7zlncTPnBd53FZDsvbcxrca6jxpjTJSzPQJjjIlyFgTGGBPlLAiMMSbKWRAYY0yUsyAwxpgoZ0FgjDFRzoLAGGOi3P8DzxpvmRd4IFoAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "labels = ['1D', '2D', '3D']\n", "\n", "dres_fig2, dres_ax2 = plt.subplots()\n", "for data, label in zip(dres_conv2.T, labels):\n", " plt.plot(data, label=label)\n", "dres_ax2.set_xlabel('Window')\n", "dres_ax2.set_ylabel('Jensen-Shannon divergence')\n", "plt.legend()" ] }, { "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.\n", "\n", "[4] Matteo Tiberti, Elena Papaleo, Tone Bengtsen, Wouter Boomsma, and Kresten Lindorff-Larsen.\n", "ENCORE: Software for Quantitative Ensemble Comparison.\n", "PLOS Computational Biology, 11(10):e1004415, October 2015.\n", "00031.\n", "URL: https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1004415, doi:10.1371/journal.pcbi.1004415." ] } ], "metadata": { "kernelspec": { "display_name": "Python (mdanalysis)", "language": "python", "name": "mdanalysis" }, "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.7.3" }, "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": false } }, "nbformat": 4, "nbformat_minor": 2 }