doc/plotting-compare.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Building plots from `sourmash compare` output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "### Running this notebook.\n",
    "\n",
    "You can run this notebook interactively via mybinder; click on this button:\n",
    "[![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/dib-lab/sourmash/master?filepath=doc%2Fplotting-compare.ipynb)\n",
    "\n",
    "A rendered version of this notebook is available at [sourmash.readthedocs.io](https://sourmash.readthedocs.io) under \"Tutorials and notebooks\".\n",
    "\n",
    "You can also get this notebook from the [doc/ subdirectory of the sourmash github repository](https://github.com/dib-lab/sourmash/tree/master/doc). See [binder/environment.yaml](https://github.com/dib-lab/sourmash/blob/master/binder/environment.yml) for installation dependencies.\n",
    "\n",
    "### What is this?\n",
    "\n",
    "This is a Jupyter Notebook using Python 3. If you are running this via [binder](https://mybinder.org), you can use Shift-ENTER to run cells, and double click on code cells to edit them.\n",
    "\n",
    "Contact: C. Titus Brown, ctbrown@ucdavis.edu. Please [file issues on GitHub](https://github.com/dib-lab/sourmash/issues/) if you have any questions or comments!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running `sourmash compare` and generating figures in Python\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First, we need to generate a similarity matrix with compare. (If you want to generate this programmatically, it's just a `numpy` matrix.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[K\n",
      "== This is sourmash version 3.3.2.dev9+g462bc387. ==\n",
      "\u001b[K== Please cite Brown and Irber (2016), doi:10.21105/joss.00027. ==\n",
      "\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR2060939_1.sig'g'\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR2060939_2.sig'g'\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR2241509_1.sig'g'\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR2255622_1.sig'g'\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR453566_1.sig'g'\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR453569_1.sig'g'\n",
      "\u001b[Kloaded 1 sigs from '../tests/test-data/demo/SRR453570_1.sig'g'\n",
      "\u001b[Kloaded 7 signatures total.                                                     \n",
      "\u001b[K\n",
      "0-SRR2060939_1.fa...\t[1.    0.356 0.078 0.086 0.    0.    0.   ]\n",
      "1-SRR2060939_2.fa...\t[0.356 1.    0.072 0.078 0.    0.    0.   ]\n",
      "2-SRR2241509_1.fa...\t[0.078 0.072 1.    0.074 0.    0.    0.   ]\n",
      "3-SRR2255622_1.fa...\t[0.086 0.078 0.074 1.    0.    0.    0.   ]\n",
      "4-SRR453566_1.fas...\t[0.    0.    0.    0.    1.    0.382 0.364]\n",
      "5-SRR453569_1.fas...\t[0.    0.    0.    0.    0.382 1.    0.386]\n",
      "6-SRR453570_1.fas...\t[0.    0.    0.    0.    0.364 0.386 1.   ]\n",
      "min similarity in matrix: 0.000\n",
      "\u001b[Ksaving labels to: compare-demo.labels.txt\n",
      "\u001b[Ksaving comparison matrix to: compare-demo\n"
     ]
    }
   ],
   "source": [
    "!sourmash compare ../tests/test-data/demo/*.sig -o compare-demo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline\n",
    "# import the `fig` module from sourmash:\n",
    "from sourmash import fig"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `sourmash.fig` module contains code to load the similarity matrix and associated labels:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "matrix, labels = fig.load_matrix_and_labels('compare-demo')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, `matrix` is a numpy matrix and `labels` is a list of labels (by default, filenames)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "matrix:\n",
      " [[1.    0.356 0.078 0.086 0.    0.    0.   ]\n",
      " [0.356 1.    0.072 0.078 0.    0.    0.   ]\n",
      " [0.078 0.072 1.    0.074 0.    0.    0.   ]\n",
      " [0.086 0.078 0.074 1.    0.    0.    0.   ]\n",
      " [0.    0.    0.    0.    1.    0.382 0.364]\n",
      " [0.    0.    0.    0.    0.382 1.    0.386]\n",
      " [0.    0.    0.    0.    0.364 0.386 1.   ]]\n",
      "labels: ['SRR2060939_1.fastq.gz', 'SRR2060939_2.fastq.gz', 'SRR2241509_1.fastq.gz', 'SRR2255622_1.fastq.gz', 'SRR453566_1.fastq.gz', 'SRR453569_1.fastq.gz', 'SRR453570_1.fastq.gz']\n"
     ]
    }
   ],
   "source": [
    "print('matrix:\\n', matrix)\n",
    "print('labels:', labels)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `plot_composite_matrix` function returns a generated plot, along with the labels and matrix as re-ordered by the clustering:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAr4AAAFxCAYAAACV0auJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzde7ycZXnv/8+Xc0vAIFGUgyAnUamkAspW1HRbFfzJtoj1gKdQRa1lV7dH1CJYe6BWLFZFRFEogohnVBSsSlUQJdagRMEigkQJEM7BgBCu3x/PvWQc1lpZa0iYNVmf9+v1vFjP3Ifnngkrueaa67knVYUkSZK0rltv2AuQJEmS7g8GvpIkSZoVDHwlSZI0Kxj4SpIkaVYw8JUkSdKsYOArSZKkWcHAV5IkSTNOko8luTbJxRO0J8m/J7ksyY+TPHZ1cxr4SpIkaSY6Cdhvkvb9gV3a8UrgQ6ub0MBXkiRJM05VfRu4YZIuzwb+ozoXAHOTPHSyOTeYrHHevHm1ww47THuh0n31wx/+cHlVPWjY69D9a968LWqHHbYZ9jK0Dru77hz2EtagDHsBGsfin1w77CWsMXXnjZP+W7z+Hz20atUdg869BLi956ETquqEaU6zDXBVz/nS9tjVEw2YNPDdYYcdWLRo0TTXIN13Sa4c9hp0/9thh21YtOhzw16G1mG33rl02EtYY9ab/J/w0ZJ1J4ift+Nxw17CGnP7VadP+m9xrbqDjR/yjEHnvr2q9hpo8D3G+x+nJhuwDv3WSJIk6X6TkAy1anYpsF3P+bbAbyYbYI2vJEmSpi1AWG+gYw05E3hp291hH+DmqpqwzAHM+EqSJGkGSvJJYAEwL8lS4EhgQ4CqOh44C3gmcBnwW+CQ1c1p4CtJkqQBrN1Sh6p64WraC/ib6cxp4CtJkqSBDLnGd9oMfCVJkjSQjNiOHKMVpkuSJEkDMuMrSZKkAYRRy6Ea+EqSJGkg1vhKkiRpnZeMXuA7WquVJEmSBmTGV5IkSQPImvwWtvuFga8kSZIGYqmDJEmSNAOZ8ZUkSdIA1u5XFq8NBr6SJEkaiIGvJEmS1nmhu71tlIxWmC5JkiQNyIyvJEmSBjB6Nb6jtVpJQ5XkwCRXJVmR5E+HvR5J0hC1b24b5BgWA1+NpCT7Jjk/yc1JbkhyXpK9kyxMsqoFZrckuSjJs3rG7ZCkWvuKJFckObynfeMkJya5MsmtSX6UZP+e9n2SfL1d87okn07y0J72o5Lc2TP/iiQ79rRfkWRlT9s5fc/rQUlOS3JTkhuTnNrT9p4k/9PWdUmSl/a07Zrki21NNyQ5O8kjpvA67t76Lk9SU3jp3wMcVlVzqupHU+g/0XUryc6DjpckzQwGvtJalmRz4MvA+4EHAtsA7wTuaF2+V1VzgLnAccDpSeb2TTO39XkucESSp7XHNwCuAp4CPAA4AjgjyQ6tfQvgBGAHYHvgVuDjfXN/qgWGY8flfe0H9LQ9va/tc8CyNveD6QLNMbcBB7R1vQx4X5InjD0f4EzgEcBWwA+AL7J6dwJnAC+fQl/aupZMsa8kSTPKjKnxPeEEOO20Ya9CI2JXgKr6ZDtfCZwDkOSxY52q6u4kpwDHA7sAF/ZPVFWLkiwB5gNfr6rbgKN6unw5yS+BPYErquqrveOTfAD4rzXxpJI8HdgOWFBVq9rDv8+qVtWRPd2/n+Q7wP8Czq+qH9AFu2Nz/Rvwd0m2rKrrJ7pmVV0KXLq67GuSjYHrgfWBi5Isq6qdWrb8ULog/Srg7VX1+TZmZ+BEutf2TuAbVfX8JN9u017Usswvr6pPJXnTnns+evIXSZI0g4RRy6HOmNWedhosXjzsVWhE/BxYleTkJPsn2WK8TknWBw6hC7qunKDPPsDuwGUTtG9FF2hPlOV88jhtB7RygyVJ/nqcMae2koRzkuzR8/g+wKXAyUmuT3JhkqdMsK4/AvZezbqWTRb0TkdV3dEy5AB7VNVO7edfAE+iy0K/E/hET+nHu+jekGwBbEuXoaeqntwzz5wW9O4HvHFNrFWSdP8ZtVKHGZPxBZg/H849d9ir0EyQSbYFrKpbkuwLvAX4CPCQJGfRZR4B9klyE7ApcBfw4qq6tm+a5S2LuQlwDPCFe68hGwKnAidX1SXjtD8GeAfw7J6Hz6ArhbgGeDzw2SQ39WSnXwT8N93b5NcCZyfZrapuogsOnw68gi5gPwj4YpKdq2p53+WPBy4Czh5nXdsCHwRe39+2plXVp3tOP5XkrcDj6Mos7qQrjdi6qpYC351kqufRlYy8ZW2tVZK0ZsVdHaT7R1X9rKoWVtW2dBnbrYFjW/MFVTWXLtN4Jl1Gst88YA5dlnEBsGFvY7rf5FOA3wGH9Q9uH+N/FXhtVX2nZ10/rarfVNWqqjofeB9dHfFY+3lVtbKqfltV/wzc1LO+lXTlFCdW1Z1VdTpd+cAT+679r+05P6+qqq/tQXRZ1uN6gu21JslLkyxuN+Pd1NY1rzW/mS7A/0HLfv/VJFNtTfdcJUlaawx8NfJaNvYkuqCr9/EVwGuAl4y39VYLTo8Bbm/9AEgSutrUrYCDqurO3nFJtgf+E3hXVZ2yuuXBpF9r09v+43Y+oSTvBPYHnl5Vt/S1bUEX9J5ZVf+4mnXdZ+11+AjdG4Mt25uNi2nPp6qWVdWhVbU18CrguElqia+mq2+WJI2QsN5Ax7AY+GrkJNktyRvaR/ok2Q54IXBBf99W4/pRupKEiRwNvDnJJu38Q8Aj6XZfWNl37W2AbwIfrKrjx1nbs5Nskc7jgL+l7a6Q5GFJnphkoySbJHkTXXb0vDb888AWSV6WZP0kz6XbseK8Nv6twMHA0/prd9tOF2cD51XV4UxRW+cmwEbtfJNWAjIVm9IF6te1sYfQ8+YjyV+O/RkBN7a+YzftXQPseM9UnAEsnOq6JUkzw6jV+Br4ahTdSlc/+/0kt9EFvBcDb5ig/7HAM1tN7ni+QheYHdqymK+i24lgWe7Zb/dFre8r6AK2I3vaVvTM9QK6G+VuBf4D+JeqOrm1bUYXVN8I/BrYD9h/LIitqhuA/0NXfnEzcDjw7J763n8CHgb8T8+139baDqS72e2Q/OEewg+b+GUEuhrcldxzk9xKuhvsAEjy1Z5r/IGq+ildffT36ALZP+GeIJ62nu+31+dMurKQX7a2o+hu4rspyfPabhnHIkkaHQkZ8BjakvtKBP/AXnvtVYsWLbpfFrJgQfdfb24TQJIfVtVew16H7l977bV7LVr0uWEvQ+uwW+9cOuwlrDHrzaz70++bIQZCa9q8HY8b9hLWmNuvOn3Sf4s3+uNtaqudXzXQ3Et/cuRQ/p1fh35rJEmSdH9yVwdJM0orV1gxzjFuCYMkSVMRRu/mNjO+0jquqvYf9hokSesi9/GVJEmSZiQzvpIkSRrIqGV8DXwlSZI0gAy1XncQBr6SJEkazIhlfEdrtZIkSdKAzPhKkiRp+mKNryRJkmaBMNyvHx7EaIXpkiRJ0oDM+EqSJGkg7uogSZKkWcEaX0mSJM0CAWt8JUmSpJnHjK8kSZKmL4xcCtXAV5IkSYOx1GH6TjgBFi8e9iokSZI0Lclgx5DMiIzvaafBzTcPexWSpHXdejPjn701Y8QybZOZt+Nxw17CGrP88tcMewlrzJwNTx/2Eta4GfM3wAMeMOwVSJIkaVpmRO3A1M2YwFeSJEkjJFAj9smDga8kSZIGM1px76glqCVJkqTBmPGVJEnSYNYbrZSvga8kSZIGMHpfWWzgK0mSpOkL1vhKkiRJM5EZX0mSJA3GGl9JkiTNCtb4SpIkaVYYrbjXGl9JkiTNDmZ8JUmSNH3BGl9JkiTNEqMV9xr4SpIkaTA1Yje3WeMrSZKkGSfJfkkuTXJZksPHaX9Aki8luSjJkiSHrG5OM76SJEmavmSt1fgmWR/4IPA0YClwYZIzq+qnPd3+BvhpVR2Q5EHApUlOrarfTTSvGV9JkiQNJgMeq/c44LKqurwFsqcDz+7rU8BmSQLMAW4A7ppsUgNfSZIk3d/mJVnUc7yyr30b4Kqe86XtsV4fAB4J/Ab4CfDaqrp7sota6iBpypIcCPw7sAXwpKr60ZCXJEkapsFvblteVXtNNvM4j1Xf+TOAxcD/BnYCvp7kO1V1y0STmvHVSEqyb5Lzk9yc5IYk5yXZO8nCJKuSrEhySyt4f1bPuB2SVGtfkeSK3oL5JBsnOTHJlUluTfKjJPv3tO+T5Ovtmtcl+XSSh46zvo2SXJJkad/j70rykyR3JTmqr21Bkrt71rYiycv61vax9ryWJXl93/gDklzcxp2f5FFTeB13T3J2kuVJ+v9CGc97gMOqas59CXrbn8HOg46XJM0Q62WwY/WWAtv1nG9Ll9ntdQjwuepcBvwS2G3S5U7jqUkzQpLNgS8D7wceSPfRxzuBO1qX71XVHGAucBxwepK5fdPMbX2eCxyR5Gnt8Q3oPlp5CvAA4AjgjCQ7tPYtgBOAHYDtgVuBj4+zzDcB147z+GXAm4GvTPD0ftOCyrHj5J62o4Bd2nX/DHhzkv3aa7ILcCrw6va8vwScmWR1n+rcCZwBvHw1/cZsDyyZYl9J0rps0PreqSWJLwR2SfLwJBsBLwDO7OvzK+CpAEm2Ah4BXD7ZpDOq1GHxYliwYNir0AjYFaCqPtnOVwLnACR57Finqro7ySnA8XQB44X9E1XVoiRLgPnA16vqNroAc8yXk/wS2BO4oqq+2js+yQeA/+p77OHAi4HXAx/pu97Jrc+LpveUAXgpcEhV3QjcmOQjwELga3Qf93ynqr7b5v8X4B10Afw3Jpqwqi6luwt20uxrko2B64H1gYuSLKuqnVq2/FDgwXRvGN5eVZ9vY3YGTqR7be8EvlFVz0/y7TbtRS3L/PKq+lSSN+2556On/6pIktY5VXVXksOAs+n+7flYVS1J8urWfjzwLuCkJD+hC6ffUlXLJ5t3xgS+W20FD73XB8bSuH4OrEpyMt1dnhe0YPAPtK1QDqELuq4cb6Ik+wC7A/88QftWdIH2RFnOJ4/T9n7gbXQB+XQ9OMk1wG+BLwB/V1W3JdkC2Bq4qKfvRcBfjC2VP3wPPXa+O5MEvlNVVXcAc1qgukf7SAngF8CTgGXAXwKfSLJzVV1N9xfSOXTZ6Y2AvdpcT+6fp2Wu33hf1ylJup+txS+wqKqzgLP6Hju+5+ffAE+fzpwzJvB96EPh3HOHvQrNFJP9HlXVLUn2Bd5Cl1F9SJKz6DKPAPskuQnYlG5bkxdXVX/ZwfKWxdwEOIYuyOxbQzakKx84uaouGaf9MXRZ1Wf3PHYgsEFVfT7Jgqk929+7hC47egldScHJwHuBV9Ft0wJwc0//m4HN2s9fB45u1zyf7rXZCPjjaa5hWqrq0z2nn0ryVrotaL5I94Zje2DrqloKfHeSqZ5HVzLylrW1VknSmpa1GviuDdb4aiRV1c+qamFVbUuX1dwaOLY1X1BVc+nqcc+ky0j2m0cXTL4RWABs2NuYZD3gFOB3wGH9g9vH+F+l2zrlO+2xTYF3A/93wOe0rKp+WlV3V9Uv6WqBn9uaV7T/bt4zZHO6GmNaYP4yuq1drm7P76d0NwesNUlemmRxkpvam43d27Vp6w/wg3TfqPNXk0y1NX+4bY0kaaYLXSQ5yDEkBr4aeS3oO4ku6Op9fAXwGuAlSf50nHGrquoY4PbWD4AkoatN3Qo4qKru7B2XZHvgP4F3VdUpPU270N309p0ky4DPAQ9tOzDsMMhTo5UvtFKOq4E9etr3oKfMoqo+U1W7V9WWwJF02dZ71TWvKe11+AjdG4Mt25uNi3vWvKyqDq2qremy1sdNUkt8NX94964kSWucga9GTpLdkrwhybbtfDvghcAF/X2r6nrgo3QlCRM5mm6HhE3a+YfoNsQ+oKr+oE43yTbAN4EP9tYZNRfTBW/z2/EK4Jr281Vt/IbtOusBGyTZpNUij21n9rB0tmvr+mLP/P8B/F2SLZLsRlfacVLP2vZMsn66r238MPCl8Uo0+p5P2no2auebtBKQqdiULji/ro09hJ43H0n+cuzPCLix9V3Vzq8BduyZ6wy6G/UkSaMkGewYEgNfjaJbgccD309yG13AezHwhgn6Hws8s9XkjucrdIHZoS2L+Sq6YHVZz366Y7swvIIuYDuyd79d6O5AbVnOZVW1jO6rE+9u52MB30fobnp7IfD29vNLWttjge8Bt9HV6V4M/G3POo+ku5nsSrqdJP61qr7W0/4+4Cbg0vbfQ1m97dsaxjLHK9t4AJJ8NcnbxhvYvi/9mLbma4A/Ac7r6bI33Z/RCrqSk9e2Eg7ods44uZVIPK/tlnEskqTRsva2M1srUjXxnvV77bVXLVq0aK0vYmwLM29u05gkP1zNN7poHbTXXrvXokWfG/YytA677c5lw17CmjNiNxVNZt6Oxw17CWvM8stfs/pOI2LOhk+e9N/ijR+4Y237jH8YaO7LT3/RUP6dN+MrSZKkWcHAV1rHtXKFFeMc45YwSJI0JWHkanxnzD6+ktaOqtp/2GuQJK2jRqzixsBXkiRJg1lvtCJfSx0kSZI0K5jxlSRJ0gBG7yuLDXwlSZI0fUPek3cQBr6SJEkajDW+kiRJ0sxjxleSJEmDGbGMr4GvJEmSpi9QoxX3GvhKkiRpQCOW8bXGV5IkSbOCGV9JkiQNxn18JUmStM5LLHWQJEmSZiIzvpIkSRrMiKVQBw58TzgBTjttzSxi8WKYP3/NzCVJkqT7yWyp8T3tNANWSdJo+eDPbh/2EtaYI5/18WEvYY1Zfvlrhr2ENadq2Cu4/4SRq/G9T6UO8+fDuefe90UsWHDf55AkSZImY42vJEmSBlKzpdRBkiRJs9xsublNkiRJs9gI1viOWJwuSZIkDcaMryRJkgaQ2bOdmSRJkma5ESt1MPCVJEnSYEYr7rXGV5IkSbODGV9JkiRNX6AsdZAkSdKsYOArSZKkWWHEdnWwxleSJEmzghlfSZIkTV8YuRSqga8kSZIG4BdYSJIkabYYsZvbRixBLUmSJA3GjK8kSZKmL4xcxtfAV5IkSQMpa3wlSZI0K4xY0eyILVfSMCU5MMlVSVYk+dNhr0eSpOkw8NVISrJvkvOT3JzkhiTnJdk7ycIkq1pgdkuSi5I8q2fcDkmqta9IckWSw3vaN05yYpIrk9ya5EdJ9u+79lOTXJLkt0m+lWT7vvbHJvl2m/+aJK/tu/632thLkvx5T1uSvD3Jr9raT0+yeU/7u1vQeUtb39v7rntAkovbdc9P8qgpvI4vS/LDNufSdo3JPgl6D3BYVc2pqh+tbv5JrltJdh50vCRpBgjddmaDHENi4KuR04LBLwPvBx4IbAO8E7ijdfleVc0B5gLHAacnmds3zdzW57nAEUme1h7fALgKeArwAOAI4IwkO7RrzwM+1x5/ILAI+FTP2uYBXwM+DGwJ7Ayc03PdTwI/am1vBz6T5EGt7aXAS4AnAlsDf9Se45gTgd2qanPgCcDBSZ7TrrsLcCrw6va8vwScuZogFuCPgdcB84DHA08F3jhJ/+2BJauZU5I0K6S7uW2QY0is8Z2qE06A004b9irU2RWgqj7ZzlfSgsskjx3rVFV3JzkFOB7YBbiwf6KqWpRkCTAf+HpV3QYc1dPly0l+CewJXAE8B1hSVZ9u1zsKWJ5kt6q6BHg9cHZVndrG3wH8rPXdFXgs8PSqWgl8NsnrgIPaGg8ATqyqq1r/fwG+meSvq+q3VXVp3/LvpgusAZ4BfKeqvtsz9h10Afw3Jnohq+pDPae/TnIq8Gf9/ZJsDFwPrA9clGRZVe3UsuWHAg+me8Pw9qr6fBuzM12wPh+4E/hGVT0/ybfbtBclKeDlVfWpJG/ac89HT7RUSZLuMzO+U3XaabB48bBXoc7PgVVJTk6yf5ItxuuUZH3gELqg68oJ+uwD7A5cNkH7VnSB9liW89HARWPtLVD+RXscYB/ghlZqcG2SLyV5WM/Yy6vq1p5LXNQzNu2g53xjuqB9bD2HJ1kBLAU2BU6bZGzac5uOJzNORreq7mgZcoA9qmqn9vMvgCfRZcffCXwiyUNb27vo3pBsAWxLy15X1ZN75pnTgt79mDzTLEmaicz4rsPmz4dzzx32KmaHSep/quqWJPsCbwE+AjwkyVl0mUeAfZLcRBcY3gW8uKqu7ZtmectibgIcA3zh3kvIhnTlAye3bC7AHOC6vq43A5u1n7ely+o+DfgJ8G668oYntrE3jzN2m/bzV4E3JzkDuLE9P+jKEcae+9Etmzsf+Iue+b4OHJ1kAXB+G7tR79jVSXIIsBfwiqmOGct8N59K8lbgccAX6d5wbA9sXVVLge9OMtXzgI9zz3OWJM10/SmXEWDGVyOpqn5WVQuralu6rObWwLGt+YKqmkuXaTyTLiPZbx5dIPpGYAGwYW9jkvWAU4DfAYf1NK0ANucPbQ6MZXFXAp+vqgur6na6LOgTkjxgCmM/Rhckn0uXdf1We3xp33OvdmPZyjY/LTB/GfAB4Or2/H7aP3YiSf4COBrYv6qWT2VMG/fSJIuT3NTebOzerg3wZrq/En+QZEmSv5pkqq3pSiUkSSOigFovAx3DYuCrkdeCvpPo+1i/qlYArwFeknG23qqqVVV1DHB76wd0uyvQ1aZuBRxUVXf2DFsC7NHTd1NgJ+4pD/gx3d8Fv7/MWNfWZ8ckm/W07zE2tqrurqojq2qHFtAvAX7djvFs0K499nw+U1W7V9WWwJF02dZ71TX3a2UGHwEOqKqfrK5/z7jt27jDgC3bm42L23OlqpZV1aFVtTXwKuC4SXZyuBrYbqrXliRpEAa+GjlJdkvyhiTbtvPtgBcCF/T3rarrgY/S3eg1kaPpSgw2aecfAh5JFwiu7Ov7eWD3JAe1/u8AftxTCvFx4MAk81upxBHAd6vqpqr6ObAYODLJJkkOBB4DfLY9jwcm2alta/Yo4L3A37eb9NZL8qokW7T2xwF/Q8+Na0n2TLJ+2yXiw8CXetY10Wv5v+nKOQ6qqh9M1nccm9IF9te1uQ6h581Hkr8c+zOiK90oYFU7vwbYsWeuM4CF07y+JGnY3M5MWutupdt66/tJbqMLeC8G3jBB/2OBZyZ5zATtX6ELzA5tWcxX0dXQLss9+/2+CKCqrqPbheEf25jHAy8Ym6iqvgm8rc15Ld2uCwf3XOsFdHW0N9IF3M9tc0JXInAWcBtdve/HquqEnrEH0t1MdivwCbqbxXq3O3sfcBNwafvvoazeEXQ3pp3V81y/OtaY5KtJ3jbewKr6KV199PfoAtk/Ac7r6bI33Z/RCrqSk9dW1S9b21HAya1E4nlV9VXuKVWRJI2C4M1t0tpWVb+muxlqPCe1o7f/UrrdEcakr724Z2eFe7WPc/3/BHabpP1DdFnj8dquoKspHq/t58AjJmi7G9hvNevad7L2Ccbca+uyvvb9+877X7u30+1HPN7YN9PV+Y7XdjzdFm69jx291167//MUli1Jmim8uU2SJEmaeQx8pXVcK1dYMc4xbgmDJElTEWC99QY7hsVSB2kd11+uIEnSmjLE+9QGYuArSZKk6RvuBg0DsdRBkiRJM06S/ZJcmuSyJIdP0GdB+yKlJUn+a3VzmvGVJEnSAELWUso3yfrAB4Gn0X0L6YVJzmxbaY71mQscB+xXVb9K8uDVzWvGV5IkSdMW1ur3VzwOuKyqLq+q3wGnA8/u63Mw8Lmq+hVAVV27ukkNfCVJkjSQ+xD4zkuyqOd4Zd/U2wBX9ZwvbY/12hXYIsm5SX6Y5KWrW6+lDpIkSbq/La+qvSZpHy8vXH3nGwB7Ak8F/gj4XpIL2hdCjcvAV5IkSdMXyNqrHVgKbNdzvi3wm3H6LK+q24Dbknwb2AOYMPC11EGSJEkDWYs1vhcCuyR5eJKNgBcAZ/b1+SLwpCQbJPlj4PHAzyab1IyvJEmSpi3AemtpH9+quivJYcDZwPrAx6pqSZJXt/bjq+pnSb4G/Bi4G/hoVV082bwGvpIkSZpxquos4Ky+x47vO/9X4F+nOqeBryRJkgYyat/cZuArSZKkgRj4SpIkad0X1to3t60tMybwXbwYFiwY9iomsfhYzp3/umGvQpIkSQOaEYHvwQcPewWSpNngyGd9fNhLWGPe+eVDhr2ENaf6v5dgdN3NXcNewv1qLe7ju1bMiMD3la/sjhltgdleSZKkMWH0anxHLE6XJEmSBjMjMr6SJEkaMVP/FrYZw8BXkiRJAzHwlSRJ0qywtr6yeG2xxleSJEmzghlfSZIkTdso7upg4CtJkqSBGPhKkiRp3RfIiBX5WuMrSZKkWcGMryRJkgZiqYMkSZJmBQNfSZIkrfNGcVcHa3wlSZI0K5jxlSRJ0vRl9L65zcBXkiRJAxm1UgcDX0mSJA0kI1Y0O2LLlSRJkgZjxleSJEnTNoq7Ohj4SpIkafoCGbHI11IHSVOW5MAkVyVZkeRPh70eSdJwJYMdw2Lgq5GUZN8k5ye5OckNSc5LsneShUlWtcDsliQXJXlWz7gdklRrX5HkiiSH97RvnOTEJFcmuTXJj5Ls33ftpya5JMlvk3wryfZ97Y9N8u02/zVJXtt3/W+1sZck+fOetiR5e5JftbWfnmTznvZ3t6Dzlra+t/dd94AkF7frnp/kUVN4HXdPcnaS5UlqCi/9e4DDqmpOVf1oCv0num4l2XnQ8ZIkDcLAVyOnBYNfBt4PPBDYBngncEfr8r2qmgPMBY4DTk8yt2+aua3Pc4EjkjytPb4BcBXwFOABwBHAGUl2aNeeB3yuPf5AYBHwqZ61zQO+BnwY2BLYGTin57qfBH7U2t4OfCbJg1rbS4GXAE8Etgb+qD3HMScCu1XV5sATgIOTPKdddxfgVODV7Xl/CTgzyerKme4EzgBevpp+Y7YHlkyxryRpHTdqGV9rfKdj8WJYsGDYqxDsClBVn2znK2nBZZLHjnWqqruTnAIcD+wCXNg/UVUtSrIEmA98vapuA47q6fLlJL8E9gSuAJ4DLKmqT7frHQUsT7JbVV0CvB44u6pObePvAH7W+u4KPBZ4elWtBD6b5HXAQW2NBwAnVtVVrf+/AN9M8tdV9duqurRv+eWHJzAAACAASURBVHfTBdYAzwC+U1Xf7Rn7DroA/hsTvZBtzktXl31NsjFwPbA+cFGSZVW1U8uWHwo8mO4Nw9ur6vNtzM50wfp8ugD7G1X1/CTfbtNe1LLML6+qTyV50557PnqyZUiSZpgRK/E14ztlBx8M8+cPexXq/BxYleTkJPsn2WK8TknWBw6hC7qunKDPPsDuwGUTtG9FF2iPZTkfDVw01t4C5V+0xwH2AW5opQbXJvlSkof1jL28qm7tucRFPWPTDnrON6YL2sfWc3iSFcBSYFPgtEnGpj23+6yq7mgZcoA9qmqn9vMvgCfRZcffCXwiyUNb27vo3pBsAWxLy15X1ZN75pnTgt79gDeuibVKku4fofvmtkGOYTHjO1WvfGV36P4xyVvIqrolyb7AW4CPAA9JchZd5hFgnyQ30QWGdwEvrqpr+6ZZ3rKYmwDHAF+49xKyIV35wMktmwswB7iur+vNwGbt523psrpPA34CvJuuvOGJbezN44zdpv38VeDNSc4AbmzPD+CPe5770S2bOx/4i575vg4cnWQBcH4bu1Hv2LVhLPPdfCrJW4HHAV+ke8OxPbB1VS0FvjvJVM8DPs49z1mSpDXOjK9GUlX9rKoWVtW2dFnNrYFjW/MFVTWXLtN4Jl1Gst88ukD0jcACYMPexiTrAacAvwMO62laAWzOH9ocGMvirgQ+X1UXVtXtdFnQJyR5wBTGfowuSD6XLsP8rfb40r7nXu3GspVtflpg/jLgA8DV7fn9tH/smpbkpUkWJ7mpvdnYvV0b4M10CYEfJFmS5K8mmWprulIJSdKoGDDbO8yMr4GvRl4L+k6i72P9qloBvAZ4yXhbb1XVqqo6Bri99QO63RXoalO3Ag6qqjt7hi0B9ujpuymwE/eUQvwY6N0dYezntD47Jtmsp32PsbFVdXdVHVlVO7SAfgnw63aMZ4N27bHn85mq2r2qtgSOpMu23quueU1pu1l8hO6NwZbtzcbFtJKLqlpWVYdW1dbAq4DjJqklvhrYbm2tVZK05nWlDjXQMSwGvho5SXZL8oYk27bz7YAXAhf0962q64GP0t3oNZGj6UoMNmnnHwIeCRzQbkLr9Xlg9yQHtf7vAH7cUwrxceDAJPNbqcQRwHer6qaq+jmwGDgyySZJDgQeA3y2PY8HJtmpbWv2KOC9wN+3m/TWS/KqJFu09scBf0PPjWtJ9kyyftsl4sPAl3rWNdFrmfY8Nmrnm7QSkKnYlC6wv66NPYSeNx9J/nLsz4iudKOAVe38GmDHnrnOABZO8bqSJA3EwFej6Fbg8cD3k9xGF/BeDLxhgv7HAs9M8pgJ2r9CF5gd2rKYr6KroV2We/b7fRFAVV1HtwvDP7YxjwdeMDZRVX0TeFub81q6XRcO7rnWC4C92tijgee2OaErETgLuI2u3vdjVXVCz9gD6W4muxX4BN3NYr3bnb0PuAm4tP33UFZve7qSibGM9co2HoAkX03ytvEGVtVP6eqjv0cXyP4JcF5Pl73p/oxW0JWcvLaqftnajgJObiUSz6uqr3JPqYokaUSMWqmDN7dp5FTVr+luhhrPSe3o7b+UbneEMelrL+7ZWeFe7eNc/z+B3SZp/xBd1ni8tivoaorHa/s58IgJ2u4G9lvNuvadrH2S9Uz4fKtq/77z/tfu7XT7EY839s10db7jtR1Pt4Vb72NH77XX7v88pYVLkmaEUcugGvhKkiRp2sZqfEfJqAXqkqaplSusGOcYt4RBkqR1lRlfaR3XX64gSdKaMsx63UEY+EqSJGnawuiVDhj4SpIkafqGvEPDIEYtUJckSZIGYsZXkiRJA8mI7epg4CtJkqRp67YzG/YqpsfAV5IkSQMZtZrZUVuvJEmSNBAzvpIkSZq2UCP3zW0GvpIkSRqINb6SJEmaFUatZnbU1itJkiQNxIyvJEmSpi0j+M1tBr6SJEkaiDe3SZIkaZ03il9gYY2vJEmSZgUzvpJmjKt/+zv+cfGVw17GGvHbu0YsDTKJ9z7nxGEvYY1Z+at3DnsJ0jpl1DKoBr6SJEkaiDW+kiRJWudZ4ytJkiTNUGZ8JUmSNH3u4ytJkqTZIIxe6cCorVeSJEkaiIGvJEmSBrJeaqBjKpLsl+TSJJclOXySfnsnWZXkuaub01IHSZIkDWRt1fgmWR/4IPA0YClwYZIzq+qn4/T7F+DsqcxrxleSJEnTNlbjO8gxBY8DLquqy6vqd8DpwLPH6fd/gc8C105lUgNfSZIk3d/mJVnUc7yyr30b4Kqe86Xtsd9Lsg1wIHD8VC9qqYMkSZIGch9KHZZX1V6TtI83c39x8LHAW6pqVTK1hRj4SpIkaSBZe19ZvBTYrud8W+A3fX32Ak5vQe884JlJ7qqqL0w0qYGvJEmSpi1r9wssLgR2SfJw4NfAC4CDeztU1cPvWUtOAr48WdALBr6SJEmaYarqriSH0e3WsD7wsapakuTVrX3Kdb29DHwlSZI0kLW5S0JVnQWc1ffYuAFvVS2cypwGvpIkSZq2MPUvo5gpDHwlSZI0kLVY47tWuI+vJEmSZgUzvpIkSRqIGV9JIyvJVkm+neTWJMcMez2SpJkrdNstDHIMi4GvRlKSfZOcn+TmJDckOS/J3kkWJlmVZEWSW5JclORZPeN2SFKtfUWSK5IcPsE1dklye5JPTDJ+RZIjetpfl+Tydu3fJPm3JBv0tF+RZGXP2HN62t7WN+/KJHcnmdfaN07ysTb3siSvn+JrdUKSS9tcC1fT/ZXAcmDzqnrDVOaf4JonJfmHQcdLkkbDeqmBjqGtd2hXlgaUZHPgy8D7gQfSfXf3O4E7WpfvVdUcYC5wHN23usztm2Zu6/Nc4IgkTxvnUh+k20B7PHOrak473tXz+JeAx1bV5sDuwB7A3/aNPaBn7NPHHqyqf+p5fA7wL8C5VbW8dTkK2AXYHvgz4M1J9ptgfb0uAl4D/PcU+m4P/LSqRus2XUmSpsAa3wGd8MMTOO0npw17GbPVrgBV9cl2vhI4ByDJY8c6VdXdSU4BjqcLGO8VxFbVoiRLgPnA18ceT/IC4CbgfGDnqS6sqn7Rcxrg7umM77l+gJcAf9/z8EuBQ6rqRuDGJB8BFgJfW82aPtjmvH011zwJeBFQSV4H/AVwC/A+4JF0r/NngddX1e/aGt/bxmwMXEn3rTpP6JvnW1V1QJI/BU6k+7M4i+471y+rqr+bymsiSZpZ1vI3t60VZnwHdNpPTmPxssXDXsZs9XNgVZKTk+yfZIvxOiVZHzgEuJMuKBuvzz50mdnLeh7bnC7gnOyj/iuTLE3y8bFShJ7xBye5ha5kYA/gw31jT01yXZJzkuwxwfxPAraiCzRpz3FruuztmIuAR0+yxmlpm3+fCry7ZZ3/E1gF/D+670D/X8BT6bLHAE8Hnkz3RmQu8Hzg+qo6oW+eA5JsBHwBOIUuS/9p4KA1tXZJ0nCsl8GOYTHjex/Mf8h8zl147rCXsU7KIRP/VlTVLUn2Bd4CfAR4SJKzgENbl32S3ARsCtwFvLiqru2bZnmSjYFNgGPogrIx7wJOrKqruqTmH44D9gYWA1vSlUOcCjyjZ32nAacl2YUuS3tNz/gX0ZUcBHgtcHaS3arqpr7rvAz4TFWtaOdz2n9v7ulzM7BZ/wLXpKr6Yc/pFUk+DDwFOJbuDcVmwG7AD6rqZ5NMtQ+wIXBsK6P4zFRrlCVJM1OA9c34SmtfVf2sqhZW1bZ0Gdut6YIxgAuqai6wBXAmXfa03zy6YPKNwAK6oIwk84E/B/5tguuuqKpFVXVXVV0DHAY8vWWJ+/v+D7CErs547LHzqmplVf22qv6ZrpziD9aX5I+AvwRO7nl4LADuvc7mwK3jrXNNSbJrki+3m+luAf6J7rWjqr4JfIAu+L+m3UR3r9eh2Rr4dV/t8LhZeEmS1hYDX428qroEOIkuAO59fAXdx/IvafWl/eNWVdUxwO3c8/H9AmAH4FdJltEFxgclmejGsLFAbqL3vBsAO022/HHGPge4ATi3Z603AlfTlU6M2YMusF6bPgRcAuzSbth7Gz3rrap/r6o96UoudgXeNNbUN8/VwDb5wxT6w9baqiVJ94tRK3Uw8NXISbJbkjck2badbwe8ELigv29VXQ98FHjHJFMeTbdDwibACXSB6vx2HA98hVbKkOTxSR6RZL0kWwL/Trfzws2t/RVJHtx+fhTwVuAb7fxhSZ6YZKMkmyR5E1329Ly+9bwM+I9xdlb4D+DvkmyRZDe60o6TpvB6bdSeW4AN27Wn+ru/Gd0NbivaNf+6Z9692+uxIXAb3RuIVa35GmDHnnm+R1d28rdJNkjyHOBxU1yDJGkGCm5nJt0fbgUeD3w/yW10Ae/FTHwz2rHAM5M8ZoL2rwA3Aoe2EoRlYwddicHtVXVd67sj3S4Kt7Zr3kEXdI95IvCTtq6z2vG21rYZXQb1RuDXwH7A/i04ByDJNsD/pgty+x0J/IKuROC/gH+tqkl3dGjOoduR4Ql0gf1KupvSSPKitqvFRN5It1PDrXT11J/qadu8PXZjW9P1wHta24nAo5LclOQLVfU7ukz2wtb/+cDnep73w5KsQJI0OgbM9npzmzQNVfVr4HkTNJ9EXxa0qpbSbbc1Jn3txQS7I1TVUX3nnwQ+OV7f1n7IJG1LgImC77E+v2aC38uqugP4q3ZMWVUtmKTtVLqb88bOF/a1f5vu5rVe72ht32CC59Pqm+f3PbYI+H3JSds+baztV8CcbR61q/sHS5LWGgNfSZIkTdvYVxaPEksdpBHXyhVWjHOs7RvfJEmznKUOku5X/eUKo6K/rEKSpLXNwFeSJEkDGeYODYMw8JUkSdK0jeI3txn4SpIkaSDDrNcdhDe3SZIkaVYw4ytJkqRpy5B3aBiEga8kSZIGYuArSZKkdV53c9to7epgja8kSZJmBTO+kiRJGsioZVANfCVJkjRtwRpfSZIkzRKjFviOWoZakiRJGogZX0mSJE1bUiO3q4OBryRJkgYyaqUOBr6SJEmaNm9um0UWL1s87CVI65wCbl81Yn+LTuC9zzlx2EtYY17/uZcPewmStEYY+EqSJGkgZnwlSZK0zuu+snjYq5geA19JkiRNX2C9EdvVwX18JUmSNCuY8ZUkSdJARi2DauArSZKkaXM7M0mSJM0ao3Zz26hlqCVJkqSBmPGVJEnStHWlDqO1q4OBryRJkgYyajW+ljpIkiRpVjDjK0mSpGlLRi/ja+ArSZKkgYxa6YCBryRJkgaSEcv4jlqgLkmSJA3EjK8kSZIGMmIJXwNfSZIkTV8YvVIHA19JkiQNZNRqZkdtvZIkSdJAzPhKkiRpIBmxryw24yvp95JsleTbSW5Ncsyw1yNJmtky4DEsBr4aSUn2TXJ+kpuT3JDkvCR7J1mYZFWSFUluSXJRkmf1jNshSbX2FUmuSHL4BNfYJcntST4xyfgVSY7oaX9dksvbtX+T5N+SbNDTfkWSlT1jz+m75oOSnJbkpiQ3Jjm1p22bJF9sz3dpkldP8bU6IcmlSe5OsnA13V8JLAc2r6o3TGX+Ca55UpJ/GHS8JGnmG7u5bZBjWAx8NXKSbA58GXg/8EBgG+CdwB2ty/eqag4wFzgOOD3J3L5p5rY+zwWOSPK0cS71QeDCCZYxt6rmtONdPY9/CXhsVW0O7A7sAfxt39gDesY+va/tc8AyYHvgwcB7eto+AfwS2Ar4/4B/SvJnE6yv10XAa4D/nkLf7YGfVtVofXYlSdIUWON7H6z43QoWnLRg2MuYjXYFqKpPtvOVwDkASR471qmq7k5yCnA8sAvjBLFVtSjJEmA+8PWxx5O8ALgJOB/YeaoLq6pf9JwGuHuq45M8HdgOWFBVq9rDP2ptc4AFwPOq6k7goiSfAf4K+NZq1vTBNsftq7n+ScCLgEryOuAvgFuA9wGPpHudPwu8vqp+lyTAe9uYjYErgYOBJ/TN862qOiDJnwIn0v1ZnAUUcFlV/d1UXh9J0swzYruZmfEd1FabbsWcjeYMexmz1c+BVUlOTrJ/ki3G65RkfeAQ4E66oGy8PvvQZWYv63lsc+Dvgck+6r+ylRt8PMm8vjkPTnILXcnAHsCH+8aemuS6JOck2aPn8X2AS4GTk1yf5MIkTxmbtu+/Yz/vPskap6WqFgKnAu9u2ej/BFYB/w+YB/wv4Kl02WOApwNPpnsjMhd4PnB9VZ3QN88BSTYCvgCcQpel/zRw0JpauyRpCALrDXgMixnfAT10s4fy0M0eyrkLzx32UtZJOWTi34qquiXJvsBbgI8AD0lyFnBo67JPkpuATYG7gBdX1bV90yxPsjGwCXAMXVA25l3AiVV1Ve5diLQc2BtYDGxJVw5xKvCMnvWdBpyWZBfgpcA1PeNfRFdyEOC1wNlJdquqm4Bt6YLJV9AF7AcBX0yyc1UtT3IeXVnGm4BHtfbrJnyh1oCq+mHP6RVJPgw8BTiW7g3FZsBuwA+q6meTTLUPsCFwbCuj+EyS16+lZUuS7idmfKX7QVX9rKoWVtW2dFnPremCMYALqmousAVwJvCkcaaYB8wB3khXQrAhQJL5wJ8D/zbBdVdU1aKququqrgEOA57essT9ff8HWEJXZzz22HlVtbKqfltV/0xXTjG2vpXAFVV1YlXdWVWnA1cBT2ztLwIe3h77EF3AvXQ1L9V9kmTXJF9Osqxlsf+J7rWjqr4JfIAu+L+m3UR3r9eh2Rr4dV/t8LhZeEmSAJLs127Ovmy8G9GTvCjJj9txft+nqOMy8NXIq6pLgJPo+9i/qlbQfSz/klZf2j9uVVUdA9zOPR/fLwB2AH6VZBldYHxQkoluDBsL5CZ607sBsNNky+8Z++Oe+e7dserKqnpWVT2oqh5Pl3H+wSRzrwkfAi4Bdmk37L2tZ71U1b9X1Z7Ao+lKHt401tQ3z9XANvnDFPrD1tqqJUlr3drc1aGVK34Q2J/uU84XJnlUX7dfAk+pqsfQfVp7wurmNfDVyEmyW5I3JNm2nW8HvBC4oL9vVV0PfBR4xyRTHg28OckmdL80O9Hd7Daf7sa4r9BKGZI8PskjkqyXZEvg34Fzq+rm1v6KJA9uPz8KeCvwjXb+sCRPTLJRkk1aycI84Ly2js8DWyR5WZL1kzyXbseK89r4RybZrI1/MV1ZxHun8Hpt1J5bgA3btaf6u78Z3Q1uK5LsBvx1z7x7t9djQ+A2ujcQYzflXQPs2DPP9+jKTv42yQZJngM8boprkCTNUGtxH9/H0d0AfXlV/Q44HXh2b4eqOr+qbmynF9CVDE7KwFej6Fbg8cD3k9xG9z/7xUx8M9qxwDOTPGaC9q8ANwKHthKEZWMHsAK4varGaml3BL7W1nAx3RZqL+yZ64nAT9q6zmrH21rbZnQZ1BuBXwP7Afu34JyqugH4P3RZ5puBw4FnV9XyNv4ZwOVt/KuB/XrWNZlz6MoonkAX2K+kuylt7GOiJZOMfSPdTg230tVTf6qnbfP22I10ZQvXc8/2aycCj2r7EX+h/aX1HGBh6/98uq3baOt4WJIVU3gukqQZ5D4EvvOSLOo5Xtk39TZ0pX1jlrbHJvJy4KurW683t2nkVNWvgedN0HxSO3r7L6XbbmtM+tqL7qP68a51VN/5J4FPjte3tR8ySdsSYKLge6zPd4A/maDtWO6pY56yqlowSdupdLXCY+cL+9q/TXfzWq93tLZvMMHzafXN8/seWwT8vuSkbZ821vYrYM7Wj9rV/YMlaXZYXlV7TdI+XmJ43H8j2p72Lwf2Xd1FDXwlSZI0kLW4NdlSur3tx2wL/Ka/U/s096P0fII6GUsdpBHXyhVWjHNMVsIgSdJ9MmiZwxRj5QuBXZI8vO0F/wK6nZruuX7yMLqyuZdU1c+nMqkZX2nE9ZcrjIr+sgpJ0uhJ1k6FWlXdleQw4GxgfeBjVbUkyatb+/F0pXdbAse1TYPuWk35hIGvJEmSZp6qGrtJvPex43t+fgXdlz5N2X0KfBcvhgUL7ssMo2vxsmOZf/jrhr0MSZKkoRm1b24bOPA9+OA1uQxJkiSNlCl+GcVMMnDg+8pXdsdsteAks72SJEmjxBpfSZIkTVsYve3BDHwlSZI0kFlT6iBJkqTZbcTi3pHLUEuSJEkDMeMrSZKkgVjqIEmSpFlhxOJeA19JkiRNX4D1RizytcZXkiRJs4IZX0mSJA1kxBK+Br6SJEkaRJHUsBcxLQa+kiRJGogZX0ka0LLLlvOeA08c9jLWiDd+/uXDXsIas8n6o5XRkaSJGPhKkiRp2hL38ZUkSdIsMWJxr4GvJEmSBjNq++KO2nolSZKkgZjxlSRJ0kCs8ZUkSdIsEEatytfAV5IkSdPWhb2jFfha4ytJkqRZwYyvJEmSBpKMVg7VwFeSJEkDGq1SBwNfSZIkDSDW+EqSJEkzkRlfSZIkDWi0Mr4GvpIkSRrIqN3cNlqrlSRJkgZkxleSJEkDstRBkiRJ67iM4K4OBr6SJEkayKgFvtb4SpIkaVYw8JX0e0m2SvLtJLcmOWbY65EkzXTrDXgMh4GvRlKSfZOcn+TmJDckOS/J3kkWJlmVZEWSW5JclORZPeN2SFKtfUWSK5IcPsE1dklye5JPTDJ+RZIjetpfl+Tydu3fJPm3JBv0tF+RZGXP2HP6rvmgJKcluSnJjUlO7Wv/8yT/neS2JFcled4UXqsTklya5O4kC1fT/ZXAcmDzqnrD6uae5JonJfmHQcdLkkZD8v+3d+/BVpV1GMe/jygypYRIjSlKmSjSRRIUUktmahScMRtxRsW8MClTjWNlmmZOf9SMWZNmlqkUChqmZWmG5CW7zOQt0cxLXqJG8wYpiXoUVPDpj/Weabk7Z5/jDthnu5/PmTWw1vuud/32ZvbZv/3u33pRS1u7pMY3Oo6kkcBi4DPAT4HhwIeBl0uXW23vo2pxweOAyyWNtb2qNswo22slTQH+IOlO2zc2XOo84I5+whhle20fx38FLLC9StJo4ErgBODsWp8Dbf+mn3F/Ua45DngJeF/tcU8ELgOOBm4E3gaM6mecur8AVwDfHETfccBfbXsQfSMioquJrOrQRe5efjfTF0xvdxjdaGcA2z8p+6uBGwAk7d7byfZrki4FLgDG00cSa3uppPuBSVTJJGWcw4BVwC3AToMNzPbfa7sCXhvs+ZL2A7YHptteVw7/udbldOBC278u+yvLNlBM55Xx1wxw/QXAEYAlfR74BPA88F1gV6rn+efAibZfUfWR/exyzubAo8BsYK+GcX5n+0BJHwTmU/1bLAEMLLN9+kCPISIiYn1IqUOLZr9/NpO2mdTuMLrVw8A6SQslzZS0VV+dJA0D5gCvUiVlffWZRjWruqx2bCTwNaDZV/2PSnpc0sWSxjSMOVvS81QlA7sBFzacu0jS05JukLRb7fg04CFgoaSVku6QtG9DO5LulfSUpB+XWeX1wvYxwCLgW7a3KLPS64AvAGOADwEfBT5bTtkP+AjVB5FRwKHAStvzGsY5UNJw4GrgUmA08DNg1vqKPSIi2kMt/rRLZnxbNHfyXOZOntvuMN60NKf/F4Xt5yXtA5wC/BDYRtISqrIGgGmSVgFvBdYCn7T9r4ZhnpG0OTACOIsqKev1dWC+7cf6qEN6BtgDuBvYmqocYhGwfy2+y4DLJI0HjgJW1M4/AriLajb4c8D1kiaUMoyxVMnksVQJ+yzgl5J2sv1MaT+y9HkSWAh8r4y5Qdi+s7b7iKQLgX2Bc6g+UGwJTAD+ZPuBJkNNAzYDzillFFdKOnEDhR0RERtNZ82hdla0EYXtB2wfY3ss1YzttlTJGMBttkcBWwHXUNX/NhoDbAGcBEynSsqQNAn4GPCdfq7bY3up7bW2VwDHA/uVWeLGvn8D7gd+UDt2s+3Vtl+y/Q2qcore+FYDj9ieb/tV25cDjwF719ovtv2w7R7gDOCAAZ+s/4OknSUtlrS8zGKfQfXcYfu3wPepkv8V5Sa6/3keim2BJxpqh/uchY+IiM7RaTO+SXyj49l+EFhA7UawcryH6mv5I0t9aeN562yfBazhv1/fTwfeBfxT0nKqxHiWpLv6u3z5s79X8abAe5qFXzv3ntp4fRmofUM4H3gQGG97JHAatcdq+1zbk4H3UpU8nNzb1DDOU8B2ev0U+g4bLOqIiIg+JPGNjiNpgqQvShpb9rcHDgdua+xreyXwI+CrTYY8E/iSpBHAPKpEdVLZLgCupZQySJoqaRdJm0jaGjgX+L3t50r7sZLeUf4+EfgycFPZ30HS3pKGSxoh6WSq2dObSxxXAVtJOlrSMEmHANvV2i8G5kjaUdJbqEo9Fg/i+RpeHpuAzcq1B/va35LqBrceSROoVtLoHXeP8nxsBrxI9QGi96a8FcCOtXFupSo7OUHSppIOBvYcZAwRETEEtbqUWTuXM0viG53oBWAqcLukF6kS3vvo/2a0c4ADJH2gn/ZrgWeB40oJwvLeDegB1th+uvTdEbiuxHAf1RJqh9fG2hu4t8S1pGynlbYtqWZQnwWeAGYAM0tyju1/Ax+nmmV+DjgVOKjU92L7IuAS4HaqMoGXqZZKG8gNVGUSe1El9qupbkpD0hFlVYv+nES1UsMLVPXUV9TaRpZjz5Z4VgLfLm3zgYmq1iO+2vYrwMHAMaX/oVRLt1Hi2EFSzyAeS0REDClqcWsPNVuuc8qUKV66dOlGDCeiUtbVndLuOGLDKcunPV5fzmyT4aO9+Tb7939SBznpqk+1O4T1ZsSwN8+yzl+ZNK7dIUR0DGmXpu/Fk6fs6ltuX9jS2CM2ndqW9/nM+EZEREREV0jiG9HhSrlCTx9bsxKGiIiI9aCzSh2yjm9Eh7O9iGot4Y5S/sOMiIjoWO29Ua0VSXwjIiIiokWdlfim1CEiIiIiukJmfCMiIiKiJeqwOdQkvhERERHRos4qdUjiGxEREREtEOqwxLez5qcjIiIi1LP+1QAAAcJJREFUIlqUGd+IiIiIeMMEWc4sIiIiIrpFZxUPdFa0EREREREtyoxvRERERLSk025uS+IbERERES0QWc4sIiIiIrpCp93clhrfiIiIiBhyJM2Q9JCkZZJO7aNdks4t7fdI2n2gMZP4RkRERESLNmlxa07SMOA8YCYwEThc0sSGbjOB8WWbC5w/mGgjIiIiIt4wtfgzCHsCy2z/w/YrwOXAQQ19DgIuceU2YJSkdzaN13b/jdLTwKODiS5iPRtn++3tDiI2rvzOiYgYUpq+F0u6DhjT4tgjgDW1/Xm259XGPgSYYfvYsn8kMNX28bU+i4Ezbf+x7N8EnGJ7aX8XbXpzWxKPiNiY8jsnIqJz2J6xAYfva1q4cbZ2MH1eJ6UOERERETHUPA5sX9sfCzzZQp/XSeIbEREREUPNHcB4Se+WNBw4DLimoc81wFFldYdpwHO2n2o2aNbxjYiIiIghxfZaSccD1wPDgIts3y/p06X9AmAJcACwDHgJmDPQuE1vbouIiIiIeLNIqUNEREREdIUkvhERERHRFZL4RkRERERXSOIbEREREV0hiW9EREREdIUkvhERERHRFZL4RkRERERX+A9JrtF27nuhGgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 792x576 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "f, reordered_labels, reordered_matrix = fig.plot_composite_matrix(matrix, labels)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "reordered matrix:\n",
      " [[1.    0.382 0.364 0.    0.    0.    0.   ]\n",
      " [0.382 1.    0.386 0.    0.    0.    0.   ]\n",
      " [0.364 0.386 1.    0.    0.    0.    0.   ]\n",
      " [0.    0.    0.    1.    0.356 0.078 0.086]\n",
      " [0.    0.    0.    0.356 1.    0.072 0.078]\n",
      " [0.    0.    0.    0.078 0.072 1.    0.074]\n",
      " [0.    0.    0.    0.086 0.078 0.074 1.   ]]\n",
      "reordered labels: ['SRR2255622_1.fastq.gz', 'SRR2241509_1.fastq.gz', 'SRR2060939_2.fastq.gz', 'SRR2060939_1.fastq.gz', 'SRR453570_1.fastq.gz', 'SRR453569_1.fastq.gz', 'SRR453566_1.fastq.gz']\n"
     ]
    }
   ],
   "source": [
    "print('reordered matrix:\\n', reordered_matrix)\n",
    "print('reordered labels:', reordered_labels)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Customizing plots"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you want to customize the plots, please see the code for `plot_composite_matrix` in [sourmash/fig.py](https://github.com/dib-lab/sourmash/blob/master/sourmash/fig.py), which is reproduced below; you can modify the code in place to (for example) [use custom dendrogram colors](https://stackoverflow.com/questions/38153829/custom-cluster-colors-of-scipy-dendrogram-in-python-link-color-func)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "import scipy.cluster.hierarchy as sch\n",
    "\n",
    "def plot_composite_matrix(D, labeltext, show_labels=True, show_indices=True,\n",
    "                          vmax=1.0, vmin=0.0, force=False):\n",
    "    \"\"\"Build a composite plot showing dendrogram + distance matrix/heatmap.\n",
    "    Returns a matplotlib figure.\"\"\"\n",
    "    if D.max() > 1.0 or D.min() < 0.0:\n",
    "        error('This matrix doesn\\'t look like a distance matrix - min value {}, max value {}', D.min(), D.max())\n",
    "        if not force:\n",
    "            raise ValueError(\"not a distance matrix\")\n",
    "        else:\n",
    "            notify('force is set; scaling to [0, 1]')\n",
    "            D -= D.min()\n",
    "            D /= D.max()\n",
    "\n",
    "    if show_labels:\n",
    "        show_indices = True\n",
    "\n",
    "    fig = pylab.figure(figsize=(11, 8))\n",
    "    ax1 = fig.add_axes([0.09, 0.1, 0.2, 0.6])\n",
    "\n",
    "    # plot dendrogram\n",
    "    Y = sch.linkage(D, method='single')  # centroid\n",
    "\n",
    "    dendrolabels = labeltext\n",
    "    if not show_labels:\n",
    "        dendrolabels = [str(i) for i in range(len(labeltext))]\n",
    "\n",
    "    Z1 = sch.dendrogram(Y, orientation='left', labels=dendrolabels,\n",
    "                        no_labels=not show_indices)\n",
    "    ax1.set_xticks([])\n",
    "\n",
    "    xstart = 0.45\n",
    "    width = 0.45\n",
    "    if not show_labels:\n",
    "        xstart = 0.315\n",
    "    scale_xstart = xstart + width + 0.01\n",
    "\n",
    "    # plot matrix\n",
    "    axmatrix = fig.add_axes([xstart, 0.1, width, 0.6])\n",
    "\n",
    "    # (this reorders D by the clustering in Z1)\n",
    "    idx1 = Z1['leaves']\n",
    "    D = D[idx1, :]\n",
    "    D = D[:, idx1]\n",
    "\n",
    "    # show matrix\n",
    "    im = axmatrix.matshow(D, aspect='auto', origin='lower',\n",
    "                          cmap=pylab.cm.YlGnBu, vmin=vmin, vmax=vmax)\n",
    "    axmatrix.set_xticks([])\n",
    "    axmatrix.set_yticks([])\n",
    "\n",
    "    # Plot colorbar.\n",
    "    axcolor = fig.add_axes([scale_xstart, 0.1, 0.02, 0.6])\n",
    "    pylab.colorbar(im, cax=axcolor)\n",
    "\n",
    "    return fig"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAr4AAAFxCAYAAACV0auJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzde7ycZXnv/8+Xc0vAIFGUgyAnUamkAspW1HRbFfzJtoj1gKdQRa1lV7dH1CJYe6BWLFZFRFEogohnVBSsSlUQJdagRMEigkQJEM7BgBCu3x/PvWQc1lpZa0iYNVmf9+v1vFjP3Ifnngkrueaa67knVYUkSZK0rltv2AuQJEmS7g8GvpIkSZoVDHwlSZI0Kxj4SpIkaVYw8JUkSdKsYOArSZKkWcHAV5IkSTNOko8luTbJxRO0J8m/J7ksyY+TPHZ1cxr4SpIkaSY6Cdhvkvb9gV3a8UrgQ6ub0MBXkiRJM05VfRu4YZIuzwb+ozoXAHOTPHSyOTeYrHHevHm1ww47THuh0n31wx/+cHlVPWjY69D9a968LWqHHbYZ9jK0Dru77hz2EtagDHsBGsfin1w77CWsMXXnjZP+W7z+Hz20atUdg869BLi956ETquqEaU6zDXBVz/nS9tjVEw2YNPDdYYcdWLRo0TTXIN13Sa4c9hp0/9thh21YtOhzw16G1mG33rl02EtYY9ab/J/w0ZJ1J4ift+Nxw17CGnP7VadP+m9xrbqDjR/yjEHnvr2q9hpo8D3G+x+nJhuwDv3WSJIk6X6TkAy1anYpsF3P+bbAbyYbYI2vJEmSpi1AWG+gYw05E3hp291hH+DmqpqwzAHM+EqSJGkGSvJJYAEwL8lS4EhgQ4CqOh44C3gmcBnwW+CQ1c1p4CtJkqQBrN1Sh6p64WraC/ib6cxp4CtJkqSBDLnGd9oMfCVJkjSQjNiOHKMVpkuSJEkDMuMrSZKkAYRRy6Ea+EqSJGkg1vhKkiRpnZeMXuA7WquVJEmSBmTGV5IkSQPImvwWtvuFga8kSZIGYqmDJEmSNAOZ8ZUkSdIA1u5XFq8NBr6SJEkaiIGvJEmS1nmhu71tlIxWmC5JkiQNyIyvJEmSBjB6Nb6jtVpJQ5XkwCRXJVmR5E+HvR5J0hC1b24b5BgWA1+NpCT7Jjk/yc1JbkhyXpK9kyxMsqoFZrckuSjJs3rG7ZCkWvuKJFckObynfeMkJya5MsmtSX6UZP+e9n2SfL1d87okn07y0J72o5Lc2TP/iiQ79rRfkWRlT9s5fc/rQUlOS3JTkhuTnNrT9p4k/9PWdUmSl/a07Zrki21NNyQ5O8kjpvA67t76Lk9SU3jp3wMcVlVzqupHU+g/0XUryc6DjpckzQwGvtJalmRz4MvA+4EHAtsA7wTuaF2+V1VzgLnAccDpSeb2TTO39XkucESSp7XHNwCuAp4CPAA4AjgjyQ6tfQvgBGAHYHvgVuDjfXN/qgWGY8flfe0H9LQ9va/tc8CyNveD6QLNMbcBB7R1vQx4X5InjD0f4EzgEcBWwA+AL7J6dwJnAC+fQl/aupZMsa8kSTPKjKnxPeEEOO20Ya9CI2JXgKr6ZDtfCZwDkOSxY52q6u4kpwDHA7sAF/ZPVFWLkiwB5gNfr6rbgKN6unw5yS+BPYErquqrveOTfAD4rzXxpJI8HdgOWFBVq9rDv8+qVtWRPd2/n+Q7wP8Czq+qH9AFu2Nz/Rvwd0m2rKrrJ7pmVV0KXLq67GuSjYHrgfWBi5Isq6qdWrb8ULog/Srg7VX1+TZmZ+BEutf2TuAbVfX8JN9u017Usswvr6pPJXnTnns+evIXSZI0g4RRy6HOmNWedhosXjzsVWhE/BxYleTkJPsn2WK8TknWBw6hC7qunKDPPsDuwGUTtG9FF2hPlOV88jhtB7RygyVJ/nqcMae2koRzkuzR8/g+wKXAyUmuT3JhkqdMsK4/AvZezbqWTRb0TkdV3dEy5AB7VNVO7edfAE+iy0K/E/hET+nHu+jekGwBbEuXoaeqntwzz5wW9O4HvHFNrFWSdP8ZtVKHGZPxBZg/H849d9ir0EyQSbYFrKpbkuwLvAX4CPCQJGfRZR4B9klyE7ApcBfw4qq6tm+a5S2LuQlwDPCFe68hGwKnAidX1SXjtD8GeAfw7J6Hz6ArhbgGeDzw2SQ39WSnXwT8N93b5NcCZyfZrapuogsOnw68gi5gPwj4YpKdq2p53+WPBy4Czh5nXdsCHwRe39+2plXVp3tOP5XkrcDj6Mos7qQrjdi6qpYC351kqufRlYy8ZW2tVZK0ZsVdHaT7R1X9rKoWVtW2dBnbrYFjW/MFVTWXLtN4Jl1Gst88YA5dlnEBsGFvY7rf5FOA3wGH9Q9uH+N/FXhtVX2nZ10/rarfVNWqqjofeB9dHfFY+3lVtbKqfltV/wzc1LO+lXTlFCdW1Z1VdTpd+cAT+679r+05P6+qqq/tQXRZ1uN6gu21JslLkyxuN+Pd1NY1rzW/mS7A/0HLfv/VJFNtTfdcJUlaawx8NfJaNvYkuqCr9/EVwGuAl4y39VYLTo8Bbm/9AEgSutrUrYCDqurO3nFJtgf+E3hXVZ2yuuXBpF9r09v+43Y+oSTvBPYHnl5Vt/S1bUEX9J5ZVf+4mnXdZ+11+AjdG4Mt25uNi2nPp6qWVdWhVbU18CrguElqia+mq2+WJI2QsN5Ax7AY+GrkJNktyRvaR/ok2Q54IXBBf99W4/pRupKEiRwNvDnJJu38Q8Aj6XZfWNl37W2AbwIfrKrjx1nbs5Nskc7jgL+l7a6Q5GFJnphkoySbJHkTXXb0vDb888AWSV6WZP0kz6XbseK8Nv6twMHA0/prd9tOF2cD51XV4UxRW+cmwEbtfJNWAjIVm9IF6te1sYfQ8+YjyV+O/RkBN7a+YzftXQPseM9UnAEsnOq6JUkzw6jV+Br4ahTdSlc/+/0kt9EFvBcDb5ig/7HAM1tN7ni+QheYHdqymK+i24lgWe7Zb/dFre8r6AK2I3vaVvTM9QK6G+VuBf4D+JeqOrm1bUYXVN8I/BrYD9h/LIitqhuA/0NXfnEzcDjw7J763n8CHgb8T8+139baDqS72e2Q/OEewg+b+GUEuhrcldxzk9xKuhvsAEjy1Z5r/IGq+ildffT36ALZP+GeIJ62nu+31+dMurKQX7a2o+hu4rspyfPabhnHIkkaHQkZ8BjakvtKBP/AXnvtVYsWLbpfFrJgQfdfb24TQJIfVtVew16H7l977bV7LVr0uWEvQ+uwW+9cOuwlrDHrzaz70++bIQZCa9q8HY8b9hLWmNuvOn3Sf4s3+uNtaqudXzXQ3Et/cuRQ/p1fh35rJEmSdH9yVwdJM0orV1gxzjFuCYMkSVMRRu/mNjO+0jquqvYf9hokSesi9/GVJEmSZiQzvpIkSRrIqGV8DXwlSZI0gAy1XncQBr6SJEkazIhlfEdrtZIkSdKAzPhKkiRp+mKNryRJkmaBMNyvHx7EaIXpkiRJ0oDM+EqSJGkg7uogSZKkWcEaX0mSJM0CAWt8JUmSpJnHjK8kSZKmL4xcCtXAV5IkSYOx1GH6TjgBFi8e9iokSZI0Lclgx5DMiIzvaafBzTcPexWSpHXdejPjn701Y8QybZOZt+Nxw17CGrP88tcMewlrzJwNTx/2Eta4GfM3wAMeMOwVSJIkaVpmRO3A1M2YwFeSJEkjJFAj9smDga8kSZIGM1px76glqCVJkqTBmPGVJEnSYNYbrZSvga8kSZIGMHpfWWzgK0mSpOkL1vhKkiRJM5EZX0mSJA3GGl9JkiTNCtb4SpIkaVYYrbjXGl9JkiTNDmZ8JUmSNH3BGl9JkiTNEqMV9xr4SpIkaTA1Yje3WeMrSZKkGSfJfkkuTXJZksPHaX9Aki8luSjJkiSHrG5OM76SJEmavmSt1fgmWR/4IPA0YClwYZIzq+qnPd3+BvhpVR2Q5EHApUlOrarfTTSvGV9JkiQNJgMeq/c44LKqurwFsqcDz+7rU8BmSQLMAW4A7ppsUgNfSZIk3d/mJVnUc7yyr30b4Kqe86XtsV4fAB4J/Ab4CfDaqrp7sota6iBpypIcCPw7sAXwpKr60ZCXJEkapsFvblteVXtNNvM4j1Xf+TOAxcD/BnYCvp7kO1V1y0STmvHVSEqyb5Lzk9yc5IYk5yXZO8nCJKuSrEhySyt4f1bPuB2SVGtfkeSK3oL5JBsnOTHJlUluTfKjJPv3tO+T5Ovtmtcl+XSSh46zvo2SXJJkad/j70rykyR3JTmqr21Bkrt71rYiycv61vax9ryWJXl93/gDklzcxp2f5FFTeB13T3J2kuVJ+v9CGc97gMOqas59CXrbn8HOg46XJM0Q62WwY/WWAtv1nG9Ll9ntdQjwuepcBvwS2G3S5U7jqUkzQpLNgS8D7wceSPfRxzuBO1qX71XVHGAucBxwepK5fdPMbX2eCxyR5Gnt8Q3oPlp5CvAA4AjgjCQ7tPYtgBOAHYDtgVuBj4+zzDcB147z+GXAm4GvTPD0ftOCyrHj5J62o4Bd2nX/DHhzkv3aa7ILcCrw6va8vwScmWR1n+rcCZwBvHw1/cZsDyyZYl9J0rps0PreqSWJLwR2SfLwJBsBLwDO7OvzK+CpAEm2Ah4BXD7ZpDOq1GHxYliwYNir0AjYFaCqPtnOVwLnACR57Finqro7ySnA8XQB44X9E1XVoiRLgPnA16vqNroAc8yXk/wS2BO4oqq+2js+yQeA/+p77OHAi4HXAx/pu97Jrc+LpveUAXgpcEhV3QjcmOQjwELga3Qf93ynqr7b5v8X4B10Afw3Jpqwqi6luwt20uxrko2B64H1gYuSLKuqnVq2/FDgwXRvGN5eVZ9vY3YGTqR7be8EvlFVz0/y7TbtRS3L/PKq+lSSN+2556On/6pIktY5VXVXksOAs+n+7flYVS1J8urWfjzwLuCkJD+hC6ffUlXLJ5t3xgS+W20FD73XB8bSuH4OrEpyMt1dnhe0YPAPtK1QDqELuq4cb6Ik+wC7A/88QftWdIH2RFnOJ4/T9n7gbXQB+XQ9OMk1wG+BLwB/V1W3JdkC2Bq4qKfvRcBfjC2VP3wPPXa+O5MEvlNVVXcAc1qgukf7SAngF8CTgGXAXwKfSLJzVV1N9xfSOXTZ6Y2AvdpcT+6fp2Wu33hf1ylJup+txS+wqKqzgLP6Hju+5+ffAE+fzpwzJvB96EPh3HOHvQrNFJP9HlXVLUn2Bd5Cl1F9SJKz6DKPAPskuQnYlG5bkxdXVX/ZwfKWxdwEOIYuyOxbQzakKx84uaouGaf9MXRZ1Wf3PHYgsEFVfT7Jgqk929+7hC47egldScHJwHuBV9Ft0wJwc0//m4HN2s9fB45u1zyf7rXZCPjjaa5hWqrq0z2nn0ryVrotaL5I94Zje2DrqloKfHeSqZ5HVzLylrW1VknSmpa1GviuDdb4aiRV1c+qamFVbUuX1dwaOLY1X1BVc+nqcc+ky0j2m0cXTL4RWABs2NuYZD3gFOB3wGH9g9vH+F+l2zrlO+2xTYF3A/93wOe0rKp+WlV3V9Uv6WqBn9uaV7T/bt4zZHO6GmNaYP4yuq1drm7P76d0NwesNUlemmRxkpvam43d27Vp6w/wg3TfqPNXk0y1NX+4bY0kaaYLXSQ5yDEkBr4aeS3oO4ku6Op9fAXwGuAlSf50nHGrquoY4PbWD4AkoatN3Qo4qKru7B2XZHvgP4F3VdUpPU270N309p0ky4DPAQ9tOzDsMMhTo5UvtFKOq4E9etr3oKfMoqo+U1W7V9WWwJF02dZ71TWvKe11+AjdG4Mt25uNi3vWvKyqDq2qremy1sdNUkt8NX94964kSWucga9GTpLdkrwhybbtfDvghcAF/X2r6nrgo3QlCRM5mm6HhE3a+YfoNsQ+oKr+oE43yTbAN4EP9tYZNRfTBW/z2/EK4Jr281Vt/IbtOusBGyTZpNUij21n9rB0tmvr+mLP/P8B/F2SLZLsRlfacVLP2vZMsn66r238MPCl8Uo0+p5P2no2auebtBKQqdiULji/ro09hJ43H0n+cuzPCLix9V3Vzq8BduyZ6wy6G/UkSaMkGewYEgNfjaJbgccD309yG13AezHwhgn6Hws8s9XkjucrdIHZoS2L+Sq6YHVZz366Y7swvIIuYDuyd79d6O5AbVnOZVW1jO6rE+9u52MB30fobnp7IfD29vNLWttjge8Bt9HV6V4M/G3POo+ku5nsSrqdJP61qr7W0/4+4Cbg0vbfQ1m97dsaxjLHK9t4AJJ8NcnbxhvYvi/9mLbma4A/Ac7r6bI33Z/RCrqSk9e2Eg7ods44uZVIPK/tlnEskqTRsva2M1srUjXxnvV77bVXLVq0aK0vYmwLM29u05gkP1zNN7poHbTXXrvXokWfG/YytA677c5lw17CmjNiNxVNZt6Oxw17CWvM8stfs/pOI2LOhk+e9N/ijR+4Y237jH8YaO7LT3/RUP6dN+MrSZKkWcHAV1rHtXKFFeMc45YwSJI0JWHkanxnzD6+ktaOqtp/2GuQJK2jRqzixsBXkiRJg1lvtCJfSx0kSZI0K5jxlSRJ0gBG7yuLDXwlSZI0fUPek3cQBr6SJEkajDW+kiRJ0sxjxleSJEmDGbGMr4GvJEmSpi9QoxX3GvhKkiRpQCOW8bXGV5IkSbOCGV9JkiQNxn18JUmStM5LLHWQJEmSZiIzvpIkSRrMiKVQBw58TzgBTjttzSxi8WKYP3/NzCVJkqT7yWyp8T3tNANWSdJo+eDPbh/2EtaYI5/18WEvYY1Zfvlrhr2ENadq2Cu4/4SRq/G9T6UO8+fDuefe90UsWHDf55AkSZImY42vJEmSBlKzpdRBkiRJs9xsublNkiRJs9gI1viOWJwuSZIkDcaMryRJkgaQ2bOdmSRJkma5ESt1MPCVJEnSYEYr7rXGV5IkSbODGV9JkiRNX6AsdZAkSdKsYOArSZKkWWHEdnWwxleSJEmzghlfSZIkTV8YuRSqga8kSZIG4BdYSJIkabYYsZvbRixBLUmSJA3GjK8kSZKmL4xcxtfAV5IkSQMpa3wlSZI0K4xY0eyILVfSMCU5MMlVSVYk+dNhr0eSpOkw8NVISrJvkvOT3JzkhiTnJdk7ycIkq1pgdkuSi5I8q2fcDkmqta9IckWSw3vaN05yYpIrk9ya5EdJ9u+79lOTXJLkt0m+lWT7vvbHJvl2m/+aJK/tu/632thLkvx5T1uSvD3Jr9raT0+yeU/7u1vQeUtb39v7rntAkovbdc9P8qgpvI4vS/LDNufSdo3JPgl6D3BYVc2pqh+tbv5JrltJdh50vCRpBgjddmaDHENi4KuR04LBLwPvBx4IbAO8E7ijdfleVc0B5gLHAacnmds3zdzW57nAEUme1h7fALgKeArwAOAI4IwkO7RrzwM+1x5/ILAI+FTP2uYBXwM+DGwJ7Ayc03PdTwI/am1vBz6T5EGt7aXAS4AnAlsDf9Se45gTgd2qanPgCcDBSZ7TrrsLcCrw6va8vwScuZogFuCPgdcB84DHA08F3jhJ/+2BJauZU5I0K6S7uW2QY0is8Z2qE06A004b9irU2RWgqj7ZzlfSgsskjx3rVFV3JzkFOB7YBbiwf6KqWpRkCTAf+HpV3QYc1dPly0l+CewJXAE8B1hSVZ9u1zsKWJ5kt6q6BHg9cHZVndrG3wH8rPXdFXgs8PSqWgl8NsnrgIPaGg8ATqyqq1r/fwG+meSvq+q3VXVp3/LvpgusAZ4BfKeqvtsz9h10Afw3Jnohq+pDPae/TnIq8Gf9/ZJsDFwPrA9clGRZVe3UsuWHAg+me8Pw9qr6fBuzM12wPh+4E/hGVT0/ybfbtBclKeDlVfWpJG/ac89HT7RUSZLuMzO+U3XaabB48bBXoc7PgVVJTk6yf5ItxuuUZH3gELqg68oJ+uwD7A5cNkH7VnSB9liW89HARWPtLVD+RXscYB/ghlZqcG2SLyV5WM/Yy6vq1p5LXNQzNu2g53xjuqB9bD2HJ1kBLAU2BU6bZGzac5uOJzNORreq7mgZcoA9qmqn9vMvgCfRZcffCXwiyUNb27vo3pBsAWxLy15X1ZN75pnTgt79mDzTLEmaicz4rsPmz4dzzx32KmaHSep/quqWJPsCbwE+AjwkyVl0mUeAfZLcRBcY3gW8uKqu7ZtmectibgIcA3zh3kvIhnTlAye3bC7AHOC6vq43A5u1n7ely+o+DfgJ8G668oYntrE3jzN2m/bzV4E3JzkDuLE9P+jKEcae+9Etmzsf+Iue+b4OHJ1kAXB+G7tR79jVSXIIsBfwiqmOGct8N59K8lbgccAX6d5wbA9sXVVLge9OMtXzgI9zz3OWJM10/SmXEWDGVyOpqn5WVQuralu6rObWwLGt+YKqmkuXaTyTLiPZbx5dIPpGYAGwYW9jkvWAU4DfAYf1NK0ANucPbQ6MZXFXAp+vqgur6na6LOgTkjxgCmM/Rhckn0uXdf1We3xp33OvdmPZyjY/LTB/GfAB4Or2/H7aP3YiSf4COBrYv6qWT2VMG/fSJIuT3NTebOzerg3wZrq/En+QZEmSv5pkqq3pSiUkSSOigFovAx3DYuCrkdeCvpPo+1i/qlYArwFeknG23qqqVVV1DHB76wd0uyvQ1aZuBRxUVXf2DFsC7NHTd1NgJ+4pD/gx3d8Fv7/MWNfWZ8ckm/W07zE2tqrurqojq2qHFtAvAX7djvFs0K499nw+U1W7V9WWwJF02dZ71TX3a2UGHwEOqKqfrK5/z7jt27jDgC3bm42L23OlqpZV1aFVtTXwKuC4SXZyuBrYbqrXliRpEAa+GjlJdkvyhiTbtvPtgBcCF/T3rarrgY/S3eg1kaPpSgw2aecfAh5JFwiu7Ov7eWD3JAe1/u8AftxTCvFx4MAk81upxBHAd6vqpqr6ObAYODLJJkkOBB4DfLY9jwcm2alta/Yo4L3A37eb9NZL8qokW7T2xwF/Q8+Na0n2TLJ+2yXiw8CXetY10Wv5v+nKOQ6qqh9M1nccm9IF9te1uQ6h581Hkr8c+zOiK90oYFU7vwbYsWeuM4CF07y+JGnY3M5MWutupdt66/tJbqMLeC8G3jBB/2OBZyZ5zATtX6ELzA5tWcxX0dXQLss9+/2+CKCqrqPbheEf25jHAy8Ym6iqvgm8rc15Ld2uCwf3XOsFdHW0N9IF3M9tc0JXInAWcBtdve/HquqEnrEH0t1MdivwCbqbxXq3O3sfcBNwafvvoazeEXQ3pp3V81y/OtaY5KtJ3jbewKr6KV199PfoAtk/Ac7r6bI33Z/RCrqSk9dW1S9b21HAya1E4nlV9VXuKVWRJI2C4M1t0tpWVb+muxlqPCe1o7f/UrrdEcakr724Z2eFe7WPc/3/BHabpP1DdFnj8dquoKspHq/t58AjJmi7G9hvNevad7L2Ccbca+uyvvb9+877X7u30+1HPN7YN9PV+Y7XdjzdFm69jx291167//MUli1Jmim8uU2SJEmaeQx8pXVcK1dYMc4xbgmDJElTEWC99QY7hsVSB2kd11+uIEnSmjLE+9QGYuArSZKk6RvuBg0DsdRBkiRJM06S/ZJcmuSyJIdP0GdB+yKlJUn+a3VzmvGVJEnSAELWUso3yfrAB4Gn0X0L6YVJzmxbaY71mQscB+xXVb9K8uDVzWvGV5IkSdMW1ur3VzwOuKyqLq+q3wGnA8/u63Mw8Lmq+hVAVV27ukkNfCVJkjSQ+xD4zkuyqOd4Zd/U2wBX9ZwvbY/12hXYIsm5SX6Y5KWrW6+lDpIkSbq/La+qvSZpHy8vXH3nGwB7Ak8F/gj4XpIL2hdCjcvAV5IkSdMXyNqrHVgKbNdzvi3wm3H6LK+q24Dbknwb2AOYMPC11EGSJEkDWYs1vhcCuyR5eJKNgBcAZ/b1+SLwpCQbJPlj4PHAzyab1IyvJEmSpi3AemtpH9+quivJYcDZwPrAx6pqSZJXt/bjq+pnSb4G/Bi4G/hoVV082bwGvpIkSZpxquos4Ky+x47vO/9X4F+nOqeBryRJkgYyat/cZuArSZKkgRj4SpIkad0X1to3t60tMybwXbwYFiwY9iomsfhYzp3/umGvQpIkSQOaEYHvwQcPewWSpNngyGd9fNhLWGPe+eVDhr2ENaf6v5dgdN3NXcNewv1qLe7ju1bMiMD3la/sjhltgdleSZKkMWH0anxHLE6XJEmSBjMjMr6SJEkaMVP/FrYZw8BXkiRJAzHwlSRJ0qywtr6yeG2xxleSJEmzghlfSZIkTdso7upg4CtJkqSBGPhKkiRp3RfIiBX5WuMrSZKkWcGMryRJkgZiqYMkSZJmBQNfSZIkrfNGcVcHa3wlSZI0K5jxlSRJ0vRl9L65zcBXkiRJAxm1UgcDX0mSJA0kI1Y0O2LLlSRJkgZjxleSJEnTNoq7Ohj4SpIkafoCGbHI11IHSVOW5MAkVyVZkeRPh70eSdJwJYMdw2Lgq5GUZN8k5ye5OckNSc5LsneShUlWtcDsliQXJXlWz7gdklRrX5HkiiSH97RvnOTEJFcmuTXJj5Ls33ftpya5JMlvk3wryfZ97Y9N8u02/zVJXtt3/W+1sZck+fOetiR5e5JftbWfnmTznvZ3t6Dzlra+t/dd94AkF7frnp/kUVN4HXdPcnaS5UlqCi/9e4DDqmpOVf1oCv0num4l2XnQ8ZIkDcLAVyOnBYNfBt4PPBDYBngncEfr8r2qmgPMBY4DTk8yt2+aua3Pc4EjkjytPb4BcBXwFOABwBHAGUl2aNeeB3yuPf5AYBHwqZ61zQO+BnwY2BLYGTin57qfBH7U2t4OfCbJg1rbS4GXAE8Etgb+qD3HMScCu1XV5sATgIOTPKdddxfgVODV7Xl/CTgzyerKme4EzgBevpp+Y7YHlkyxryRpHTdqGV9rfKdj8WJYsGDYqxDsClBVn2znK2nBZZLHjnWqqruTnAIcD+wCXNg/UVUtSrIEmA98vapuA47q6fLlJL8E9gSuAJ4DLKmqT7frHQUsT7JbVV0CvB44u6pObePvAH7W+u4KPBZ4elWtBD6b5HXAQW2NBwAnVtVVrf+/AN9M8tdV9duqurRv+eWHJzAAACAASURBVHfTBdYAzwC+U1Xf7Rn7DroA/hsTvZBtzktXl31NsjFwPbA+cFGSZVW1U8uWHwo8mO4Nw9ur6vNtzM50wfp8ugD7G1X1/CTfbtNe1LLML6+qTyV50557PnqyZUiSZpgRK/E14ztlBx8M8+cPexXq/BxYleTkJPsn2WK8TknWBw6hC7qunKDPPsDuwGUTtG9FF2iPZTkfDVw01t4C5V+0xwH2AW5opQbXJvlSkof1jL28qm7tucRFPWPTDnrON6YL2sfWc3iSFcBSYFPgtEnGpj23+6yq7mgZcoA9qmqn9vMvgCfRZcffCXwiyUNb27vo3pBsAWxLy15X1ZN75pnTgt79gDeuibVKku4fofvmtkGOYTHjO1WvfGV36P4xyVvIqrolyb7AW4CPAA9JchZd5hFgnyQ30QWGdwEvrqpr+6ZZ3rKYmwDHAF+49xKyIV35wMktmwswB7iur+vNwGbt523psrpPA34CvJuuvOGJbezN44zdpv38VeDNSc4AbmzPD+CPe5770S2bOx/4i575vg4cnWQBcH4bu1Hv2LVhLPPdfCrJW4HHAV+ke8OxPbB1VS0FvjvJVM8DPs49z1mSpDXOjK9GUlX9rKoWVtW2dFnNrYFjW/MFVTWXLtN4Jl1Gst88ukD0jcACYMPexiTrAacAvwMO62laAWzOH9ocGMvirgQ+X1UXVtXtdFnQJyR5wBTGfowuSD6XLsP8rfb40r7nXu3GspVtflpg/jLgA8DV7fn9tH/smpbkpUkWJ7mpvdnYvV0b4M10CYEfJFmS5K8mmWprulIJSdKoGDDbO8yMr4GvRl4L+k6i72P9qloBvAZ4yXhbb1XVqqo6Bri99QO63RXoalO3Ag6qqjt7hi0B9ujpuymwE/eUQvwY6N0dYezntD47Jtmsp32PsbFVdXdVHVlVO7SAfgnw63aMZ4N27bHn85mq2r2qtgSOpMu23quueU1pu1l8hO6NwZbtzcbFtJKLqlpWVYdW1dbAq4DjJqklvhrYbm2tVZK05nWlDjXQMSwGvho5SXZL8oYk27bz7YAXAhf0962q64GP0t3oNZGj6UoMNmnnHwIeCRzQbkLr9Xlg9yQHtf7vAH7cUwrxceDAJPNbqcQRwHer6qaq+jmwGDgyySZJDgQeA3y2PY8HJtmpbWv2KOC9wN+3m/TWS/KqJFu09scBf0PPjWtJ9kyyftsl4sPAl3rWNdFrmfY8Nmrnm7QSkKnYlC6wv66NPYSeNx9J/nLsz4iudKOAVe38GmDHnrnOABZO8bqSJA3EwFej6Fbg8cD3k9xGF/BeDLxhgv7HAs9M8pgJ2r9CF5gd2rKYr6KroV2We/b7fRFAVV1HtwvDP7YxjwdeMDZRVX0TeFub81q6XRcO7rnWC4C92tijgee2OaErETgLuI2u3vdjVXVCz9gD6W4muxX4BN3NYr3bnb0PuAm4tP33UFZve7qSibGM9co2HoAkX03ytvEGVtVP6eqjv0cXyP4JcF5Pl73p/oxW0JWcvLaqftnajgJObiUSz6uqr3JPqYokaUSMWqmDN7dp5FTVr+luhhrPSe3o7b+UbneEMelrL+7ZWeFe7eNc/z+B3SZp/xBd1ni8tivoaorHa/s58IgJ2u4G9lvNuvadrH2S9Uz4fKtq/77z/tfu7XT7EY839s10db7jtR1Pt4Vb72NH77XX7v88pYVLkmaEUcugGvhKkiRp2sZqfEfJqAXqkqaplSusGOcYt4RBkqR1lRlfaR3XX64gSdKaMsx63UEY+EqSJGnawuiVDhj4SpIkafqGvEPDIEYtUJckSZIGYsZXkiRJA8mI7epg4CtJkqRp67YzG/YqpsfAV5IkSQMZtZrZUVuvJEmSNBAzvpIkSZq2UCP3zW0GvpIkSRqINb6SJEmaFUatZnbU1itJkiQNxIyvJEmSpi0j+M1tBr6SJEkaiDe3SZIkaZ03il9gYY2vJEmSZgUzvpJmjKt/+zv+cfGVw17GGvHbu0YsDTKJ9z7nxGEvYY1Z+at3DnsJ0jpl1DKoBr6SJEkaiDW+kiRJWudZ4ytJkiTNUGZ8JUmSNH3u4ytJkqTZIIxe6cCorVeSJEkaiIGvJEmSBrJeaqBjKpLsl+TSJJclOXySfnsnWZXkuaub01IHSZIkDWRt1fgmWR/4IPA0YClwYZIzq+qn4/T7F+DsqcxrxleSJEnTNlbjO8gxBY8DLquqy6vqd8DpwLPH6fd/gc8C105lUgNfSZIk3d/mJVnUc7yyr30b4Kqe86Xtsd9Lsg1wIHD8VC9qqYMkSZIGch9KHZZX1V6TtI83c39x8LHAW6pqVTK1hRj4SpIkaSBZe19ZvBTYrud8W+A3fX32Ak5vQe884JlJ7qqqL0w0qYGvJEmSpi1r9wssLgR2SfJw4NfAC4CDeztU1cPvWUtOAr48WdALBr6SJEmaYarqriSH0e3WsD7wsapakuTVrX3Kdb29DHwlSZI0kLW5S0JVnQWc1ffYuAFvVS2cypwGvpIkSZq2MPUvo5gpDHwlSZI0kLVY47tWuI+vJEmSZgUzvpIkSRqIGV9JIyvJVkm+neTWJMcMez2SpJkrdNstDHIMi4GvRlKSfZOcn+TmJDckOS/J3kkWJlmVZEWSW5JclORZPeN2SFKtfUWSK5IcPsE1dklye5JPTDJ+RZIjetpfl+Tydu3fJPm3JBv0tF+RZGXP2HN62t7WN+/KJHcnmdfaN07ysTb3siSvn+JrdUKSS9tcC1fT/ZXAcmDzqnrDVOaf4JonJfmHQcdLkkbDeqmBjqGtd2hXlgaUZHPgy8D7gQfSfXf3O4E7WpfvVdUcYC5wHN23usztm2Zu6/Nc4IgkTxvnUh+k20B7PHOrak473tXz+JeAx1bV5sDuwB7A3/aNPaBn7NPHHqyqf+p5fA7wL8C5VbW8dTkK2AXYHvgz4M1J9ptgfb0uAl4D/PcU+m4P/LSqRus2XUmSpsAa3wGd8MMTOO0npw17GbPVrgBV9cl2vhI4ByDJY8c6VdXdSU4BjqcLGO8VxFbVoiRLgPnA18ceT/IC4CbgfGDnqS6sqn7Rcxrg7umM77l+gJcAf9/z8EuBQ6rqRuDGJB8BFgJfW82aPtjmvH011zwJeBFQSV4H/AVwC/A+4JF0r/NngddX1e/aGt/bxmwMXEn3rTpP6JvnW1V1QJI/BU6k+7M4i+471y+rqr+bymsiSZpZ1vI3t60VZnwHdNpPTmPxssXDXsZs9XNgVZKTk+yfZIvxOiVZHzgEuJMuKBuvzz50mdnLeh7bnC7gnOyj/iuTLE3y8bFShJ7xBye5ha5kYA/gw31jT01yXZJzkuwxwfxPAraiCzRpz3FruuztmIuAR0+yxmlpm3+fCry7ZZ3/E1gF/D+670D/X8BT6bLHAE8Hnkz3RmQu8Hzg+qo6oW+eA5JsBHwBOIUuS/9p4KA1tXZJ0nCsl8GOYTHjex/Mf8h8zl147rCXsU7KIRP/VlTVLUn2Bd4CfAR4SJKzgENbl32S3ARsCtwFvLiqru2bZnmSjYFNgGPogrIx7wJOrKqruqTmH44D9gYWA1vSlUOcCjyjZ32nAacl2YUuS3tNz/gX0ZUcBHgtcHaS3arqpr7rvAz4TFWtaOdz2n9v7ulzM7BZ/wLXpKr6Yc/pFUk+DDwFOJbuDcVmwG7AD6rqZ5NMtQ+wIXBsK6P4zFRrlCVJM1OA9c34SmtfVf2sqhZW1bZ0Gdut6YIxgAuqai6wBXAmXfa03zy6YPKNwAK6oIwk84E/B/5tguuuqKpFVXVXVV0DHAY8vWWJ+/v+D7CErs547LHzqmplVf22qv6ZrpziD9aX5I+AvwRO7nl4LADuvc7mwK3jrXNNSbJrki+3m+luAf6J7rWjqr4JfIAu+L+m3UR3r9eh2Rr4dV/t8LhZeEmS1hYDX428qroEOIkuAO59fAXdx/IvafWl/eNWVdUxwO3c8/H9AmAH4FdJltEFxgclmejGsLFAbqL3vBsAO022/HHGPge4ATi3Z603AlfTlU6M2YMusF6bPgRcAuzSbth7Gz3rrap/r6o96UoudgXeNNbUN8/VwDb5wxT6w9baqiVJ94tRK3Uw8NXISbJbkjck2badbwe8ELigv29VXQ98FHjHJFMeTbdDwibACXSB6vx2HA98hVbKkOTxSR6RZL0kWwL/Trfzws2t/RVJHtx+fhTwVuAb7fxhSZ6YZKMkmyR5E1329Ly+9bwM+I9xdlb4D+DvkmyRZDe60o6TpvB6bdSeW4AN27Wn+ru/Gd0NbivaNf+6Z9692+uxIXAb3RuIVa35GmDHnnm+R1d28rdJNkjyHOBxU1yDJGkGCm5nJt0fbgUeD3w/yW10Ae/FTHwz2rHAM5M8ZoL2rwA3Aoe2EoRlYwddicHtVXVd67sj3S4Kt7Zr3kEXdI95IvCTtq6z2vG21rYZXQb1RuDXwH7A/i04ByDJNsD/pgty+x0J/IKuROC/gH+tqkl3dGjOoduR4Ql0gf1KupvSSPKitqvFRN5It1PDrXT11J/qadu8PXZjW9P1wHta24nAo5LclOQLVfU7ukz2wtb/+cDnep73w5KsQJI0OgbM9npzmzQNVfVr4HkTNJ9EXxa0qpbSbbc1Jn3txQS7I1TVUX3nnwQ+OV7f1n7IJG1LgImC77E+v2aC38uqugP4q3ZMWVUtmKTtVLqb88bOF/a1f5vu5rVe72ht32CC59Pqm+f3PbYI+H3JSds+baztV8CcbR61q/sHS5LWGgNfSZIkTdvYVxaPEksdpBHXyhVWjHOs7RvfJEmznKUOku5X/eUKo6K/rEKSpLXNwFeSJEkDGeYODYMw8JUkSdK0jeI3txn4SpIkaSDDrNcdhDe3SZIkaVYw4ytJkqRpy5B3aBiEga8kSZIGYuArSZKkdV53c9to7epgja8kSZJmBTO+kiRJGsioZVANfCVJkjRtwRpfSZIkzRKjFviOWoZakiRJGogZX0mSJE1bUiO3q4OBryRJkgYyaqUOBr6SJEmaNm9um0UWL1s87CVI65wCbl81Yn+LTuC9zzlx2EtYY17/uZcPewmStEYY+EqSJGkgZnwlSZK0zuu+snjYq5geA19JkiRNX2C9EdvVwX18JUmSNCuY8ZUkSdJARi2DauArSZKkaXM7M0mSJM0ao3Zz26hlqCVJkqSBmPGVJEnStHWlDqO1q4OBryRJkgYyajW+ljpIkiRpVjDjK0mSpGlLRi/ja+ArSZKkgYxa6YCBryRJkgaSEcv4jlqgLkmSJA3EjK8kSZIGMmIJXwNfSZIkTV8YvVIHA19JkiQNZNRqZkdtvZIkSdJAzPhKkiRpIBmxryw24yvp95JsleTbSW5Ncsyw1yNJmtky4DEsBr4aSUn2TXJ+kpuT3JDkvCR7J1mYZFWSFUluSXJRkmf1jNshSbX2FUmuSHL4BNfYJcntST4xyfgVSY7oaX9dksvbtX+T5N+SbNDTfkWSlT1jz+m75oOSnJbkpiQ3Jjm1p22bJF9sz3dpkldP8bU6IcmlSe5OsnA13V8JLAc2r6o3TGX+Ca55UpJ/GHS8JGnmG7u5bZBjWAx8NXKSbA58GXg/8EBgG+CdwB2ty/eqag4wFzgOOD3J3L5p5rY+zwWOSPK0cS71QeDCCZYxt6rmtONdPY9/CXhsVW0O7A7sAfxt39gDesY+va/tc8AyYHvgwcB7eto+AfwS2Ar4/4B/SvJnE6yv10XAa4D/nkLf7YGfVtVofXYlSdIUWON7H6z43QoWnLRg2MuYjXYFqKpPtvOVwDkASR471qmq7k5yCnA8sAvjBLFVtSjJEmA+8PWxx5O8ALgJOB/YeaoLq6pf9JwGuHuq45M8HdgOWFBVq9rDP2ptc4AFwPOq6k7goiSfAf4K+NZq1vTBNsftq7n+ScCLgEryOuAvgFuA9wGPpHudPwu8vqp+lyTAe9uYjYErgYOBJ/TN862qOiDJnwIn0v1ZnAUUcFlV/d1UXh9J0swzYruZmfEd1FabbsWcjeYMexmz1c+BVUlOTrJ/ki3G65RkfeAQ4E66oGy8PvvQZWYv63lsc+Dvgck+6r+ylRt8PMm8vjkPTnILXcnAHsCH+8aemuS6JOck2aPn8X2AS4GTk1yf5MIkTxmbtu+/Yz/vPskap6WqFgKnAu9u2ej/BFYB/w+YB/wv4Kl02WOApwNPpnsjMhd4PnB9VZ3QN88BSTYCvgCcQpel/zRw0JpauyRpCALrDXgMixnfAT10s4fy0M0eyrkLzx32UtZJOWTi34qquiXJvsBbgI8AD0lyFnBo67JPkpuATYG7gBdX1bV90yxPsjGwCXAMXVA25l3AiVV1Ve5diLQc2BtYDGxJVw5xKvCMnvWdBpyWZBfgpcA1PeNfRFdyEOC1wNlJdquqm4Bt6YLJV9AF7AcBX0yyc1UtT3IeXVnGm4BHtfbrJnyh1oCq+mHP6RVJPgw8BTiW7g3FZsBuwA+q6meTTLUPsCFwbCuj+EyS16+lZUuS7idmfKX7QVX9rKoWVtW2dFnPremCMYALqmousAVwJvCkcaaYB8wB3khXQrAhQJL5wJ8D/zbBdVdU1aKququqrgEOA57essT9ff8HWEJXZzz22HlVtbKqfltV/0xXTjG2vpXAFVV1YlXdWVWnA1cBT2ztLwIe3h77EF3AvXQ1L9V9kmTXJF9Osqxlsf+J7rWjqr4JfIAu+L+m3UR3r9eh2Rr4dV/t8LhZeEmSAJLs127Ovmy8G9GTvCjJj9txft+nqOMy8NXIq6pLgJPo+9i/qlbQfSz/klZf2j9uVVUdA9zOPR/fLwB2AH6VZBldYHxQkoluDBsL5CZ607sBsNNky+8Z++Oe+e7dserKqnpWVT2oqh5Pl3H+wSRzrwkfAi4Bdmk37L2tZ71U1b9X1Z7Ao+lKHt401tQ3z9XANvnDFPrD1tqqJUlr3drc1aGVK34Q2J/uU84XJnlUX7dfAk+pqsfQfVp7wurmNfDVyEmyW5I3JNm2nW8HvBC4oL9vVV0PfBR4xyRTHg28OckmdL80O9Hd7Daf7sa4r9BKGZI8PskjkqyXZEvg34Fzq+rm1v6KJA9uPz8KeCvwjXb+sCRPTLJRkk1aycI84Ly2js8DWyR5WZL1kzyXbseK89r4RybZrI1/MV1ZxHun8Hpt1J5bgA3btaf6u78Z3Q1uK5LsBvx1z7x7t9djQ+A2ujcQYzflXQPs2DPP9+jKTv42yQZJngM8boprkCTNUGtxH9/H0d0AfXlV/Q44HXh2b4eqOr+qbmynF9CVDE7KwFej6Fbg8cD3k9xG9z/7xUx8M9qxwDOTPGaC9q8ANwKHthKEZWMHsAK4varGaml3BL7W1nAx3RZqL+yZ64nAT9q6zmrH21rbZnQZ1BuBXwP7Afu34JyqugH4P3RZ5puBw4FnV9XyNv4ZwOVt/KuB/XrWNZlz6MoonkAX2K+kuylt7GOiJZOMfSPdTg230tVTf6qnbfP22I10ZQvXc8/2aycCj2r7EX+h/aX1HGBh6/98uq3baOt4WJIVU3gukqQZ5D4EvvOSLOo5Xtk39TZ0pX1jlrbHJvJy4KurW683t2nkVNWvgedN0HxSO3r7L6XbbmtM+tqL7qP68a51VN/5J4FPjte3tR8ySdsSYKLge6zPd4A/maDtWO6pY56yqlowSdupdLXCY+cL+9q/TXfzWq93tLZvMMHzafXN8/seWwT8vuSkbZ821vYrYM7Wj9rV/YMlaXZYXlV7TdI+XmJ43H8j2p72Lwf2Xd1FDXwlSZI0kLW4NdlSur3tx2wL/Ka/U/s096P0fII6GUsdpBHXyhVWjHNMVsIgSdJ9MmiZwxRj5QuBXZI8vO0F/wK6nZruuX7yMLqyuZdU1c+nMqkZX2nE9ZcrjIr+sgpJ0uhJ1k6FWlXdleQw4GxgfeBjVbUkyatb+/F0pXdbAse1TYPuWk35hIGvJEmSZp6qGrtJvPex43t+fgXdlz5N2X0KfBcvhgUL7ssMo2vxsmOZf/jrhr0MSZKkoRm1b24bOPA9+OA1uQxJkiSNlCl+GcVMMnDg+8pXdsdsteAks72SJEmjxBpfSZIkTVsYve3BDHwlSZI0kFlT6iBJkqTZbcTi3pHLUEuSJEkDMeMrSZKkgVjqIEmSpFlhxOJeA19JkiRNX4D1RizytcZXkiRJs4IZX0mSJA1kxBK+Br6SJEkaRJHUsBcxLQa+kiRJGogZX0ka0LLLlvOeA08c9jLWiDd+/uXDXsIas8n6o5XRkaSJGPhKkiRp2hL38ZUkSdIsMWJxr4GvJEmSBjNq++KO2nolSZKkgZjxlSRJ0kCs8ZUkSdIsEEatytfAV5IkSdPWhb2jFfha4ytJkqRZwYyvJEmSBpKMVg7VwFeSJEkDGq1SBwNfSZIkDSDW+EqSJEkzkRlfSZIkDWi0Mr4GvpIkSRrIqN3cNlqrlSRJkgZkxleSJEkDstRBkiRJ67iM4K4OBr6SJEkayKgFvtb4SpIkaVYw8JX0e0m2SvLtJLcmOWbY65EkzXTrDXgMh4GvRlKSfZOcn+TmJDckOS/J3kkWJlmVZEWSW5JclORZPeN2SFKtfUWSK5IcPsE1dklye5JPTDJ+RZIjetpfl+Tydu3fJPm3JBv0tF+RZGXP2HP6rvmgJKcluSnJjUlO7Wv/8yT/neS2JFcled4UXqsTklya5O4kC1fT/ZXAcmDzqnrD6uae5JonJfmHQcdLkkZD8v+3d+/BVpV1GMe/jygypYRIjSlKmSjSRRIUUktmahScMRtxRsW8MClTjWNlmmZOf9SMWZNmlqkUChqmZWmG5CW7zOQt0cxLXqJG8wYpiXoUVPDpj/Weabk7Z5/jDthnu5/PmTWw1vuud/32ZvbZv/3u33pRS1u7pMY3Oo6kkcBi4DPAT4HhwIeBl0uXW23vo2pxweOAyyWNtb2qNswo22slTQH+IOlO2zc2XOo84I5+whhle20fx38FLLC9StJo4ErgBODsWp8Dbf+mn3F/Ua45DngJeF/tcU8ELgOOBm4E3gaM6mecur8AVwDfHETfccBfbXsQfSMioquJrOrQRe5efjfTF0xvdxjdaGcA2z8p+6uBGwAk7d7byfZrki4FLgDG00cSa3uppPuBSVTJJGWcw4BVwC3AToMNzPbfa7sCXhvs+ZL2A7YHptteVw7/udbldOBC278u+yvLNlBM55Xx1wxw/QXAEYAlfR74BPA88F1gV6rn+efAibZfUfWR/exyzubAo8BsYK+GcX5n+0BJHwTmU/1bLAEMLLN9+kCPISIiYn1IqUOLZr9/NpO2mdTuMLrVw8A6SQslzZS0VV+dJA0D5gCvUiVlffWZRjWruqx2bCTwNaDZV/2PSnpc0sWSxjSMOVvS81QlA7sBFzacu0jS05JukLRb7fg04CFgoaSVku6QtG9DO5LulfSUpB+XWeX1wvYxwCLgW7a3KLPS64AvAGOADwEfBT5bTtkP+AjVB5FRwKHAStvzGsY5UNJw4GrgUmA08DNg1vqKPSIi2kMt/rRLZnxbNHfyXOZOntvuMN60NKf/F4Xt5yXtA5wC/BDYRtISqrIGgGmSVgFvBdYCn7T9r4ZhnpG0OTACOIsqKev1dWC+7cf6qEN6BtgDuBvYmqocYhGwfy2+y4DLJI0HjgJW1M4/AriLajb4c8D1kiaUMoyxVMnksVQJ+yzgl5J2sv1MaT+y9HkSWAh8r4y5Qdi+s7b7iKQLgX2Bc6g+UGwJTAD+ZPuBJkNNAzYDzillFFdKOnEDhR0RERtNZ82hdla0EYXtB2wfY3ss1YzttlTJGMBttkcBWwHXUNX/NhoDbAGcBEynSsqQNAn4GPCdfq7bY3up7bW2VwDHA/uVWeLGvn8D7gd+UDt2s+3Vtl+y/Q2qcore+FYDj9ieb/tV25cDjwF719ovtv2w7R7gDOCAAZ+s/4OknSUtlrS8zGKfQfXcYfu3wPepkv8V5Sa6/3keim2BJxpqh/uchY+IiM7RaTO+SXyj49l+EFhA7UawcryH6mv5I0t9aeN562yfBazhv1/fTwfeBfxT0nKqxHiWpLv6u3z5s79X8abAe5qFXzv3ntp4fRmofUM4H3gQGG97JHAatcdq+1zbk4H3UpU8nNzb1DDOU8B2ev0U+g4bLOqIiIg+JPGNjiNpgqQvShpb9rcHDgdua+xreyXwI+CrTYY8E/iSpBHAPKpEdVLZLgCupZQySJoqaRdJm0jaGjgX+L3t50r7sZLeUf4+EfgycFPZ30HS3pKGSxoh6WSq2dObSxxXAVtJOlrSMEmHANvV2i8G5kjaUdJbqEo9Fg/i+RpeHpuAzcq1B/va35LqBrceSROoVtLoHXeP8nxsBrxI9QGi96a8FcCOtXFupSo7OUHSppIOBvYcZAwRETEEtbqUWTuXM0viG53oBWAqcLukF6kS3vvo/2a0c4ADJH2gn/ZrgWeB40oJwvLeDegB1th+uvTdEbiuxHAf1RJqh9fG2hu4t8S1pGynlbYtqWZQnwWeAGYAM0tyju1/Ax+nmmV+DjgVOKjU92L7IuAS4HaqMoGXqZZKG8gNVGUSe1El9qupbkpD0hFlVYv+nES1UsMLVPXUV9TaRpZjz5Z4VgLfLm3zgYmq1iO+2vYrwMHAMaX/oVRLt1Hi2EFSzyAeS0REDClqcWsPNVuuc8qUKV66dOlGDCeiUtbVndLuOGLDKcunPV5fzmyT4aO9+Tb7939SBznpqk+1O4T1ZsSwN8+yzl+ZNK7dIUR0DGmXpu/Fk6fs6ltuX9jS2CM2ndqW9/nM+EZEREREV0jiG9HhSrlCTx9bsxKGiIiI9aCzSh2yjm9Eh7O9iGot4Y5S/sOMiIjoWO29Ua0VSXwjIiIiokWdlfim1CEiIiIiukJmfCMiIiKiJeqwOdQkvhERERHRos4qdUjiGxEREREtEOqwxLez5qcjIiIi1LP+1QAAAcJJREFUIlqUGd+IiIiIeMMEWc4sIiIiIrpFZxUPdFa0EREREREtyoxvRERERLSk025uS+IbERERES0QWc4sIiIiIrpCp93clhrfiIiIiBhyJM2Q9JCkZZJO7aNdks4t7fdI2n2gMZP4RkRERESLNmlxa07SMOA8YCYwEThc0sSGbjOB8WWbC5w/mGgjIiIiIt4wtfgzCHsCy2z/w/YrwOXAQQ19DgIuceU2YJSkdzaN13b/jdLTwKODiS5iPRtn++3tDiI2rvzOiYgYUpq+F0u6DhjT4tgjgDW1/Xm259XGPgSYYfvYsn8kMNX28bU+i4Ezbf+x7N8EnGJ7aX8XbXpzWxKPiNiY8jsnIqJz2J6xAYfva1q4cbZ2MH1eJ6UOERERETHUPA5sX9sfCzzZQp/XSeIbEREREUPNHcB4Se+WNBw4DLimoc81wFFldYdpwHO2n2o2aNbxjYiIiIghxfZaSccD1wPDgIts3y/p06X9AmAJcACwDHgJmDPQuE1vbouIiIiIeLNIqUNEREREdIUkvhERERHRFZL4RkRERERXSOIbEREREV0hiW9EREREdIUkvhERERHRFZL4RkRERERX+A9JrtF27nuhGgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 792x576 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "_ = plot_composite_matrix(matrix, labels)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python (myenv)",
   "language": "python",
   "name": "myenv"
  },
  "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}