SMG-epi-HSP47-mean-intensity-plotting.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "mkdir: ../jupyter_figures: File exists\n"
     ]
    }
   ],
   "source": [
    "import sys, os, glob, warnings, logging\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "from scipy import stats\n",
    "from sw_plotting import change_bar_width, plotCountBar\n",
    "from sw_utilities import tukeyTest\n",
    "\n",
    "logging.basicConfig(stream=sys.stdout, format='%(asctime)s - %(levelname)s - %(message)s', level=logging.ERROR)\n",
    "# logging.basicConfig(stream=sys.stdout, format='%(asctime)s - %(levelname)s - %(message)s', level=logging.INFO)\n",
    "# logging.basicConfig(stream=sys.stdout, format='%(asctime)s - %(levelname)s - %(message)s', level=logging.DEBUG)\n",
    "\n",
    "# ignore warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "# plotting configuration\n",
    "font = {'family' : 'Arial',\n",
    "        'size'   : 7}\n",
    "matplotlib.rc('font', **font)\n",
    "plt.rcParams['svg.fonttype'] = 'none'\n",
    "\n",
    "# Make a folder if it is not already there to store exported figures\n",
    "!mkdir ../jupyter_figures\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>index</th>\n",
       "      <th>file_name</th>\n",
       "      <th>mean_intensity</th>\n",
       "      <th>groups</th>\n",
       "      <th>gland_id</th>\n",
       "      <th>sample_id</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>epi1-1.nd2</td>\n",
       "      <td>1.256417</td>\n",
       "      <td>sg-Control</td>\n",
       "      <td>1</td>\n",
       "      <td>sg-Control_1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>epi1-2.nd2</td>\n",
       "      <td>0.953481</td>\n",
       "      <td>sg-Control</td>\n",
       "      <td>1</td>\n",
       "      <td>sg-Control_1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>epi1-3.nd2</td>\n",
       "      <td>0.941935</td>\n",
       "      <td>sg-Control</td>\n",
       "      <td>1</td>\n",
       "      <td>sg-Control_1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>epi2-1.nd2</td>\n",
       "      <td>1.059278</td>\n",
       "      <td>sg-Control</td>\n",
       "      <td>2</td>\n",
       "      <td>sg-Control_2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>epi2-2.nd2</td>\n",
       "      <td>0.841262</td>\n",
       "      <td>sg-Control</td>\n",
       "      <td>2</td>\n",
       "      <td>sg-Control_2</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   index   file_name  mean_intensity      groups gland_id     sample_id\n",
       "0      0  epi1-1.nd2        1.256417  sg-Control        1  sg-Control_1\n",
       "1      1  epi1-2.nd2        0.953481  sg-Control        1  sg-Control_1\n",
       "2      2  epi1-3.nd2        0.941935  sg-Control        1  sg-Control_1\n",
       "3      3  epi2-1.nd2        1.059278  sg-Control        2  sg-Control_2\n",
       "4      4  epi2-2.nd2        0.841262  sg-Control        2  sg-Control_2"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "input_folder = '../data/190523-CRISPR-KO-SMG-IF-HSP47/'\n",
    "# Read in 1st dataset\n",
    "\n",
    "df0 = pd.read_csv(input_folder+'190523-LV-188-control-anti-Hsp47-647-anti-Col4-red-mean-intensity-5478.txt')\n",
    "\n",
    "# average measurements on different z slices for the same image stack\n",
    "df = df0.groupby('file_name').describe()['mean_intensity']['mean'].reset_index()\n",
    "df.rename(columns={'mean': 'mean_intensity'}, inplace=True)\n",
    "\n",
    "# annotate groups and samples (2 epi buds, 3 and 4 FOVs each, 21 z slices per image stack)\n",
    "df['groups'] = ['sg-Control']*len(df)\n",
    "df['gland_id'] = ['1']*3 + ['2']*4\n",
    "df['sample_id'] = [df.groups[i] + '_' + df.gland_id[i] for i in range(len(df))]\n",
    "\n",
    "# pass to new variable\n",
    "df1 = df\n",
    "\n",
    "# Read in 2nd dataset\n",
    "df0 = pd.read_csv(input_folder+'190523-LV-228-anti-Hsp47-647-anti-Col4-red-mean-intensity-4974.txt')\n",
    "\n",
    "# average measurements on different z slices for the same image stack\n",
    "df = df0.groupby('file_name').describe()['mean_intensity']['mean'].reset_index()\n",
    "df.rename(columns={'mean': 'mean_intensity'}, inplace=True)\n",
    "\n",
    "# annotate groups and samples (2 epi buds, 3 and 4 FOVs each, 21 z slices per image stack)\n",
    "df['groups'] = ['sg1-Serpinh1']*len(df)\n",
    "df['gland_id'] = ['1']*3 + ['2']*4\n",
    "df['sample_id'] = [df.groups[i] + '_' + df.gland_id[i] for i in range(len(df))]\n",
    "\n",
    "# pass to new variable\n",
    "df2 = df\n",
    "\n",
    "# Read in 3rd dataset\n",
    "df0 = pd.read_csv(input_folder+'190523-LV-229-anti-Hsp47-647-anti-Col4-red-mean-intensity-7337.txt')\n",
    "\n",
    "# average measurements on different z slices for the same image stack\n",
    "df = df0.groupby('file_name').describe()['mean_intensity']['mean'].reset_index()\n",
    "df.rename(columns={'mean': 'mean_intensity'}, inplace=True)\n",
    "\n",
    "# annotate groups and samples (2 epi buds, 2 FOVs each, 21 z slices per image stack)\n",
    "df['groups'] = ['sg2-Serpinh1']*len(df)\n",
    "df['gland_id'] = ['1']*2 + ['2']*2\n",
    "df['sample_id'] = [df.groups[i] + '_' + df.gland_id[i] for i in range(len(df))]\n",
    "\n",
    "# pass to new variable\n",
    "df3 = df\n",
    "\n",
    "# Read in 4th dataset\n",
    "df0 = pd.read_csv(input_folder+'190523-LV-230-anti-Hsp47-647-anti-Col4-red-mean-intensity-4549.txt')\n",
    "\n",
    "# average measurements on different z slices for the same image stack\n",
    "df = df0.groupby('file_name').describe()['mean_intensity']['mean'].reset_index()\n",
    "df.rename(columns={'mean': 'mean_intensity'}, inplace=True)\n",
    "\n",
    "# annotate groups and samples (2 epi buds, 3 FOVs each)\n",
    "df['groups'] = ['sg3-Serpinh1']*len(df)\n",
    "df['gland_id'] = ['1']*3 + ['2']*3\n",
    "df['sample_id'] = [df.groups[i] + '_' + df.gland_id[i] for i in range(len(df))]\n",
    "\n",
    "# pass to new variable\n",
    "df4 = df\n",
    "\n",
    "# Concatenate the 2 dataframes and reset the index\n",
    "df = pd.concat([df1, df2, df3, df4])\n",
    "df.reset_index(inplace=True)\n",
    "\n",
    "# Normalize intensities by the average of sg-Control\n",
    "df['mean_intensity'] = df['mean_intensity'] / df1.mean_intensity.mean()\n",
    "\n",
    "# Check the header\n",
    "df.head()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<AxesSubplot:xlabel='groups', ylabel='mean_intensity'>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXkAAAD/CAYAAAAUnaZMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAAAY9ElEQVR4nO3de3hV5Zn38e+dEHLgaCBIQBQRq2LFEQNlUAGxVuzVqnXmhb5S0daKp1doa6X24LRFK7Uz08NceKKiI1WrteLUlrdqFYotHjBIVUSpQgGhHEI5ExIScs8faxNCCMkirJ2198rvc137Staz19rrzkPy2w/PXgdzd0REJJly4i5ARETSRyEvIpJgCnkRkQRTyIuIJJhCXkQkwRTyIiIJ1iHuAhrq2bOn9+/fP+4yRESyyuLFize7e0lTz2VUyPfv35/y8vK4yxARySpmtvpwz2m6RkQkwRTyIiIJppAXEUkwhbyISIIp5EVEEiyjjq7JdO+t38FT5WvpnJ/LhOEncGzXgrhLEhFplkI+pPfW7+DSexayt7YOgKcWr+XFr42iU766UEQyl6ZrQnp68dr6gAdYv72K+cs3xViRiEjLFPIhdS44dMTeWaN4EclwCvmQrvjE8fTtXli/PHxAMeed3ORZxCIiGUND0ZB6dSngxa+NYv7yTXTO78C5A3uSk2NxlyUi0iyF/BEo7JjLp88ojbsMEZHQNF0jIpJgCnkRkQRTyIuIJFjkIW9mY8zswQbLE8xslpk9bWZ9o96fiIgcXqQhb2YDgSFAw/P9d7r7NcCjwIgo9yciIs2LNOTd/UN3/49Gbc+aWQnwOeB3jbcxs0lmVm5m5RUVFVGWIyLS7qV9Tt7MTgR+AEx29z2Nn3f3me5e5u5lJSU6uUhEJEppC3kzu9LMBgGPAMXADDMbna79iYjIodJyMpS7f6HB4sh07ENERFqmQyhFRBJMIS8ikmAKeRGRBFPIi4gkmEJeRCTBFPIiIgmmkBcRSTCFvIhIginkRUQSTCEvIpJgCnkRkQRTyIuIJJhCXkQkwRTyIiIJppAXEUkwhbyISIIp5EVEEiwtd4ZKqpUVu5jz5jo65Xdg/NB+FHfqGHdJIiLNUsiH9OGmnVwyYyGVe/cB8MtFa3j+KyMp7Jgbc2UiIoen6ZqQflW+tj7gAdZsqWT+8k0xViQi0jKFfEgFHQ7tqoI8dZ+IZDalVEhXfOIEenXJr18+6/jujDy5JMaKRERapjn5kHp3K2DOjSO4Z/6HdC3MY/KYk+mQq/dIEclsCvmQNmyv4l/ue4WNO6oBKF+1lScnDVfQi0hGU0KF9Pjrq+sDHmDx6q28/EFFjBWJiLRMIR/Snpp9h7btrYuhEhGR8BTyIY0r60dh3oFj4vsVFzLm1F4xViQi0jLNyYd08rFd+O3N5zLnzbV0yu/A54f204lQIpLxIg95MxsDXOHuX04tXw6MJfhfw/XuXhv1PtvKwF6dmTr21LjLEBEJLdLpGjMbCAwBCho0j3P3ScDLwAVR7k9ERJoXaci7+4fu/h+NmmtSX9cDpY23MbNJZlZuZuUVFTpaRUQkSm3xwev+w1JKgQ2Nn3T3me5e5u5lJSU6g1REJEpp++DVzK4EFgNPmdnPAQOuS9f+RETkUGkJeXf/QoPFZcDcdOxHRESap+PkRUQSTMfJH6G1Wyv50XPL+evGnYw5tRdf+eTH6NjEZYhFRDKBQv4IXTt7Me+t3wHA+xt2Uudw28U6dl5EMpOGoM34x65qJv9yCf88/SVueuxN3l67rT7g93vxvY0xVSci0jKN5Jtx25x3+MOyIMTnvrOebZV76VaYx/Y9NfXrnFTSKa7yRERapJF8M/7U6FLCr6z8B3d97uN0LQjeG08q6cQ3Lz4tjtJERELRSL4Zp5V2ZcmabfXLp/buygWnHcvC20rYsnsvxxcXYWbxFSgi0gKN5Jsx/fIzGNirMwD9exRxXPdCPv7d5xkxfR5/WLZRAS8iGc/cPe4a6pWVlXl5eXncZRxiy+69PL90Pd98ZulB7b+fch6nlXaNqSoRkYCZLXb3sqae00g+hOJOHXm30VE1AEvXbY+hGhGR8BTyIZ1zUs+DljvkGMMH9IipGhGRcBTyIV18Rim3XXwqxx1TyKm9u3DPhCH0Ky6KuywRkWbp6JojcP2ok7h+1ElxlyEiEppG8iIiCaaQFxFJMIW8iEiCKeRFRBJMIS8ikmAKeRGRBFPIH6G6usy5DISISEt0nHxIdXXOHXOX8fjra+iU34FbPvUxJnzihLjLEhFplkbyIc1Zso6HF66iuraOLbv38u1nlvLBxp1xlyUi0iyFfEhL1mw9pO0vH21r+0JERI6AQj6kYScWH7ScYzC0f/Fh1hYRyQyakw/pkjP7sKJiN4+9trp+Tr5/T93fVUQym24aIiKS5XTTEBGRdirSkDezEWb2SOrRPdX2GTP7qZk9YWZ9o9yfiIg0L+qR/CTgOmAWMD7Vdi4wECgGNke8PxERaUbUIZ/r7lXAeqA01fYscAkwBxjeeAMzm2Rm5WZWXlFREXE5IiLtW4shb2bfNrN+IV+v0szyCQJ+Q6rta+5eB1QAh9wU1d1nunuZu5eVlJSErVtEREIIcwjlbGCimZUCC4Ffu3vNYdadCTwA5AHPmdkg4CUzexioA26KoGYREQkpTMj35sAIPAd4DBjX1Iruvhi4ulHzMuC+VtYnIiJHIUzInwPc6e5bzKw38Nc01yQiIhE5bMinpmduBEYBxWZmwLnufn5bFSciIkfnsCHv7uvN7OfAe8CfAQd+1FaFiYjI0WtpumYicAowFjCCoJ+Y7qJERCQazYa8u99pZnnAPmAw8FGbVCUiIpEIczLUwwQfvk4FfpLeckREJEphQn4ncBFwF8GIXkREskSYkJ8PbCWYk38+veWIiEiUwhwn/xHBh6+XE3zwKiIiWSJMyH8RuBuoTXMtIiISsTDTNVuBQiA39RARkSwRZiRfSTBVs9+0NNUiIiIRCzOSX0Qwkn8bWJDeckREJEphQv4ygg9cfw9cn9ZqREQkUmFCfh+QD/Qn3PSOiIhkiDAh/58EQX8t8K30liMiIlE6bMib2bFmdg3BpQzeJ7ga5T1tVZiIiBy95qZfKoAPgNeADwnm5W9oi6JERCQazV1Pvg542cx6AlcRjPod+FIb1SYiIkcpzAepI4FJ7q6Lk4mIZJkwIV8EXGxmlQDuPi+9JYmISFTChPxCoDj1EBGRLNJsyJvZdwiuQLmfA7PTWpGIiESmxdv/NdVuZhe5u64tLyKS4cKcDNWUsyOtQkRE0qK1IS8iIlmgtSFvkVYhIiJp0eLRNWZmQF9Sbwjuvgb4Y3rLEhGRKIQ5hPIRYDXB7f8cmObuC5ta0cxGANelFqe4+zYzGwWMI7ir1PfcfcPRly0iImGECfn33f2ukK83iSDkhwHjgQcILomwkeByxZtaU6SIiLROmDn5MWb2azObbWYtHSOf6+5VwHqgNNXWH/gusAwY23gDM5tkZuVmVl5RUXEEpYuISEvCjOTHAn0IPmzt1sK6lWaWTxDw+6dl1rn7XjPbQhNvKu4+E5gJUFZW5mELFxGRloUJ+V8BBcBKYBAwppl1ZxJM0eQBz5nZIOAJM5sF7AW+enTliojIkQgT8uuA3wLVwCXNrejui4GrGzUvA+a2pjgRETk6YebkVwNbgClAr/SWIyIiUQozkn8KGA7cQTCqFxGRLBFmJP994FygErgvveWIiEiUwoT8DmCnuy8nCHoREckSYUJ+OTDFzHYAb6W5HhERiVCYOfl/Bj4L7CO4rIGIiGSJMCP59wkuSZBLuDcFEWljSzcvZeLvJ3L+r87nztfupHpfddwlSYYIE9pDgL+jkbxIRqrZV8PN825m857NADy5/Em6duzK5CGTY65MMkGYkP8QKAd2pbkWEWmFFdtX1Af8fos2LIqpGsk0YUJ+N3AZwbVrHJiWzoJE5Mgc3+V4Oud1ZlfNgXHYoB6DYqxIMkmLIe/u32+LQkSkdYryiph+3nTueO0ONlVu4py+53DTP90Ud1mSIfRBqkgCjO43mpHHjaSqtoqivKK4y5EMoht5iyREjuUo4OUQCnkRkQRTyIuIJJhCXkQkwRTyIiIJppAXEUkwhbyISIIp5EVEEkwhLyKSYAp5EZEEU8iLiCSYQl5EpJEts2fzwajRfDBqNFsefSzuco6KLlAmItLA7kWL2HjX9PrljXfeScGg0ygaMiTGqlpPI3kRkQYqy8sPbXvj0LZsoZAXEWmgcPCZh7adOTiGSqIRacib2QgzeyT16N6g/TQzeyHKfYmIpEPnc8+h583/j5xOncjp3JmSKZPpNHx43GW1WtRz8pOA64BhwHjgATPLB64FKiLel4hIWpTcdBM9b7gBzDCzuMs5KlFP1+S6exWwHihNtX0D+DHB/WEPYWaTzKzczMorKvQ+ICKZwXJysj7gIfqQr0yN3EuBDWZWBAwHvgMMM7MxjTdw95nuXubuZSUlJRGXIyLSvkU9XTMTeADIA54D+rv7pwHM7FF3nxfx/kREpBmRhry7LwauPsxzX4hyXyIi0jKdDCWSxV5e+zKz3plFrddy5aArGdt/bNwlSYZRyItkqZXbVzJl3hRqvRaAqQumUtqplDNLDj3OW9ovnQwlkqUWrltYH/AAjvPy2pdjrEgykUJeJEsN6DYgVJu0bwp5kSw1os8Ixp8ynlzLxTBGHTeK/Jx8dtfsjrs0ySDm3uQ5SrEoKyvz8iYuDiQih7e1aiv3v3U/j7//OADd8rsx61OzOKX4lJgry167Fy1iy+zZmBnFEydSNHRo3CU1y8wWu3tZU89pJC+S5fb5Pp5c/mT98vbq7Tzw9gMxVpTdqlesYM01X2bXiy+x8w8vsuZL11D9t7/FXVarKeQl/Sq3wPPfhsf+D7zxIOz/3+PGZfDafbDqz/HWl+V2VO9gn+87qG1b9bZ4ikmAnS++BDU19cteU8Oul16KsaKjo0MoJf2euALWvBp8/8ELUL0Tuh8PT38ZvC5oP+/rcMHt8dWYxQZ0H8DgksG8XfF2fdtlAy+Lr6Asl9enz6FtffvGUEk0FPKSXtvXHgj4/d75dfB1f8ADvDoDRn4d8grbrrYEufeCe3nsvcdYs3MNFx5/IReccEHcJWWtrmMvYsfcuez64x8B6DxmDF0uvDDeoo6CQl7Sq6A75BVBTeWBtq59YOuqg9er2xc8pFW65Xfjxn+6Me4yEsHy8uh3/31Ur1gBZuQPyO7DUjUnL+mV3xkunAY5qfFEpxIYczsMbxRIZ18VrCuSIfJPOinrAx40kpe2MOxaOO2z8I8V0PdsyCuA0sHQYyCsmAe9z4BBl8VdpUgiKeSlbXTpHTz221UBr94DK+cHIV9yChx7enz1iSSUpmskHs99A/76e6itgrVvwFNfjLsikURSyEs8Vjc64mbz8uB4ehGJlEJe4tGv0WniPQZC4THx1CKSYAp5icfF/w4njQHLCebk//UhSMBNk0UyjT54lXh0ORaufCbuKkQSTyN5EZEEU8iLiCSYQl5EJMEU8iIJ8Pr61/nsM5/lrNlnceuCW3V3KKmnkBfJclW1Vdyy4BZW7VhFrdfy3KrnuPcv98ZdlmQIhbxIllu1YxXbq7cf1PZWxVsxVZMcXlfHrgUL2PLII/V3hqpZt47VE6/ivUGns2r856leuTLmKlumkBfJcid2O5Hu+d0Pajur11nxFJMg62+/nY+uu56N03/IyksuZfcrr7D+u9+jctEiqKtjz1tv8fep34i7zBYp5EWyXH5uPj8e/WMGdh9IQW4BnxnwGW4484a4y8oqdZWVVNx7L2snT2HrE0+yd+NGts9pcB5HTQ3/eHAWe5YsOWi7qqVL8draNq72yOhkKJEEGNp7KM9cqpPLWmvd1KnsejG4j+vOF17gmL9ddeBexCleW0vhkCHs/tOf6tsKBg/GOmR2jEY6kjezEWb2SOrRPdU2wcxmmdnTZpa9N0oUkUTat3Mnu16ad1Dbrnnz6HLRRQcacnIovmoipXdMo9OIEVheHoVnn03fH93dxtUeuajfgiYB1wHDgPHAA8BOd7/GzD4HjACeinifIiKtlpOfT06XLtTt2FHf1qFHD/rc/UN2nD+amjVr6Hz+GArP+DgAxz80K6ZKWyfqOflcd68C1gOlAO7+rJmVAJ8Dftd4AzObZGblZlZeUVERcTkiIs2zjh3p9fVbIDcXgJyiIjqecAJ/HXEOG6fdATm59QGfjaIeyVeaWT5BwG8AMLMTgW8Ck919T+MN3H0mMBOgrKzMGz8vIpJux4wbR+fzzqP6gw+o21PFuilTAHBg84wZFA05i04jRsRbZCtFPZKfSTBFcx2w28wGAY8AxcAMMxsd8f4km9VWw7o3oWp7y+uKpFleaSmdR46kesWHhzy3552lMVQUjUhH8u6+GLi6UfPIKPchCfH3JfDYONi9CfKK4NJ74OOXx12VCJ2GDWNzo7aioUObXDcb6Dh5icfz3wkCHqCmEv7/rbCvJt6aRICisjKOvf07dOhTSl6fPvT+7r9RNCR7Ty7L7AM8Jbm2rTl4uXIz7N0Nhd1jKUekoeIJEyieMCHuMiKhkbzEY9AlBy8PGK2AF0kDjeQlHp/8XhDqK+YH93gdlfnXABHJRgp5iUduHoy8NXiISNpoukZEJMEU8iIiCaaQFxFJMIW8iEiCKeRFRBJMIS9tp2o7vPs/8NEbB7fX1cH8u+Cng2HWp2D1K7GUJ5JEOoRS2kbFcnhoLOzZEiyfdSVcOiP4vnwWLEjdfGHbanh8PHz1XSjoGk+tIgmikby0jYU/OxDwAEt+AZtTV/tbMf/gdat3wNpGo30RaRWFvLSNpi4nvL+t9xkHt+d0gF6D0l+TSDugkJe2MWQiYAeWew2CV2fA3SfCqj/DiaOC5wu6w2d+Al1LYypUJFk0Jy9t42MXwVXPwtKnoWtf2LgM3p0TPLf6z9DrdLhtNXQohA4d461VJEEU8tJ2ThwZPAB+3Gg6ZtO7ULdPAS8SMU3XSDz6NLoJwzH9g6kaEYmUQl7icfHd0G948H2PgXD5g5CjX0eRqGm6RuLR7Ti45nmo2QN5hXFXI5JYGjpJvBTwImmlkBcRSTCFvIhIginkRUQSTCEvIpJgCnkRkQQzd4+7hnpmVgGsjruOEHoCm+MuIkHUn9FRX0YrW/rzBHcvaeqJjAr5bGFm5e5eFncdSaH+jI76MlpJ6E9N14iIJJhCXkQkwRTyrTMz7gISRv0ZHfVltLK+PzUnLyKSYBrJi4gkmEI+QmbWzcw6xV2HiMh+CvnDMLMCM7vHzB4ys7lmdkqIzb4CHPbmpGZ2tZl9MrIiM4yZ5ZrZvY1/RjMbYWa/NLMHzewbR/iatzfRNtrMvtyorZ+ZLWhd5ZnJzKab2QNm9rCZ5TRoV38eITP7VurveVqj9sT3Zbu8nryZjQfOB7oCywlOePgEMM7dV6VWuwH4jbu/YGYlwHAzKwZuJHhznAZ8Dfg7cDLwPeBcYLeZnUxw1+r/Bm4CqoC5bfGzpVtzfZdqtyY2uxS4293/ktoeM/tB6jV2A/cDjwPPAmXAX4AewDeBU8xsNHBzan/rgXeA/2tm5wJvuvt/AZOBbVH/vOnWTH+OBxa6++/M7GcE/VGR2kz92YQWfjffcve7zOzxRpslvi/b60j+eGAL8ADQ291vBt5otM6pwGsA7l7h7r8FrgO+BHyVIOzzgZ8AjwAjgIXA/6TapwJnATPc/RqCX7QkOGzfufvDwOtNbDMN+LSZ3QP0NrPTCf74dgMnAd2ARe7+Q6Da3e8E/ghc2OA1Frj7t1LbAcx196sJ/vBw91uBnRH+nG3lcP3pqYAfBux194oG26g/m9bc7+bc1MjbzazhQCTxfdkuR/LAC0AdwT/w/jtHO9DHzH5I8G68DBgO7B/Jf5GDR6n73yArgZrUcsNDlXY1amtqhJuNDtd3BzGzMwlGO/MI/lB+5O61ZjYDWAIsdvfbzGwcwR/Urv2bpr4WADsavGRl6mtt6uv+5+qO+ieK12H708zGEgwepqo/Q2muL89x9zvM7BZgsJm1m75sryP504FbCP4B15jZfwLnASvd/fPu/gOC42OvMLPZBCOD36TaZgLTgZ818borCUb7+/0SuDH1y/NEun6YNna4vqtquJK7v5Xqy5kEo6k5ZvYLIMfdXwaOMbP7CUY7tQ027WJmPwFGE/wRJl1z/flzYADwC2CD+rNFzfXlJWZ2LzAIeL899WW7P07ezL5PML2Sm/pvlYSUjr4zs0fd/QtRvFa2UX9GR315QLsPeRGRJGuv0zUiIu2CQl5EJMEU8iIiCaaQFxFJsPZ6nLwIAKmTjaYAewluPdkfeIjgvIgqYCnwLjDQ3R80s/8mOJP5q8CrBIfNrgeuITh+eoq7b2zTH0KkGRrJS3t3LXAVcG9q+adAF+AFd7+B4DjrpvwudRbkpwiuV7SXIPy3pbFWkSOmkby0d43P2m3qTOVaoEPqdPhuDdoB8oC3CUb8NxAcm/2bNNcsEppCXtq7B1OPPOBfgEcJTo+/LzWVM5/gMhffIrie0e7UduPNbDAwh2AkPwnYB7zfptWLtEAnQ0m7ZmaXElyjqAfBNUw+DLHNaFJz9OmtTuToKeRFRBJMH7yKiCSYQl5EJMEU8iIiCaaQFxFJMIW8iEiCKeRFRBLsfwEvHRTGzo/bUwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# quick visulization\n",
    "sns.swarmplot(df.groups, df.mean_intensity)\n",
    "# sns.barplot(df.groups, df.mean_intensity)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Summary of test:\n",
      "       Multiple Comparison of Means - Tukey HSD, FWER=0.05       \n",
      "================================================================\n",
      "   group1       group2    meandiff p-adj   lower   upper  reject\n",
      "----------------------------------------------------------------\n",
      "  sg-Control sg1-Serpinh1  -0.8044  0.001 -0.9773 -0.6315   True\n",
      "  sg-Control sg2-Serpinh1   -0.636  0.001 -0.8388 -0.4332   True\n",
      "  sg-Control sg3-Serpinh1  -0.7357  0.001 -0.9157 -0.5557   True\n",
      "sg1-Serpinh1 sg2-Serpinh1   0.1684 0.1256 -0.0344  0.3712  False\n",
      "sg1-Serpinh1 sg3-Serpinh1   0.0687 0.6934 -0.1112  0.2487  False\n",
      "sg2-Serpinh1 sg3-Serpinh1  -0.0997 0.5477 -0.3085  0.1091  False\n",
      "----------------------------------------------------------------\n",
      "p values of all pair-wise tests:\n",
      " [0.001      0.001      0.001      0.12562561 0.69336076 0.54766578]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<statsmodels.sandbox.stats.multicomp.TukeyHSDResults at 0x7f84c9a104c0>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# df.head()\n",
    "tukeyTest(df.mean_intensity, df.groups)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAI3CAYAAAC7y7EtAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAC4jAAAuIwF4pT92AABiOUlEQVR4nO2deZhcVZnGf193kk46K9kTCIQk7BD2fXUDhSCLiIojyqAjoqiMy6gjMm6j4r7iuIuioCCIqOwKgqyyhCUsgSTsJCF7Oumku77545xL3a70Vt1VdW91v7/nOc+tuvfcU2919X3r1Lnf+Y65O0IIIeqHhqwFCCGEKA8ZtxBC1BkybiGEqDNk3EIIUWfIuIUQos6QcQshRJ0h4xZCiDpDxi2EEHWGjFsIIeoMGbcQQtQZMm4hhKgzZNxCCFFnyLiFEKLOkHELIUSdIeMWQog6Q8YthBB1hoxbCCHqDBm3EELUGTJuIYSoM2TcQghRZ8i4hRCizpBxCyFEnSHjFkKIOkPGLYQQdYaMWwgh6gwZtxBC1BlDshYg6hMzGwscmdr1DLApIzlC1JphwIzU85vdfXWtXlzGLfrKkcAfsxYhRE44AbiqVi+moRIhhKgzZNxCCFFnaKhE9JVn0k+uvPJK5syZk5UWIWrKwoULOfHEE9O7numialWQcYu+0uFG5Jw5c9htt92y0iJE1tT0xryGSoQQos6QcQshRJ0h4xZCiDpDxi2EEHWGjFsIIeoMGbcQQtQZMm4hhKgzZNxCCFFnyLiFEKLOkHELIUSdIeMWQog6Q8YthBB1hoxbCCHqDBm3EELUGTJuIYSoM2TcQghRZ8i4hRCizpBxCyFEnSHjFkKIOkPGLYQQdYaMWwgh6oxcGbeZnWlmY7LWIYQQeSZXxg38GHjRzC41szea2ZCsBQkhRN7Im3EDDAdOAa4AXjCz75nZQRlrEkKI3JA3494V+CKwGDBgAvA+4DYze8LMPmNmszPUJ4QQmZMr43b3R939PHefDRwKfB9YTjDx2cD5wONmdruZvc/MJmQoVwghMiFXxp3G3W9393OA6cBxwG+AFoKJHwh8D3jezP5oZm82s6bs1AohRO3IrXEnuHu7u//V3f8NmAycClwELAOGAvOASwg3NX9kZodmp1YIIapP7o07jbtvcPfLgHOBjwML4yEDxgJnAreY2XwzOzUjmUIIUVXqxrjNbJKZnW1mNwEvAT8njHsbsJowlPJIfL478Fszu0whhUKIgUaujdvMxpjZu8zsWuA54LvAkcAQoB34K/A2YKq7/5u77wHsB9xIMPCTgM9mIl4IIapE7nqjZjYCeCPwVuD1wLDkUNzeB/wK+I27Ly09393vNbPjYr1dgHcA/11t3UIIUStyZdxmdjHBtJuTXXH7HGEo5CJ3f7indtx9k5ndRjDuidXQKoQQWZEr4yYMeySsJ8yevAi40d29zLYmx+0DlRAmhBB5IW/G7cBNBLO+3N1b+tHWfwJnuPvKiigTQoickDfj3tbdn6tEQ+7+VCXaEUKIvJG3qJLPm9lPzWyb3p5gZruY2S1mdkc1hQkhRF7IW4/7XYThkm8Dz/bynGbgMGBdlTQJIUSuyFuPO6FXNyLNrBk4rZxzhBCi3smkx21muxLirEtfPzHf+82MMnAUPSKEGCRk0uN290eArxDitCtRNqJJNkKIQUKWY9xfAFpLNJxP6D3/CHixh/ML8fyXCHHevR0TF0KIuiYz43b3TYTVbl7BzM6PDy909/m1V9U/zGwnwpDNcnfvdWRMJ+1sAzzTy+rvcfef9PW1hBD1R96iSs6I26czVdEH4kIOPwMqsaDDnhVoQwgxQMmVcbv7L7PW0Bdi6thLgUMq1GRi3Bt60WbdfckJIfpHVlEls5LH6RmO6f19IYvZkmY2kbACz2sq2Gxi3A+6+/0VbFcIMQDIqsf9RNx6iYaF9D0eu7StqhPTx/4A2DbuKlCZSJ3EuBXiKITYgqwm4KRD+bo7Vm6pGWZ2KXA1wbSdECXzjwq02wzsEJ/e39/2hBADj6x63F2tSlNPq9UcGLeLgXe7+41m9vcKtLs7xS9U9biFEFuQiXG7e6cG3dX+nPI88H3gO+7eWsF294pbB1aa2QXAcYT1NTcCDwMXAz92980VfF0hRJ2Qq6iSOuMwdy9Uod1kfLsduAcYkTrWRIgyOQR4t5nNc/fn+/uCZjYZmFTmabP7+7pCiL5Rd8ZtZjsARxMM7V/u/rcsdFTJtKFo3EOATcA3gWuANYRhlA8CewB7A9ea2YH9XHAC4GzCrFUhRB2QS+M2s7cSJuP8wt1/m9r/UeBLpG6qxjzcJ7r7spoLrTAWMmvNjU+XAa9294dSVe4ws4uA3xPW5twd+BTw6ZoKFUJkSu7SuprZDwhjuK8lGFOy/wBCYqpGOkaSHARcY2a5ey/lEtfV3Bk4HDiyxLSTOpuA04HVcdf7zayxdiqFEFmTqx63mR0OnBWfriTcAEw4l2DUDnySMHzwNuDjhBt6pwO/qJHUqhHHrLsdt3b31WZ2OfDvwDjCsMk9/XjZHxB68eUwG/hjP15TCNFHcmXcFHOVLAIOT268xTwgxxNM+wZ3vyDWm29m2xIM/FQGgHGXQToJ17b0w7jdfSmwtJxzysyXLoSoIHkbXjicYM4XlERLvJqwRBnAFSXn/C5uB1tipvQNyWGZqRBC1Jy89binxu2DJfuPST2+ruTYC3E7oSqKaoiZzSAMe0wGrukhx3g6fK+s3rIQor7JW487ucnWXrL/9XG7pJNEUlPidmPVVNWO1xLGjX9MmHTTHUnWQAfuraYoIUS+yJtxJylKd012xMUJdiQY1F86OefouF1cVWW1IZ3r5LSuKsVY9uTL7Hp3X1VNUUKIfJE34/4bIXLkE2Y2LYa5/W/q+OXpymZ2LPAegqnfWDOVVcLdF1J8H0eY2ZmldcxsHCGNbCPF5FZCiEFE3oz7+4RhkjnAEsIklBMJBjU/mSVpZrua2U3AVYQbcxuB72UhuFzMbKaZeSx/76TKB4G18fGPzewnZvY6MzvIzN5PSDy1Tzz+NXfvd0ZCIUR9kSvjjhNOzgLaCDdOxxF64C8B/5aq2ggcRdDfCpzu7otqqbVauPsjwDzCl5YBZxJuyN5O+HJK0sh+DfivjGQKITIkb1EluPtPzewfhNjsKcBjwK/cfUWq2uMEY/sz8BV3f6z2SquHu99iZrsA7ydMbd+R8MviBeDvwA/d/c7sFAohsiR3xg3g7o/TTW7umEZ1SlfHs8Ldj+pFncX0YtEHd38Z+FwsQgjxCrkaKhFCCNEzuexxlxKz5vWml1qtVKtCCJEbctfjNrMGM3uvmd1mZqvNrJ1ws3JzD2VTVpqFEKKW5KrHHXvWV1KcNahMRkIIUUKujJuQmnVefFwA7iNkClxPCIETQohBT96M+11xuwI4xt3/laEWIYTIJXkb455L6Fl/UaYthBCdkzfjHh63t2eqQgghckzejDvJDjgqUxVCCJFj8mbcVxMiSU7IWogQQuSVvBn3V4HlwHvN7A1ZixFCiDySq6gSd18aDftPwFVm9hfgekKK1/W9OP+mKksUQojMyZVxx1mSrzwlxHTP66J6KU7O3o8QQlSDvBld6UxJzZwUQogS8mbcZ2QtQAgh8k6ujNvdf5m1BiGEyDt5iyoRQgjRA7nqcZdiZtOBI4FZhPUnv+HuL5jZNsAOyeLBQggxmMilcZvZDODbhPUW0zcof0VYd/EQ4Ldm9ghwlrvfVnuVQgiRDbkbKjGzPYF7CbMnGwjGXRpdMjPu2xX4m5mdVEuNQgiRJbkybjMbAfwRmABsAL4IHNFJ1SsIve8kdvuXcfhECCEGPLkybuB9wLaEWZKHuvt57n5raSV3f8Ld3wmcTFiybCRwTk2VCiFERuTNuE8i9KIvdPcHeqrs7lcBPyIMm7y+ytqEECIX5M24d4nba8s4589xO6vCWoQQIpfkzbiTPNyryjhnddw2VlaKEELkk7wZ90txW07veb+4XVphLUIIkUvyZtz/JIxX9ypniZmNBT5KGBfXcmdCiEFB3oz7p3F7jJl9vLuKZjabkKt727jromoKE0KIvJCrmZPufoOZXQK8FfiSmb0ZSE9rP9rMXgccBRxDGNd24C/u/tda6xXVY1NbgbsWrWD1hs1sP3Eku04fk7UkIXJDrow78i6giRAauC+wD8GcAb6SqpfMprwZeFutxInq888nl/PdGxeyrrXtlX07TB7FJ4/dhUmjmzJUJkQ+yNtQCe6+yd3fBJwG3B13WyflMeBs4DXuvi4LraLyLFq+nguueayDaQM8sXQdn7v6Edy9izOFGDzksccNgLtfAlxiZuOBPQnT4IcAK4CH3f25LPWJ6nD1A8/TXujcnBcvX8/9z6xi7223qrEqIfJFbo07wd1X0HGcewvMbBgww92frI0qUS2eWNr9j6cnXlon4xaDnlwNlZjZIjN70sx2LuOcA4F1wHXVUyZqxcim7vsSo4bnvq8hRNXJlXED2xFStg4r87whwLSKqxE151U7Tery2NBG49A5E2uoRoh8kkn3xcymADt1U2U/MxvXi6ZGAR+Oj3WDcgDwqp0nc8sTy3jgmdVbHDvj0O0ZO2JoBqqEyBdZ/e7cDPweKO0+JXelflxmew78vZ+aRA4Y2tjA+cfvxrUPv8hNC5ayesNmZk4cyfF7TmevGeOylidELsjEuN19hZl9GLi4Qk0+CPxnhdoSGTO0sYF5c6czb+70rKUIkUsyu9Pj7r81s6dTGgy4idB7/g9gYQ9NFIBW4CV3X1I1oUIIkTMyvUVfusiv2StLS97t7vNrr0gIIfJP3mKrto/b5zNVIYQQOSZXxq0hDyGE6JlcGXcaM9sB2J0Q8jeEYlKpLnH3n1VblxBCZE3ujNvMdiXk5T6gzFMdkHELIQY8uTJuM5tESNM6nl70sIUQYjCSK+MGPkbIAujAv4BfAouA9RQn5wghxKAmb8Z9PMX1I49y97Ye6gshxKAjj0mmAL4h0xZCiM7Jm3Gvj9tnMlUhhBA5Jm/G/XDczslUhRBC5Ji8GfcvCdEkH7DU/HchhBBFcmXc7v5zwko2BwEXxfBAIYQQKXIVVWJmpwN/AOYSVnk/1cweBJZQHP/uCnf3d1ZZohBCZE6ujBv4BR3jtYcCe8fSG2TcQogBT96MG7acMamxbiGESJEr43b3XI25CyFEHsmVcYvBzaRJk2hpaXnlefpxc3Nzh7rNzc0sW7asZtqEyBMybpEbWlpaOph16TEhRCAT4zazWcljd3+qs/19Id2WqD9Ke9U99biFGKxk1eN+Im69RMNC+p4FsLQtUWeUDn2MHDmSlpYWmpubWb++p2hQIQYPWRldd5EiiiIRQohuyMq4P1vmfiGEEJFMjNvdOzXorvYLIYQoorhpIYSoM2TcQghRZygKQ+SHTS2w8AZYugCGjWTH8XC/wreF2AIZt8gHK56CP38UWl5+ZddPj2vktw8O4wf3ZydLiDwi4xbZ4w7XndfBtBPesttQHl9Ve0lC5BmNcYvsee5eWP1sl4dP3FGh/UKkkXGL7FnzXLeHtx5dIx1C1Aky7gpiZjuZ2UYz67r7WF57bzWzG81sRWz3STO70Mx2rET7uWHUlG4Pv6jZ7kJ0QMZdIcysCfgZ0FSBtoaa2e+B3wKvBraK7c4CzgLuM7NT+vs6uWGb/WH0tNQOh40r2W4szBlvjB8OPPuvrNQJkTsyMW4z29bMpvVcsz4wsyHApcAhFWryq0BizNcDJ8S2PwqsApqBX5vZvhV6vWxpaIDXng9NowGH1c/BuqUMHwItm2H6aIM//yc89IeslQqRC7LqcS8GnjWzL5lZXd95MrOJwDUEc61Ee7sD58SnlwPHuPtV7n67u38dOJhg3k3ANyrxmrlg8i7wll/DtgeDNUDTWJ5ZDc+vTSWLvONC2LAqM4lC5IUsh0oM+Dhwq5nNzlBHnzGz44B/Aa+JuwoVaPYDhM9lM3Cuu3dIc+vujwJfiE+PGDC9boAR46BtI4yZDqOn0tJWcrx9Ezz19wyECZEvsjRuJ5j3QcCDZnZeHCeuC8zsUuBqYFvCe/kC8I8KNJ303G9x92e6qHNR6vHJFXjN/NC6rofja2ujQ4gck/XNyW8DrcBw4H+Ax83snWaWta7ecGDcLgZe5+7n9bfBuALQ1Pj0lq7qufsy4LH49Mj+vm6umLRz98cn71IbHULkmKwN8ufAAYThBgO2IURmPGVm55pZniN4nycM9ezs7jdWqM20Ky3soe6iTs6pf3Y/GRq6mNA7YTZsPXBGhoToK1kbN+7+EKH3+iHCTTcDZgBfA54xs6+b2T7ZKeySw9z9q+7eWsE2p6ceP91D3WTWyvh6GmLqkQmzUxEmKSbtDK//MtT3vWwhKkIucpW4ewH4rpn9itCLfT8wGhgDfBj4sJk9DVwB/AG4w91Lb13VlKi50myVetzDYC/paSljgaV9fVEzmwxMKvO06t1Q3v4ImHEQ531lPCOslWdbhnHHD/+vai8nRL2RC+NOcPdVwKfM7EvAGcDZQDJLcDtCr/xDwCYzux+4G3gSeBZ41t3vrLXmCpPuOW/soe6GLs7rC2cD5/ezjcoyZBh/f9ppaWmjuXlY1mqEyBW5Mu4Ed18LfAf4jpkdDZwKHAckc6ObCGPjB6RPI6fvpwzSvfieVrtPjxlUo/dfFa6++upe121vb39l29vz5s2b1yddQtQTuTc6d78OuA7AzA4AjgUOA/YjDKUMJNLDI8N7qJs+XslxdiFEzsm9cadx97uAu5LnZjYD2AHYmjAmXu+kg5RH9lA3Oe6Em7r94QfA78s8Zzbwx36+rhCiD9SVcZcSJ6h0NUmlHlmSerx1D3WT48v6e6PW3ZdS5s3NOs9UIERdU9fGPQBZkHrcU9TGrE7OGVhYA807HcqoHQ/ih4/A5OGw7yTYbiD8thKiH2Rl3GfEbU+xyoMKd38+5vLeBji0q3pmNolitM0/a6Gt1rQ7jD7iXQyZuiMNDcbqTbB6EyxcA6/dGvYpN3hRiAFEJsbt7r/M4nXrhCsI2QFfZ2bT3f35Tuq8M/X4ypqoqgIFhydWw6OrYHM7bDMK5k6A5iHwyAoYNn0nCoWOwTUO/O152HmrUE+IwUjmMyfFFvwIaCeEPP7QzBrTB81sJ+C/49M74w3buqOtAJc9BVcuDsb95Fq4+QX46aOwdAM8vLKbcx0e7ea4EAOdujJuM9vFzJ4ysyez1tJXzGymmXksfy89HlMA/CA+PR64xcxONrODzexc4HZgHCHt6wdqJLvi3L0MFnWS6K+lDa5eAhvauz9/Yw/HhRjI1NuPzWHATHqenFLvfISQt+RNhJVvSlfWaQXe6e731FpYpZj/ctfHlm2EWT3cgJzaXFk9QtQTmRh3P9K2vjJsEFfO6RCTVqX8ITXH3TcDp5jZW4B3A3sTJhu9BNwAfNXdH8lQYr9Zu7n74zNHbmbs0ALW2MCmditmbwcmDYftFVkiBjFZ9bh7uGy7Jeltl8YuZz7l3d2P6kWdxZR84XRT91LCWpYDjvFNoWfdGUPbN3DSyz9l+6EruWz4yWywEQzZtJSNQycwaWQjJ22vJIFicJOV0Rkd+lBisLH3RLju2S33NxY2s2fhEWbai7S3PclOax/jkaG7MWrsLCaPKLB+x7eENSmFGMRkZdyfIGSjG0Ew8OXAX3px3njCDTun4/Jdos7YawK81AIPrOi4f2uW8YER173yfCht7NU2n9lN66EAD69fyMpROyLEYCarOO4LzOxy4P+AVwMTCZNOznL3LiNGzGxPgnHj7md0VU/kHzN4/bah571gFWwuwDYj4dQVVzC2dX2X543a+LyMWwx6MvvN6e5PuvtrCTffVhFWSp9vZp8ojV0WA5cpzXDUdHjdNrDLVuBDRnRbf3NjT7m3hBj4ZD5Y6O4/I6ybeBlh6OSLwL1mdmC3J4oBydKxe3V5rGBDWD5mt9qJESKnZG7cELLTufupwAmERXj3AG4zs+/lfMFgUWGWjt2LFaN26uSI8dTU42hrVAC3ELkw7gR3/xOh9/1DQsTJ+4BHzOykTIWJ2mENPDLj7Twx7UTmL3WWrCpw0+IC82eeyYtb7Z+1OiFyQa6MG8Dd17n72cARwGOEvNOXmdmVhNXfxUDHGnhpq/1431/bOOX3Gzjv5nbWNM/MWpUQuSF3xp3g7rcBewKfJ0zYOZ6wwrsQQgxqcmvcEKZ+u/v5wD7AndRfbhUhhKg4dWGE7v6ImR0KvBdQtMkA5e1vfzutrcV1j5PHra2tnHLKKR3qNjU1cfHFF9dUnxB5oS6MG8DdnXDT8odZaxHVobW1tYNxlx4TQgQyN24z2wE4lRACOJYw/f0O4BJ37yb5pxhoNDU1dXieNuvSY6XPhRhMZGbcZjYE+DbwH2w51n4a8BUzO8/dv1lzcSITkqGPhsJmpq+4ncmr72No23rWD5/G8+MPYsXoXTJWKEQ+yLLHfTFwCl1nCGwGvmZmY939f2qmSmSKFdrY/elfMKZlySv7xq1/knHrn2TR5GN4buLhGaoTIh9kElViZq8B3hyfrgK+ABwO7BS3XwZaCKb+aTObk4FMkQFTVt/bwbTTzFx2PUPb1tVYkRD5I6twwLfF7TLgQHf/jLvf5u5PxO2ngKOAjQTzPjMjnaLGTFo9v8tj5gUmrnmohmqEyCdZGfdBhJzaX3H3JzqrENdTvIhg3KVrLooBypBCF8vi9PK4EIOBrIx767j9Vw/1koz6O1RRi8gRa0Z0n9Wgp+NCDAayMu4k6fLaHuo9F7djq6hF5IgXxh9MwTq/Z75uxNasHjm7xoqEyB9ZGfewuG3voV7yu3h4FbWIHNHSNJkFM05j05CO2XxXN8/k4Rn/lpEqIfJF5hNwhChl5agduXuHjzJu3UKGtq9n/fCprB8+PWtZQuQGGbfIJW6NrBzd2YIKQohcZwcUQgixJTJuIYSoM7IeKjnLzF7s5vjU5IGZfaanxtz9cxVRJYQQOSZr435vL+p43J7fi7oybiHEgCdL4+4quVRf8Z6rCCFE/ZOVcb8qo9cVQoi6JxPjdvebs3hdIYQYCCiqRAgh6gwZtxBC1Bm5MG4zm2Bmh5nZuC6OzzGzH5nZU2a2MW5/bmZzayxVCCEyJ1PjNrNpZvZr4AXgZmD/TuocA9xLWExhO0KCqpnA6cC/zOyDNRMshBA5IMvFgucANwAzKIYGTiqpsy1wOWH9SYA24G/AesISZxOAb5rZWnf/eS10CyFE1mTZ4/4ZsC3BtO8A3k/odaf5JkXTXgbs6+7HuPvJwPbA1fH8r5vZxJqoFkKIjMlqseBjgMMIk2Y+4+6HuPuF7v5cqs7WwAmxjgMfdfcHk+Puvo6wduWLhIUW3lrDtyCEEJmRVY/7xLi91d2/0EWdEwj6DFgK/Lq0gruvB34R65xUcZVCCJFDsjLupLf9q27qHBO3DvzF3bua0n5b3GpNKyHEoCAr454Stws6O2hmBhyR2nVdZ/UiS0vaFEKIAU1Wxp0s/ruhi+N70XGB4Ju6aWtkD20JIcSAIivjXh6347o4/tq4deBhd1/WTVuz4ra7OkIIMWDIyrgXxe1eXRw/IfX4rz20dXzcPtkfQUIIUS9kZdx/JUSC/LuZNaYPmNnuwCGpXZd11YiZ7QnMI/TMr62CTiGEyB1ZGffFwEZgF+CSZPKMme1IMezPgXvd/e7OGohx3pcAjcAm4A/VFi2EEHkgE+N298XAlwi97pOB581sOSHKZI9YbTPwH+nzzGyomR1uZl8E5gM7EQz+a+7+TI3kCyFEpmSWq8TdP29mQ4FPRR3jU4fXA+9w9/tKTtsT+Ht8nOQ3uQr4fBWlCiFErsh0sWB3/4yZXQScQRg2ceB+4Cfu/kInp7xE0bA3EHKZfMbdCzWQK4QQuSDrVd5x94XAf/ey+ovAJ4CFwA3uvqZqwoQQIqdkbtzl4O6bgQuy1iGEEFmSixVwhBBC9B4ZtxBC1BmZDJWYWXe5R/qKu/trqtCuEELkiqzGuI8iRJBUCqtwe0IIkVuyvDlpPVcRQghRSibG7e49jq2b2V6E1d3d3Rt7qC6EEIOGPN+c1NCHEEJ0Qp6NWwghRCfIuIUQos6QcQshRJ0h4xZCiDpDxi2EEHWGjFsIIeoMGbcQQtQZMm4hhKgzskoytW0vqk1P1Z9BL6bIu/vT/dElhBD1QFa5ShbTu5mRSZ3FvaxbVwtDCCFEX1CSKSGEqDOyMu5fZvS6QghR92SVHfCMLF5XCCEGAooqEUKIOkPG3Q/MbISZfcLM7jWzdWa21szmm9lnzGxcP9rdxsy8l+XdFXxLQog6QFEYfcTMpgPXA7uWHNojljPM7Dh3f6QPze/ZX31C9JVCwbntyeXcuGApazZsZvuJIzlu7jRmTRqVtTQRkXH3ATMbAvyRYNoO/BD4PVAATgY+AMwErjKzfdx9TZkvkRj3BuCQHuoqdl1UjELB+ep1j3HrE8tf2ffE0nXc8OhSzn3tDhy10+QM1YkEGXffeA+wX3z8IXf/burYzWZ2O/BbYDbwEeD8MttPjPtBd7+/P0KFKId/LFzewbQTCgXn+39byP4zxzOySbaRNRrj7hsfjNsFwPdLD7r7JcDV8ek5sYdeDolxP9A3eUL0jRsXvNTlsY2bC9y6cEtTF7VHxl0mZrYzsHN8+ht3L3RRNYlV3wp4VRntNwM7xKf390WjEH1lVcvmbo+v7uG4qA0y7vJJjznf0k29W1OPjyyj/d0pfi7qcYuast2E5n4dF7VBxl0+u6QeL+yqkru/SLi5WHpOT+yVNAGsNLMLzOxhM9toZqvM7DYzO9vMhpalWoheMG/udBq6SEYxdexw9p85vraCRKfoLkP5JFkL24AXeqj7POEG5fQe6qVJxrfbgXuAEaljTYQe/yHAu81snrs/X0bbnWJmk4FJZZ42u7+vK/LHTlNHc/ar5vB/Nz/J5vZiHrgpY5r4zLxdaejK1UVNkXGXz1Zxu97de8pwuD5ux5bRfmLcQ4BNwDeBa4A1hGGUDxLixPcGrjWzA929pYz2O+Nsyo98EQOUY3abykGzJvCPJ5axZkMbMyc0c+CsCTTKtHNDLo3bzI4BTicY1SiCzp7+a9zdt6u2NkKvF2BjL+omQyVN3daKmJkBc+PTZcCr3f2hVJU7zOwiQsz4Gwl/n08Bn+5N+0L0lrEjhjJvbjk/FEUtyZ1xm9mFwH+kd/Xy1N7k964ESRRJb14v0d5V5EkH3N1j1Mos4GV3X9BJnU1mdjqwhNCTf7+Zne/u7b15DSFE/ZMr4zazE4D3EkzRgJcJiyisp3bG3BPr4nZ4L+omdVp723gcs+523NrdV5vZ5cC/A+MIwyb39PY1OuEHhF58OcwmzB4VQtSYXBk3wbQhDEP8m7tfkaWYLlgbt72JixoZtyuroGN+6vG29MO43X0psLScc8KojhAiC/IWDrgfoWf9lZyaNoQhCoBhZtZTJEYySNjvyI9OSN+QHFaF9oUQOSVvPe7RcXtDpiq6Jz3uPJtwE3ELzGwaxVC+LcaquzhnBmHYYzJwjbs/20319JdGWb1lITpj0qRJtLQU+wPpx83NHX9gNjc3s2xZp//6ogbkzbifJ2TVy/PkkrtSjw8F7uii3mGpx//sZduvBX4WH58F/F83dZMZnA7c28v2heiSlpaWDmZdekzkh7wNlVwbt8dkqqIb3H0h8GB8erp1Pdj7zrhdCdzcy+b/kXp8WleVzGwH4PXx6fXuvqqX7QvRJc3NzR1Kd8dKj4vakjfj/joh9vkcM9srYy3dcWHczgU+VnrQzN4CHBef/sjdexVVEr8UboxPjzCzMztpexxwCdBI6G1/oSzlQnTBsmXLWL9+/SslMefm5uYO+9evX69hkozJ1VCJuz9pZqcClwK3mtlPCOPdSyjOQuzu/KeqLDHhR4Sc3HsDXzGz3YBfEWY6ngScE+stBr6UPtHMZgKL4tOb3f2okrY/SBh+GQ382MwOJvw91gL7Ah8nRJEAfM3d/4EQYlCRK+M2s8R4CwTjOoeiCfaEU6P34+7tZjaP0DvemTDL8/SSas8B89x9dZltPxLbvoxwA/LMWDpUI/w6+a8+yBdiQLOprcBfH3qBGxYsZfWGzcyaOJLj95zOvttt1fPJdUKujJtwYzJNboOF3f15M9ub0EN+CyGH9jDgKeBKQm94RR/bvsXMdgHeT5javmNs+wXg78AP3f3O/r4HIQYam9sL/M+fHubBZ4v9pX+t38S/lqzk3Ydvzwl7bZ2husqRN+P+bNYCysHdNwIXxNLbcxbTiy8kd38Z+FwsQohecOOCpR1MO80v/rmYo3aczNjmPAet9Y5cGbe715VxCyHyxc2Pdz2loa3duXXhco6bO62GiqpD3qJKhBCiz6xv7T7X2vrWthopqS4Dxrhj9IUQYhCz09TR3R7fsYfj9UKuhkoSzGxX4ERC2NtwOv+CGUK4YTeGsDTYNHL6foQQteH4udO5YcFLtLVvmUx01qSR7LlNOWua5JfcGZ2ZnUdYjaWciBIjP2lfhRAZse2EZj75hl34zo1PsHpDcUX6naeO5hNv2HnAZLXMlXGb2WGEyJIkH3dPJGZ9A/CXaukSQtQPB2w/np+fsT/3LF4Z4rgnjWTHKQNjiCQhV8ZNmI0IwZA/Rpja/RIhfnkrYAZhBuEM4G2ECShDgSfd/Vu1FiuEyCdDGxs4ePaErGVUjbzdnDyUYNq/cfdvuPvzcUmuWwhaj3b39e7+qLufD5wSz/sPMzsoI81CCFFT8mbck+O2dEmsewlDJ4ekd7r71YRZiga8r9rihBiQrF8O914Ef/sS3P0TWPNC1opED+RtqCRZyeWZkv3JQgRz2ZJLCYmdDqiWKCEGLItugRs/D+2bivvuuxgOOzc7TaJH8tbjTnJ7lMbsLIzbXTs5J1lKbGAkIRCiVrSs2NK0AbwAt36TmQMjcm5Akjfjfixu9yvZ/2TcjjazOSXHkpRfvVl1XQiR8Pg1W5p2ghc4foeBETo3EMnbUMn1wJHAh8zsT+7+EIC7t5jZk8As4FTgf1PnnBS3y2uqVIg65Oqrr37l8Zznb2Lq6q6zDk8bGaJt29vbO5zXHfPmzeufQNEr8tbj/j9gHSEP9d1mdmnq2BWEm5DnmdknzewNZvZt4N2ESJTerusohABah43r8lhbAV4evSNjDnozTTscTA8pQESNyVWP291fNrN3Ar8DmoBXpQ5/Dfh3wtBIerkuA9qBb9RKpxADgZfG7sO2y27CvNBh//rNsGazcdOY42necTwNDcb/PQKnzILpIzMS2wuSVepHDIHTdjOOnllgq+HG4y8XuGJhIzcuLk6urvdV6vPW48bdryCMcV8BPJLav5SwjuPzBLNOykbgfe6uHrcQZbBp6Bgen/4m3Io2sKk9mPalDW9gUVtxAsuGdrhiEbQXOmspH7S0tNDe2sI3Xl3g9N1h6qgGmoYYe0xp5DOHwtt33vzKSvb1vmp9rnrcCe7+APCmTvbfGVc4nwfMBl4G/uzuCjwVtWP1s7DgT7ByCYycADsdC1N261infTM8fQdsXAUT5sDkXTKR2hPLxu7J2hEzmLLqHkZsepl/rRnL5Y378YJNphgTEFjXBo+vhl1yugJYc3Mzr9/J2H1K+CIqFIrfMg0NDbxn3+Fcu2QoyzdQ96vU59K4uyOuOnNZ1jrEIOXJm+CmL0Ihldd5wdWwz+mwf1wadMk/4eYLYMPKYp0pu8PrPgsjJ9ZWby/YOGw8SyYfDcDvHoOl3QSTrGqtkag+sGzZMvjjB+DFBwG49777KBQKNDQ0sM/eewOw5JxzYI9TumumLsjdUEkpZjbdzA41s+PMbHTc12xm9f2VKeqPDSvD7MJCJ8n4770IXngAViyC6z/T0bQBXnoIrvkEeL6TWI7pYVWv0cO6P545mzd0f7xtY210VJlcGreZNZrZB83sMcIsyluAq4DtY5XXAsvM7PtmNnAzyYh88cQNXcc9Azz6Z3jo8jBM0hnLn4Dn7q2Otgoxt5uraUQj7DSuZlL6xtQ9+ne8TsidcZvZeIJRf5OwcnpyEzLNTGAEcBZwv5ntVEuNYpCyvuv1DAFYtxSWLui+ztJHuj+eMTuMhX06Gc0ZajBvOxiaO8coYfc3wZAu5uJN2R2m7VlbPVUijx/DZcDBBLO+DfhMJ3UeBZ6KdbYGrjazETVTKAYnY7fp+fiwHuLlejqeA163DZw2BzY+eTetzz5My8M38e5dYNaYrJX1gnEz4A1fhtElCwLPOBCO+WI2mqpArozbzN4MHEWYUPNedz/c3b9QWs/drwN2BL4ed80C3lsrnWKQMue1MGxU58esAXY9IdRJaNsIm9ZBIQ6dNA6FWUdVXWYlmDEK1t91GStv+gkb5l/LmLyPbaeZvje89Tecc107/3XDRk77YzscewGMGJe1soqRK+MG3hG3l7r7j7ur6O4Fd/8Y8GdCz3uL8EEhKsqwkXD057bsNVsDHPZhmLgD7PQG2Go7WLkIVi2BNc+FG5ZrnoO93wHN4zORPuhoaOD+l+CmRe08syZrMZUnb+GA+xF625f2VDHFzwgTc3brqaIQ/WbrfeFtl8AT18HKxSG8b8c3wOgpTJo0iXENLfxiXiPTR8OYYc6QBmNjm7OqdT23feEDnHPd2UD9z9wT2ZI34066I8+Wcc5zcZv/wUMxMBg+ptNY4JaWFk7fs42mRuPlFni5BYrLojpzJ8HsURt5cGmOpx+KuiBvQyWr4nZyd5VKmB23K7qtJUSVaW5u5oBthtHQ0EBDQ8dLK9l3yMwRNDc31/3MPZEteetxzwdeAxwP/LWX55xN6NbMr5YoIXrDsmXL4Mr3h8k2dD5z74fvP4cfDoCZeyJb8tbj/h3hRuOZZnZEdxXNbIiZ/ZCwwDDA5dUWJ0SPzH5V18esAbY/snZaxIAlb8b9c+BBYChwjZl9zcyOTR3fysz2NbOPAA8D74n7FxJuUgqRLTsdG5JKdcZep8GoSbXVIwYkuTJud28nZP5bTFiK7FzgTxTv8NwE3AVcAMwh9M6XAm90904SSAhRY4Y1w/Hfhj3fxsqNTsFh4QqHoz4JB7yn5/OF6AV5G+PG3Z8xs70I5nw6Xa8l2U4YHvmwu79YI3lC9EzTKDjoLN74mo/Q0tJCc3Mzp/7u9Vmr6pG3v/3ttLYW0/8lj1tbWznllI7j8k1NTVx88cU11Qf0egk1CEuuJduBtvRa7owbwN3XAGeZ2ScJMyn3AiYQ9K4gLLBwk7s/n5VGIQYara2tHYx7j8kNjBtuLFxR4IV1Oc7nOgjJpXEnuPtKwko4V2StRYiBTlNTEwBzJxufOrSRbUaH/e5w23POF29tZ+2mjnVFNuTauIUQtePiiy9m+KaX2fup79NY6Ji+docd4KRjZ/LgzHdnpK66bGiDxm4WkMgbmRi3mVXlpqi7a0qaEP1g+oo7tjDthLEtixnd8jRrm7etsaq+M6nZOHqHBrZe/g/WNG+3hfYFK+GOl2DpxhCpMd8e5vSDt2PWpC6SieWErHrcXWSa7xeOfkEI0S/Gtizu4fiiujHud+/VwDt2b2ZIo7H90msBWD1yex7Z5jTaG0fwwHK4JpVcowD8a8lKHnl+DRecMpeZE/ObRSOrcECrUhFC9IN2637tskIPx2vN6JYlzFj2N7ZZfgvDW5e/sn/yqvs5Y89GGkscbuz6Rez4/BW0F+AfXcSibdjcziV3P1NF1f0nyx5qEpvdSojPVqo0ITJm+ZjdGbPh6S6OGsvH5CMJZ2P7RnZ99mLGrl/0yr6ZS6/nha0O4Mmp85i+4p+81MW5E9YuYNnqNaxv63pliDueehl3xyyf/cGsjPtRYOf4eDhwNGG5ssuAP7h7D2tECSGqwYtb7ceU1fcxcuMLWxx7ZuIRbBo6NgNVWzLnhas6mHbAmbbyTjY0TWTUxhe6NG5whrWuALo2bnfHHXLq29kMlbj7rsCuwPmEKe5DgFcD3weeM7Mbzex9ZjY1C31CDFYKDcOYv92ZPDvhcDYNGQUY64dP4/Hpb2LJ5NdlLQ+AoW1rmbj2oS6PT19xO5uHdJ99ceqoRoZ14357zhhHQ0NOXZsMh0rc/VHg88DnzWwOcCpwCmGyzasIE2++Y2b/BH4PXO7uW3YDhBAVpb1xOIunHMPiKcdkLaVTmluXYd0EkA3ftIJntjoYG7YwDHUUOmbD2DBsAptGbsP+k+G2Tsa5GxqMU/ebUWnZFSUXuUrcfaG7/6+770PIr/0J4B6gETgc+DbwjJn9w8w+aGZbZyhXCJEhm4d0H6q3ttDEf790FEuG70DjyK1oHD2J5Ruh3aG9YRhPTD8JzDh0ChwxDUY0Fs/dZqsRnHfcLuy+dT6GhLoid+Fz7r6IkKfkAjObQeiFnwIcREjhegjwDTO7i5AG9nJ3z/ctYCFExWhpmsy6EVszasNzWxzbXIA/b9qLlQ0j+ejLJ3DM8AUc2fwUE4YP57b27dh6zsFsbJoIhPHrg6fA/pNg2QYY0gBnvHnfWr+dPpE7404TDfmbwDfNbDphQeA3EQz8IOBA4Otmdjfwe3f/epeNCVEByklyZOOm0zS+CV+/fMAlOcqaJ6adyB5LfsaQ9g0d9j/tk7nCXgPARh/GFevn8scNezJr3GwowJtbYVbJbP0hDTAtvyHbnZJr404TE0p9F/iumU0krJLzPsICwwcA+wMybpE5z62Ha56GsUd/gELBMYMrF8HrZ8Dwurni8s364dO4d9YHmLbyLsatf5KCDeHlMbvxtaX70NLedR6Vl1pgVtfBJHVD3f0bmdnuwInACcA+hHjw/N7+FYOKVa3wuydhU+remZnx2GrY0A5v62KNBVE+m4aOZcnk17GEEO3y9re/naZD38nQKWEZ2kLBX9k+9eSTAMy/5Aq+vvCOzNLSVoq6MG4zOwg4GTgJmJXsjtsWwvqUWrpMZM6/lnU07TRPrwu98a3r7Gd5vdDa2srmJ+5k7KRZWxwrFBwvtLH2ibvxTfWfojaXxm1mjYS47pMIPesknjsx69XA1QSzvsbdN9ZcpBCd8Mz67o8/vU7GXS2amppoffZB2l7cl2HTd3qlxw3QYMa6+65mmLVDU1Pdp6XNjXGb2Qjg9QSzngck8TiJWS8H/kgw6xvdvRqJqoToF0N6GLTr6bjoO8nQR7vDgy/DQyugpQ0mj4B9J8GMUz6ascLKkalxm9k4wk3GkwnT3pNlypJ/7+eAKwlmfYvStoq8s9M4eK6l82MN8bioLo0Ge00MZaCSVT7uswk96yOihnQ/ZBHwB0J89h0ZyBOiz+w5AR5aCUs3bHnsgMkwZljtNYmBR1Y97u/RMRrkEYpm/UBGmoToN8Ma4W2z4Y6lsPDhdTCsmbaVL3LsjNnsMSFrdWKgkPUYtwNLCQsAHwUc1Y80iu7uR1ZIlxB9ZvgQOGo6fO+Kz9Pa2kpTUxN7vPuwrGWJAUTWxg0wOZb+YBTzewshxIAmK+N+GhmtEEL0iUyM291nZvG6QggxEMjDUIkQ9YsXGNreQlvDcLyh58tp8Vp4eAVsbIdpzTB3AozK1zKOog6QcQvRB6zQxrbLb2LqynsY2t5Ce8Mwlo7dk8WTj6a9ccQW9d3hL0+HUMGEhWvgrqXw5tmaTSnKIxcLKQhRb+zy3CXMWH4LQ9vDbJvGwiamrbybPZb8HCtZcQWCYadNO6G1AH9cDAXd8RFlIOMWokzGtCxm/NpHOz02auPzTFozf4v9Dyzvur21m+GpNZVSJwYDMm4hymTC2gW9Ot641dYMm74zDSPGsLaHzDprlHlHlIHGuIUok64Wql20aBGFQoHbl7/I5wszGXP0OeGAO/c/tig+DueaNbD9rO1fOXd8fSerEzVGPW4hymTlyM5XQygUCrg1MH/CMTSMmVI8YIY3DqNh5HgKBQ+5oVPmP74Jtut+/VshOiDjFqJMVo7agbUjZmyxv6GhgWUNk7ipbXcaGjqmbrC2jdC+mYZhw2loMMzCpTdqCJw4MyxcK0Rv0VCJEOViDTy07enMfvFqJq15MA6dGOP2OIbvr38j01rHdXnq3hNCIqrWdpjaDLtuBUPVfRJlIuMWog+0N47g8a3fzFNTjmX45pVsGjKaTUPHsu7x7s/bqgn2729mHjHokXEL0Q/ahoxk3ZDi7BktpCBqgX6kCVFB9pwAk4Z3fmy/SVpIQVQGGbcQFWRYI7xtTjDpEY1h38Th8Ppt4FVbZ6tNDBw0VNIP4gLHHwJOBXYkpKpdBFwGfMfdV/Wz/bcC7wH2BpoJa3BeB3zT3XsYTRVZMWIIvGbrUAoODYoYERVGxt1HzGw6cD2wa8mhPWI5w8yOc/dH+tD2UOA3wCklh2YBZwGnm9k73f2y8pWLWiLTFtVAQyV9wMyGAH8kmLYDFwKvJiy/9h2gAMwErjKzMX14ia9SNO3rgROAQ4CPAqsIve9fm9m+fX0PQoj6RT3uvvEeYL/4+EPu/t3UsZvN7Hbgt8Bs4CPA+b1t2Mx2B+JcaS4H3uzuSe64283sz8DtwDjgG4DW2RRikKEed9/4YNwuAL5fetDdLwGujk/PiT303vIBwueyGTg3ZdpJ248CX4hPj1CvW4jBh4y7TMxsZ2Dn+PQ37l1kHIJfxu1WwKvKeIkT4vYWd3+mizoXpR6fXEbbQogBgIy7fA5JPb6lm3q3ph73ajjDzGYBU3tq292XAY+V07YQYuAg4y6fXVKPF3ZVyd1fBDZ0ck6/244sKrNtIcQAQTcny2d63LYBL/RQ93nCDcrpPdQrbRvg6R7qPhe3482syd1be/kaW2Bmk4FJZZ62c/rJwoU9fc8ElixZUubLlMfDDz9c1falv3sGi/5O/t9rOidWxl0+W8Xt+tIbh52wPm7Hltk2wLpetp20v7SXr9EZZ1NG5EtnnHjiif05XYh6ZwZwX61eTEMl5ZOsVbKxF3WToZLerm+SrtdT+xtSj7V+ihCDCBl3+SRRJL1ZlzuZN9dV5ElXbfem/fScvN62L4QYAGiopHySIYwucsB1IKnT2/Hn9PBIT+2nj/d5fDvyA+D3ZZ4zijAJaQ2wGngG2NRPHWlmE2anJpwAPFnB9quN9GdLtfUPIwyPJNxcwbZ7RMZdPmvjtrkXdZNEzSvLbDt9bk9tO2EafJ9x96X0bYz8zv68bnfYlmt5Penu1b3zVUGkP1tqpL9mY9qlaKikfJLb2sPMrKdIjCRK5Pky2wboKQlocnyZu7f1sn0hxABAxl0+C1KPZ3dVycymASM6OaffbUdmldm2EGKAIOMun7tSjw/tpt5hqcf/7E3D7v488GxPbcee/o7ltC2EGDjIuMvE3RcCD8anp1sng2mRd8btSsq7cXFF3L4u5vzurm2AK8toWwgxAJBx940L43Yu8LHSg2b2FuC4+PRHZc5q/BHQTojN/qGZNZa0vRPw3/Hpne5+F0KIQYWMu2/8iOId5a+Y2S/N7LVmdoSZfRO4OB5bDHwpfaKZzTQzj+XvpQ27+0OE8DyA44FbzOxkMzvYzM6lmIt7MyEFrBBikKFwwD7g7u1mNg+4kZCz4/RY0jwHzHP31X14iY8QIlLeRMhGeEjJ8Vbgne5+Tx/aFkLUOepx95F4I3Fv4L+Aewkx2K2EKI8vAXP7Gjfq7pvd/RTgrcANwMuEHvazwC+Afdz90v6+ByFEfaIedz9w943ABbH09pzFdJyu3l3dSwEZtBCiA+pxCyFEnSHjFkKIOkNDJSKvLAM+W/K8npD+bKl3/d1iPa8FIIQQIk9oqEQIIeoMGbcQQtQZMm4hhKgzZNxCCFFnyLiFEKLOkHELIUSdIeMWQog6Q8YthBB1hoxbCCHqDBm3EELUGTJuIYSoM2TcQghRZ8i4hRCizpBxCyFEnSHjFkKIOkPGLSqCmZ1rZgdkrUOIwYCMW/QbMzsP+DrwezPbK2M5Qgx4ZNyiX5jZaOBVQCswA/iDme2drSohBjYybtEv3H0t8A7gaqAFmAlcbmb7ZKmrt5jZmWZ2iJnV5bUg/YMT/bFEvzCzBnd/DvgwcC1F874s7+ZtZm8Ffgx8Edin3sxD+gcv+kOJfuHuhZR5f5D6Mu9FwGbgSOBrwH51Zh7SP0jRH0n0mzo275eB1UABOAK4gPoyD+kfpOgPJCpCnZr3SmAd0A44wTy+RP2Yh/QPUvTHERWjDs17K2AiMBS4m/Cz/VXUj3lI/yBFfxjRL8zMUo8bo3k31ol5DwU2xceXAt8C2qgf85D+QYr+KKJPJIbt7p7sc/f2km0uzTv1ZbMnMB5YBlwOfBv4Jjk3D+kX+oOIsok9ajezsWa2v5mdb2Y/MLOfm9l7zeywpG4ezTv1ZbNL3G4ECu7+PCE8LdfmIf0Cd1dR6XUBGuN2G+AqYDEhKiBdXgR+CowBGlL1LyfcjCoATwH7ZPg+RgLfj1p+VHJsDvAVws/4AnAjcEDyXvJQpH9wl8wFqNRPSZnwtsAT8aJ6GXgWuBO4N+7bHLd/BfYChsbzts6ZeR8O/An4UPr9xce5Nw/pH7wlcwEq9VWAyYQIgAJwF3ASMDV1/F2Ennhr6oLbPXW81LyfyNi89wGsi2O5Nw/pH5wlcwEq9VEI90MM+O94ES0ADk4uOuIQSny8FyFb4MZY92pgSOp4Yt7J5Iv1wIwavx8red6pGeTVPKR/cJu3BvxFr3D3gocr6Yi46zbggbgPd29PRZrcD/wI+Alh2ORYQqhX0tZzwDmEixDg1+7+TA3exiskulPPC13UW0gOb5hJ/yC/YZn1N4dKfRRCb3t7wmy3AvC2uL/Lng9wEHAr4YJbAOwa9yc3OGcA707Vb6yG7gq1U9rzu47Q86u4ZunPj/68liEV/h4QA4QY8tduZsPcfZO7u5kNJUSKAEyFrntK8dgdZvYH4BBgJ2AH4JHYboOHXvZP0q/XT83m8SoveQ/jCF86hwFLgafd/fZy2nb3hWb24/j0HOC1wCTgRGBJf3Sn9Ep/F9RCf12R9TeHSv4KxR7xdODTwLbx+TRC+F8r8K103U7aSMa+twNeIvSSzk0fq5DWM4FTO3nddK/+Gophi0mP7b+A7frwerMJyZAKwK+kv77112vJXIBKvgodQ/4ejhfIZ1LHr4z7NgB7p8/por0xwNPxnHfHfZX6+fzt2O5twBtT+4fE7QyKYYuJabSnnn8X2K0Pr7szcELp30z660t/PZfMBajkrxCS/yQx2Y/GHs64eOxNhPjrAnAfsEPc31VUwLzUhXpIhXVelfoS+UfJxTyeMA5aAP4F/C/weuDLwE0pTf8H7NEPDX0eY5X+bPXXc8lcgEo+CsWftg2E8euVwCPA/sDIVL1JwG8I0SJtwA3Ajsm5JW2OIkQAtAF/Scy/Alot9fiXFCf9/AM4Ke4/GFgDPBgfD4/7hwH7Ab+vlHlIf33pHwglcwEqNfywu+4VJ+OSM4FfE2K124Cfl1ykyTDKLMJMyeTCW0C4ATkiVXc28J+pOmdV+L2k48J/nTKPWwk3rt4f930DaC55nwbsCvwuQ/OT/gz113vJXIBKDT/sbia5EH7aPhkvouSn7E87qZeY947A7YSfyQXCZJurCOOWPwH+mboov5w6v5I3Jrsyj2sJvb/NwNHdnJ9X85N+le7/9lkLUKnRBx3C8RYB81L70rMdGwmTZlalLqS/Eaez03nPe/t4sT0S67elzt0MLAfO7+z1Kvi+OjOPTRRvhJ3U3Wv3ZB6V/KKR/vzpr9eSuQCVGn3QsHe8MJYBx6T2TwMOjY8bCT9tX05dSB9P1e3MvMcRJtr8hDATcgFwP/BJ4NhU/ardROrCPJIvkW90pr/k/FLzuBCYW8PPRvoz1F+PJXMBKjX6oOFAikMhKwk3HacT4mtXAPvFeo2EPCNLUxfSO1PtWGePU+eOJI5ppvZXPVyrxDx+k9K+Eji9K82p/aXmcQkwvYafj/RnqL/eSuYCVGr0QYcIj7cRQvgKhARPycSYWwlj1klkSal5twCnpdrqVYKgKr+fhtLXBoalHl+c6vndDZzYlf7U/l0pRjtcKP0DV3+9l8wFqNTwww7DGm+gOLGmEB/vlqqTDgv8GmGcOjHvt6fqZTLpIWUYic5pJb294anH6Rtmt9Axzrgr89gTOLOnetJfn/oHSslcgEqNPuhiKNbudBzDXkNqzDvWyaV5l5jF/wLXE4Z5/gqcl6o3NPU4bR6lk0S6NYW0IUl//esfSCVzASo1/LDDqto/jhfSE4RZkck45OtL6ubKvOk4Ff/B1BdPulwCjI310r3AX/XVPKR/YOgfaCVzASo1/sDDT9EvEZaNeivFMe+VlNfzPi0D7VOB+VHDovj4z4TVdJIohsuASbF+2jwu6sY8avXlI/0Z6h9IJXMBKlX+gIs9pfSstRHxcTMdb1iWa97vrtV7iOXL8bXvAo4CJsbjpxDWLkzWury8C/NI9/xuJnXDTPoHrv6BWDIXoFKlD3bLvCGWNuHU/pGU1/P+KsUx8m9UUO8Wcd6UjHECf4+vfWwn7+9wijlUemMeyeIOFYkXlv5s9Q+2krkAlSp8qEWjnQJ8FLgCuIMQonV4J/V7Y95JkqBGwmSbr1RI63sp5vtOz+RMfiHMAb4IbEP4ST4fGJ+ql/4SOqRM8/iO9Ne3/sFaMhegUuEPtGja20YjTnKJpMvZxJtIqfM6M+/XxWNDgLOAU+JzK329Pmr9THytexLzKDk+jhC1UCBkoVsK/DYeS0cupPX01jzSk0L6mo9a+jPUP5hL5gJUKvhhFntJ0ynGai8FniEs75Ssqt5G6ImPLzm/1LxXEabAf5+i6e9T+np91Dolvk5rbPeCTupMJCwSm/7SeZBi77/TWZw9mMfQktfoU8iZ9Gerf7CXzAWoVPgDhSZCop7kJtKxwGRgLPARwiy2xLw/DkwoOX8kcCqhF5bUS1Yt+WqFtc4l5DW5tJNjyZfQZEIUTGIc64CPAU3pep08LjWPPxATZkn/wNA/mEvmAlQq8CF2HJscSlhx5AXgyJKLqRk4mbDUVHfmPQI4mjAVvkDosX+hs9frh+ZX1qQs2T8dmFZSp9Q8biWsrNOYrtfJ40MIE0CSULWKRcFIf7b6B3vJXIBKPz/AYrjfdoRhjiPjRXIjHVeuSS7CJuCNKVPu1Lxj3SHAvsTlyeK+amb5m0kYnlkIbF2iezJhtl5iHtcAr+ulefyVGtwIk/5s9Q+mkrkAlQp8iGFiRHJBnQesB74bj3UW5tWjeXdxXtVuIhF+KVyeeh93AdvEY2nz+GIfzGP71OOqjKlKf7b6B1vJXIBKBT7E8PP2/ngxbYzbX6eOb3ETsTfmncH72JuQ/yKJhKmYecTnVY1ekP5s9Q+mkrkAlQp9kDCDMJU4uaBaSa1208U5pebdGs17Ugb6E2PYg7DyTkXMQ/oHh/7BVjIXoFLBD7OjebcCfwT27+GcxLxvSV2Ir62B1vTEjUbCzM5kvL4780jqdGse0j+w9Q/2krkAlTI/sE6mssdt0tuZkTLhdcBPgb16aDMx70eBi2qpv5t6fen5XQscX80en/Rnq18lfj5ZC1Dp5Qe1pWFPINyUHN1J3RmEJD4FYG0vzXs4HRdprcbCvsmXy0TgNOB7hFSgfwDeDuxdUn9uL83j8ynzeH8VPwPpz1C/SuqzyVqASi8+pOIFN5kw4/Fq4GlCj3oBIcf2sXQM/+vMvPfu5etV/CZS6j1sC9xJcRZnUtqAx0gl5I/1e2se3wF+WIPPQPoz0K9S8nlmLUClhw+o4xDI7RRnoqUvuELKwNMJgMrueVf5vWwdzaFAmCC0IL6n+SXv6SfArJQx9MY8RqVep1ohc9KfoX6V1GeZtQCVbj6c4kWxDcUV2h8Dfgu8hzAh4uclJn49MU9yPDdz8ybc+GqiuPrOnYQbWeNSdT4MXJd6L78GpqSOl4653pGYR2d/M+kfOPpVOvlMsxag0sMHFDK0/SV1wR2W7t3EOsdEQ096438BtkodT5v3KuBnwAE10N4hPIww1X4dYVgn2ZfOGHdY1JaYx/dK2kvMY108vihtPtI/sPSrdPPZZi1ApYsPpnhhHQk8C7wEHJfaX7qyzYGEn7AFwuo03yKudBOPzwBuSl2Ub6+y/kTfzKjlzcBzwIvA9M7ea3y8D8UFZgvAWSV19yAk9C8AP5L+galfpYfPN2sBKj18QGHFmQJhqKTbiTGEvBDPx/r3AwfG/clFPAN4CLi4RtonpXpnvyLcEHuwi7pp8ziakJ2wjTANexRh9Z3kS2of4EOp+lWZkSf92epX6eazzVrAYCvd/ZOXXDzJRfK7eOFd29P58fibUr2lznIsj0s9rmbI35B4sV9FcULQJkICoymdvZeS9/+/qfP27KxONd6D9GerX6V3pQFRU9y9YGaTzGy2mR1oZvuYWZOZNbq7m5nFqsln0xS3w1PnW2m7FiFkBUwm4OzQyeuvivUb3L29ku8tvod2M5sOnE+4qXoSYeXvoQQzmUUYY8XdCyXa3MyS9/1TwiIQQ4Gd03VKzqnYe5D+bPWL3iPjriFm9loz+xwh/OpuQijWPYQbh9eZ2XbJhZG6IO6K28PN7Ph4zCnBi6wijIc3ANPMbGhnWkov2v6SfBGY2baEL4//Bt4Z38dpwO9T1T9iZgf3oGsN4ac6hKiIqiL9gaz0i/KQcdcIM/sQ4abPJwnLRjWlDh8AvIpUzybVq74bWEa4iN5hZjv14uWSi281sIXJV4P4S2ArwpjoTsBTwFQzG+PubYSZeYl5zALeb2Zzk/M7+RWxDyHrIYTohaoi/dnqF+UxJGsBgwEz+xTwhfj0dkIM7D8IEyKmAacT8jw8mJyT6lVfH+sfT4i9fdTMfuLuT8e2Ld0Dj8a+K8G8/xAv2qqR+nluhP+nOYTQxNOAhe6+1syGuHubmZ0Wdb0FOAFwM/uBu99e0uZI4OBY9x+EG7PSPwD1iz6S9SD7QC/AByjeLPwfYOdO6uwIzOpkf3qqexLq9zIhsc/undQfSZgSn8TYHlzl95bcQJ0F/ImwFuEKQv6LRjre7EreSyNhLcIkbHE+cAYxRI1gGJ9M/c0+LP0DU79KPz77rAUM5EIY/ng6XgCfI5UQitA76nGWWeriPAi4N7a1mtBDfxshYZAReuOfS11wn6zRexxHcVX4f8XtL7qomzaP36a0Fgix6gvpOKX/K6V/B+kfWPpV+vi5Zy1gIBaKcdPnE0KwrgF27GebjYSx8KTnnay8vpiQjjW974JSLVV8ryOBb6S+oAqEcfktflkk7yP1fi4pMY81hGGhH5NaWJYqhpxJf7b6Vfr4uWctYKAWwvj18nhBfKIP55fGy46M21GEyRSPlFx0BcKU5v9MnVPVC47ir4FRhNSeS1JaPtbNeZ2Zx2bgAUpmdFbZ9KQ/Q/0q/fjssxYwUAtwFCGp0zJgZtxX9s9Nwiy2zxB67+m8ErsB70odewOwXep4rVZSSUwgMY9FKfN4S/K+O/ki6so8HgWOT9Wr6k906c9Wv0ofP/esBQzUAnyI4s3EGX1sYwYh93aBEJs9vpfnVSXDXC/qJOaRZDIsAG9Mzu/BPJIx182EXxPHVfL9SH+2+lUqWzIXMFAL8B/xQlgP7N+Xf35gDPBtwtj1ZuA/4v4OU+Op/jh2coEPJ9wMfRNhqbNjKYmGAUYTQh9zYx7SL/MeaCVzAQO1AO+meMPwpLivL+Z9ZOoC/E4G7yO5sKcBPyDM9Ewu7gLwMPALYHzy/qJ5VKLn90Dyt5P++tSvUp2SuYCBWgix2c/Gi+AGesjs10UbDfEiTMK8Lqnxe0gvd/Vw1LCRkNM7WTVlU9zeCrwaaI7nVOJne4EQ5TBK+utPv0r1SuYCBmohxNcmF9cS4M19bKeJkKK1APyphvqT3ts0QirYJE74bGAXwmLFrwcuJNyELRCWwppNMRyyM/M4Pt1+yWumzePPhGGmXaW//vSrVLdkLmAgltRF92ZCD6lASD6fpNPsdcQHYaJO0uP+Qbr9GryPEYSf54lp7AcMLalzYcoUPkuMfOnBPE7t5jXT5rHF0ljSXz/6VapXMhcwkAuwPfAHij9nr0od65X5RvMvEJJMzaux/vGEyRwbgbMIKT4bUsc/nzKDjxOT7RNWXUkvHjsqmkpiHhf28LoVCWWU/mz1q1SvZC5goBdCPuRFFFdjv4LU1Pcezp0Ye1zthFSdMyukqbdfGidGzWuBXUqOpU3jP4GmuP9YwkSgP5fUH0VYQut/pL++9atkXzIXMFBL+uIE3keYTpxMSb8WmEtxNmTys3Zo6pzRwAdTF+f7K6mpk2NbhBQSMhIWCLkw0utXpk3jI8kxYCvg4rj/MWBqyftrSrUxRPrrS79KforSulYJ97CaiLsX3P1CM2sEziP8/H0d8Evgz2b2c3d/Mp6zGcDMDifMhPxEbO6r7v79eKxDGtcyaTaznYFDCRf5YuBxd7/dQz7nRu+4osm6uJ1DyM98m5l9npCkH0I2uh+4+4aof6WZPRSP7UDIC/2ix+T87t6aeg99STcr/dnqFzlBxl1F4sWYmPf3zGwVYWLO/sCehAvyPWb2B8IyUSMIy029irDYAgTT/i8o5l7uixYzOxmYR8jTPCx16Gkzu8fdT+mk7ZdjmQAcFds4Nx7rYBpmNjR+8TxO+Gk/gjA229nfpewvHunPVr/IGVl3+QdDoeMNpYOA/yKkZk2GTpJtujwInJs6r883jAjT79enXmcR8EzJ63Wauxv4SSfaXvl5Huukh4XeT/Fm6hY5w6W//vSr5K9kLmCwFGBCyfNdCAmi/kqI824h5EOeT1h84aBU3f6Y9idSF/zPgXcQblRNJfT+ryTEBo8uOS8ZF92HYirZAmGWXmOqXnrcdAJwUax3FWGcvl+hi9KfrX6VfJbMBdRjKfdiIMyi/A3w+S6OTyGEDo4Etio51uc8JBTzpRSATxNvVpXU2QYY000bwwlRC4tSbX2YkokZhPUJP5yq894K/J2lP0P9KvktmQuop0JxckOy7bEnHE37ungxPUqq553qVb0yBblSPSRC1MoD8XW/yJY9us5mznX6JUGYBfp54ImUMdxLCFX8KPC/hMUikmMXdPc60p9//Sr5LpkLqLdCuFv/CDA3Pu/SvKNpXxsvptXAHjXQl3wB/DthnPNOuhjrLL2oCb2/2YQc4EcTftInIYtbEXp0/6TjuHxL6vkm4Iup9soe4pH+bPWr1EfJXEC9FEKveATF3vNLwG7x2BYXSIlpr0pMu9oXE2H2XBNwS3ztb3RRLz1OuhPhhtdjwEqKPbf7Cb3F2bFeM2Fs/kKKy6WtBzYQxm/f3ln70l8/+lXqo2QuIO+Fkp+vhHCuJHfIis7MG9iV4k/Xmpl2ic5/EnphSQ7vpi7qfQa4OWUWSUl6dC8Df6Rk1iZhXH53Qkjjdt39zaS//vSr5LtkLiDPhWLCnSnAJwk3GBfEHk6Sf2Ql8acwITFPA2EpseRYTU2b8MugmWJmwl+kjg2N2/HAIYRp9KWG8Syhp3cDYXinQFg7838IPcmkjdKf+RUZo5f+bPWr1EfJXEBeS8q0t01dhIVo2hsp5h7p0POO50wiTD3eJ91WjXQnF/CXosbngI/GfUMJk38uomO2uAIh3/PPCculjYkGtGuq3t3AMOkf2PpV6qNkLiDPhZALeUG8cO4h5A7ZlTDr8TWE4ZBkZe1S825MbzPQPi9+yRQormZyN51P9vkhIQl/ejHiJHJmf8LU6wJwjPQPDv0q+S6ZC8hjib2docDX4wUzP15Ayc/UJIxvKnAOxbCvFaSGTXLwPj6cMo4C0Joyi7WEsdV3lZyTnuXZSBgmSmb5nSz9g0e/Sn6LcpV0gru7mbUB+8ZddwEPeUwC5SEHibn7i2b2C8INyI8DuwG3mNkR7v5Qf3KL9IdUfpRvmdk64G2EL54RhDHUlYQFZR9w98fjOeaBQup5O/BSzLEyhfA+pX+A6xf5JxmPEynMrAGYRehJjwDOcvcfdVN/LPBeQk6KaYSbSofnwbzj44mEOODJhIx0q9x9fapulxkHzexthPH6VsKY/YJqa4+vK/1kp1/kG/W4OyH2qFsIkxuGxtJldj53X21mvyGsVjMNGAvcamaHZWXeqV8F7u7LCZEJTyTH02bRjWlMI0wEgXBD7dF+ppXtNdKfrX6RbxqyFpBjlhEiAoYAJ5hZU1fmGy+mZ4FvUozBHQNcb2Y7u3t77MXXlO4u8PQxM7P0Nj4eRYhZfyfhC+yGaEI1Mw3pz1a/yC+D2rjTF0oXXE8w4lnAcd2Yb9LOyvh4OSEp0BTg92a2fdIDq4DsimBmo8xs9/S+xBRisv8PAV+Nhy5w99/HY7l4D9IvBjWegzukWRSK4XrNBGM+E/g34ECK0SNHUowC+BMhDLC75afOjHX/AnyfMOutQEjFOTbr91yiNUk3+nnCElnTgZ3j/htS7zudsKhPM/L6el5e9Nf7319l4JXMBWTypoumPZ0w6eHh1IWylDBLLanzodSxywjRAUmMbToX8liKSe/fREgQdGV8/jiwVz/0VvSCJSzm8HLqfa0j5F55iWL42irgE6V/swq9fn9nJ9ZEf1d/93rRrzJwy6C7OZncKDSzbQk9410JM9zWE/JhrwWeJ4xttwOXAtsRYnJPJoxdX2pmF3lxDb9tCau5/zvhIlzp7uvM7DzCxIo5wDGEqczlaN3P3e/x1BJo/XjfOxFSyv6TEC3zbeAw4LWEyJlmwAnG8TPgene/LJ7bp5urZjaH8N4PIfxdVhP+nhuA1nJutGWhP547jTDVfDjwmEfKbTMr/WJgMqjCAROjMLPpwN8IKVrvA34E3EaICNkWuMzdV6fOm0Mw7rPjrpWEsK47CREnuwIHx2P/6+6fjmORw+PrHEDoPV1QhtYvEvKjfMDdfxD39cm8zWwH4KcEozjJ3f8Yx+snAccBEwlTrRcQ8jw/4MW1DPv6mu8jLCSwZ8mh+4DfA7909xd6Y94Z6T+V8GX7RsK9oDHAH4Bb3f27sU6vDDUL/WKAk3WXv9aF0LP5NcX8DwfQcarxsLhtIvTAd4jbSYQx7BZCrpIkeiT5udtKKhdybGMy8FQ8/rEydX4/nrceODu1v6xhEzqml20D9qzB3/hTqb/LBkIPc1l8LwXgaeA7wJSc6v9YyWe7NvV4PV2kas2LfpWBXzIXUPM3HHrUC4A1wKmE3lRpprVdCT9l7yP0rhcA74umvzthpe0bCT3uhcDXgFNT5xthqOXthOGBhZQsNdULneelzGJTX8ybLXOCV306PsXFagvAl4FXx/07EW7CrYrHFgFndvd+MtJ/bkr/xYT1P48nrDKzJhr6y8Bpefz7qwyOkrmAmr/hcOMwuQk5u+TY/oSloNI9rCR963OEnmRpfu7RnbxGI2GR19viuZfQzbqCJecmeVDeF00iyXOxmTCDs0O9btqp+UIOwH4Ul9f6fPo9p97Xgam/781d6clI/yGEKekFwpT0qSXH35X6f/hq3v7+KoOnZC6g5m+4aNzrgZ3jviHRsO+lOAzyPOEn/sWEn/YFQorNbeM5ljKj5mj6wwgrbZ+cMu3Hge37oPPfU18ey8sxb8KviptqbRqEXxibCMND+6T2W1ov8K3UF9JJOdJ/RtS1kLg0Xdw/NPX4/6KufwHDu2hnuyz0qwyeMuiiSoDbgYcIQx5/N7M7COPXB6fq3AT8EviTu68ys6MJBr494WbShe7ugFvIU/I/wLEUF1IYTTDw54Dj3X1Rb8WlbkZNjLseJCyDdWrU+V0zw91/6F1Hm7QQhntaCTlTHqx2ZEK82XYM4UvwcXe/NzkW/1akdP6L8LcCmNlJc1npPyrqWuzu85Nj7r7ZzIa4exvwQtx9KzDRzMYRPpeHgbUebipuJCQcq5l+MbgYjMa9CriKEHc9gxA1kHAfYfLDpz1mAowkY+ITCCFbaWYARxAm8SRm9CwhV/c57v5kmfqS9neN25WE9QiXE6JaujXv+Hi5me0LjHD3hbUwjagjef8TzWy4u29M10lFkNxB+EUznWBw6TDNPOifbmYjPZUIyt3bzGwEsFfcdSzwOsICvyMJ60XebGbfdvdHzWx/Qpz/EzJtUWkGnXG7e4uZfYswxn0yIWpkAyGZ/d/d/R7YItRrFsG0k6GPdHsPmdm7CClg5xLilf8CPOHuL/dBX2Lc4+J2kbtvAj4be4Vn0Y15px4/F99HQ7VNI2XIS+KuYYThgsfS4X6p97aGmLiLEHlDojFj/UkSqCEEg76t5PgMwirspLYQhld2JrznA83sne7+YK30i8HHoDNuAA/Z2r4DfCf+1MXdVyXHLSSUSibXjCYMAYwh/DzeogcdL9IHK6XPzOYSJmZAGNpJXuf8mKqiR/NOnVP1GOCUId9ImLK9O2GBiQ+kjqVJln2DooGn28tK/y1xOwd4PdG43d3jl+b3CL8QlgOXEz6bNkL8/6cIqVv3Ar5pZu9y92droV8MPgalcUOHxPWrkuRRZjYZaHP3FfH5TOAUghkB/M7dl3TXXunjPvIQYWXvYwnDN5jZUHff3Fvz7sdr94d/EFIIvBXY08xmuvviTuq1EMIkpxB656TGkLPkNsJCz9MI4aBp5hCG114mDFld6+5rk4Nm9hfgCkI0yW7Aq4BfaQKNqAaDauZkd8TxzW8RQsJ+RvhZvDvBPCEk+/lErFu1nMjJhR71fNjdv57WmPzsNrPPUjTvdsJ4+g/TbVRDX0+Y2VGEmO1b3f3hLuqMJ0TsbA38xd3ndabZzIYAze6+plbvycymElIWtJaOTZvZsYShnd+5e0vcl0QXtZvZnoRfZSOBP7v78dXWKwYnMu6ImX0Y+EZ82kYxV4kBX3b3T8d6VbnRVNJjLzWM0puPyePcmHfJl0pz2thKv+TMbAYhCmMUcIW7v6mT9oYRwvNOAs71Kq/80t3fq+SzGZ3uaafqJDdXf0kIi3wIOMLd11RTtxicDNqhkk64hpBXYzdCj3EoISHS39z9IqiOaaeGCBqBts5eo3TMN6mTh2GTlP5kMYDGxLSj3s56Bg2EG8KjCKFzpcY/DHg3YUwZ4HBCZE819TcAhS7+/p56vIVpR5K/8/MUQ0IHdb57UT30jxVx90eB9xBC+/YGdnL3M1KmXZXoAA9hZjsA881sbuy1NfZwzit13P18QkTMMoL5f9fM3huPVX3xhr7oB8YTEnBBiO6BGAbZiWl/x7tZ77O/9FF/Z+0k5j4qbpMZlkJUHBl3Cndvc/eN7v6Uuz+T7I8/lSt+EVpgBCGh1M6Epc52q4B5f8/MPhKPVW0srB/6h1A07lcm53Ri2t9y9w/H1yrbTKuoPzm/IbmxHZ9vR7gv0g78VsMkolrIuHtBpc0vudhjVMsGwgo59xGGO/7RD/P+ASFUrZEwbbwqVEB/I8WJRptim8PpaNrfdvf/jMcqOkTVH/1mNsJCuCaeCr2MN1zfRFg1aT1l5l4Xohw0xl1jUjexphDykexBGJqZSbgpOo6wQvzh3osV4r0427Dg7p81szHAJnf/ZI71J8mzhlKM4/4PQlQPBNM+N/16OdJ/IHClmV0PXEeI/T6YMAHr/Sn9V1VKsxClKKqkhqRMY1vgasLPagg5LZJUsMmvoFWEPBcP98a8OqtTRdPrl34LC1kkPdybCZNZvhMP18K0+6TfzJoJmR7nxTpJru6x8XED8DV3/3h8PcVwi6og464xFpbCuokQuXIvcBEhP8pQQmKpjwG7EOLIV1GeeVdyElDV9JvZhHjuDMJEnLGx+aqZdgX0J/H17yCkr01WQ0qGU/4MXOPu36+mfiGAwZfWNatC6NENBb5O6J3NJ6SCTVaUT1KeTiVMF38g1lsB7BaPZZYWtFL6KXYWLiCMbycrzXwn9VoVf5+V0F/S3i6E8ezXE6a5j6ymfhWVdNHNyRrh7k4YQ9037roLeMhjFkKPoXvu/iLhZtnXCJNUxhFumO3ufQxVqwSV0k+xhzqfMDRhwPfc/YNQvZ5qJfRHfUNivQXufrO7X+Pu93vMJBiPqactqoqMu0bESIbZhF4ewF0eF4RNiOaCh0keVwG/IuR/HkfG5l1B/W2xzq8Jk55+VW3TroL+Lq+bpA0hqomMu0Z4uEnVEksbMZqiKxP2sMr8bwgz8SCMA9+alXlXUn+q2gnu/s6knWr2VOv97y9EGhl3bVlGWBVnCHCChfSxnZpV/Mn9LPBNiivKjyFMEtk5mketP7+K6Y/VkpuVtcpZXe9/fyEAGXfFMetxivn1BCOYBRzXzcWftLMyPl5OWBl9CvB7M9veqzClPQv9XsGQuXr/+wvRG2TcFST+3HczazazWWZ2ppn9m5kdaDGfNiF+uIFgHGcAe3R28afMbDrBOO4jjAmvJCTCOt/MxlZyTFX6s9UvRK/JOqxloBRiCBjhQv85ISIhWaV9KWFB4aTOh1LHLiPcMBsSjzWl2hwL/CTWexMhgdGVFJdQ20v6B4Z+FZVySuYCBkJJGcK2hDzMBUKM8tr4+EnClO6mWG8qIfd3Yh7XAWcSY4pTbSUGswZ4ddy/R3xeAP5L+utfv4pKuUW5SvpJErcbp3FfT1h8+D7gR4SlsKYRTOAyj+tYuvuLZvYDwrJdZxPWl9wXONvM7iREPOxKyIEBYXLKTfEn/ULgEeAAtlxxXvrrTL8QfSLrb46BUIBm4NeEXtjdhIt6SOr4sLhtIixrtUPcTiL09FoICwok0QtJT7AV+GLJa00GnorHPyb99a9fRaXcoh53ZZhI6LGtA74K3EOqN+bum8xsV+C9hIUaZgIvEhIr/RK4E3gdIXnRKGACYSz1Lnf/HbwSLdEY600gmMefpX9A6BeiPLL+5hgIhXDjKrkJNrvk2P7ARymOtybjrwVCTPGniHkyUueM7uQ1GoF9CD//C4QsdWOkv/71q6iUWzIXMBBKyjjWAzvHfUOiYdxL8Wf484TkRRcDT1O8cbZtPCdZMRzCz//9CeOwE4CTU6bxOLC99A8M/Soq5ZbMBQyEQghBmx8v6hcJP7OTizwpNwDvAMbFc44mzOQrAO8raW8sYcbeY4SbYU+l6j5DWA9T+geIfhWVcovGuCvDKkJSorGEPM5vTB27j2Aan/aYiS6ygBBWNoEtoxNmEMZiZ1HMpvcsYQLIOe7+pPQPKP1ClIWMuwK4e4uZfYswxnoyIWphA2EB37+7+z2wRSKlWQTTSH56p9t7yMzeRbjhNpdw0+0vwBPu/rL0Dyz9QpSLVsCpAmY2DsDdV6X2NXmMIzaz0cAngU8AtwLvcPcltVfaOdIvRL5Rj7uCxMkg7u6rkuRFZjYZaHP3FfH5TOAUgmkA/K4r00jaK30s/QNTvxC9RT3uKmIhZ/O3gEOAnxHGTncHjo1VLnD3T8S6uTMG6Rcin8i4q4iZfZiQEwNC8v4hhJl5BnzZ3T8d6+VyYVnpFyKfaKikulwD7ElIA7oTIQfGpcDf3P0iyL1pSL8QOUQ97ipjZkMIX5DTgc3u/kzqWINXcBGBaiD9QuQPGXdG1PuYqvQLkR0ybiGEqDO0dJkQQtQZMm4hhKgzZNxCCFFnyLiFEKLOkHELIUSdIeMWQog6Q8YthBB1hoxbCCHqDBm3EELUGTJuIYSoM2TcQghRZ8i4hRCizpBxCyFEnSHjFkKIOkPGLYQQdYaMWwgh6gwZtxBC1BkybiGEqDNk3EIIUWfIuIUQos6QcQshRJ0h4xZCiDpDxi2EEHWGjFsIIeoMGbcQQtQZMm4hhKgzZNxCCFFnyLiFEKLOkHELIUSdIeMWQog6Q8YthBB1hoxbCCHqDBm3EELUGf8PclK6QKSYZwsAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 240x420 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Swarm plot of 4 categories\n",
    "\n",
    "outputPrefix = '190523-HSP47-epi-ROI-swarm'\n",
    "outputFigPath = \"../jupyter_figures/\" + outputPrefix + \".svg\"\n",
    "\n",
    "fig_width=0.8\n",
    "fig_height=1.4\n",
    "fig = plt.figure(figsize=(fig_width,fig_height), dpi=300)\n",
    "ax = fig.add_axes([0.1, 0.1, 0.8, 0.8])\n",
    "\n",
    "ax = sns.swarmplot(x='groups', y='mean_intensity', data=df,\n",
    "                   hue='gland_id',\n",
    "#                    palette={'1': 'red', '2': 'green', '3': 'blue', '4': yellow', '5': 'purple'},\n",
    "#                    color=\"blue\",\n",
    "                   size=2.0, alpha=.8)\n",
    "\n",
    "ax = sns.barplot(x='groups', y='mean_intensity', data=df,\n",
    "                 color=\".7\", alpha=1.0,\n",
    "                 errwidth=.7, errcolor=\"k\", capsize=.2, ci='sd')\n",
    "\n",
    "plt.ylim(0, 1.8)\n",
    "# plt.yticks([0, 4, 8])\n",
    "\n",
    "plt.xlabel(None)\n",
    "plt.ylabel(\"Mean HSP47 intensity\")\n",
    "\n",
    "# make the bar width narrower\n",
    "change_bar_width(ax, .5)\n",
    "\n",
    "# rotate x tick labels\n",
    "x_labels = ax.get_xticklabels()\n",
    "ax.set_xticklabels(labels=x_labels, rotation=45, ha=\"right\")\n",
    "\n",
    "# Remove legend (symbols representing orgnoids)\n",
    "ax.get_legend().remove()\n",
    "\n",
    "for o in fig.findobj():\n",
    "    o.set_clip_on(False)\n",
    "for o in ax.findobj():\n",
    "    o.set_clip_on(False)\n",
    "\n",
    "if outputFigPath is not None:\n",
    "    plt.savefig(outputFigPath)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>count</th>\n",
       "      <th>mean</th>\n",
       "      <th>std</th>\n",
       "      <th>min</th>\n",
       "      <th>25%</th>\n",
       "      <th>50%</th>\n",
       "      <th>75%</th>\n",
       "      <th>max</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>groups</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>sg-Control</th>\n",
       "      <td>7.0</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>0.133765</td>\n",
       "      <td>0.841262</td>\n",
       "      <td>0.930510</td>\n",
       "      <td>0.953481</td>\n",
       "      <td>1.043909</td>\n",
       "      <td>1.256417</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>sg1-Serpinh1</th>\n",
       "      <td>7.0</td>\n",
       "      <td>0.195582</td>\n",
       "      <td>0.109826</td>\n",
       "      <td>0.077105</td>\n",
       "      <td>0.099136</td>\n",
       "      <td>0.168279</td>\n",
       "      <td>0.302483</td>\n",
       "      <td>0.320450</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>sg2-Serpinh1</th>\n",
       "      <td>4.0</td>\n",
       "      <td>0.363992</td>\n",
       "      <td>0.120840</td>\n",
       "      <td>0.232612</td>\n",
       "      <td>0.310126</td>\n",
       "      <td>0.349433</td>\n",
       "      <td>0.403299</td>\n",
       "      <td>0.524491</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>sg3-Serpinh1</th>\n",
       "      <td>6.0</td>\n",
       "      <td>0.264319</td>\n",
       "      <td>0.093387</td>\n",
       "      <td>0.185096</td>\n",
       "      <td>0.209514</td>\n",
       "      <td>0.219051</td>\n",
       "      <td>0.304585</td>\n",
       "      <td>0.423285</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "              count      mean       std       min       25%       50%  \\\n",
       "groups                                                                  \n",
       "sg-Control      7.0  1.000000  0.133765  0.841262  0.930510  0.953481   \n",
       "sg1-Serpinh1    7.0  0.195582  0.109826  0.077105  0.099136  0.168279   \n",
       "sg2-Serpinh1    4.0  0.363992  0.120840  0.232612  0.310126  0.349433   \n",
       "sg3-Serpinh1    6.0  0.264319  0.093387  0.185096  0.209514  0.219051   \n",
       "\n",
       "                   75%       max  \n",
       "groups                            \n",
       "sg-Control    1.043909  1.256417  \n",
       "sg1-Serpinh1  0.302483  0.320450  \n",
       "sg2-Serpinh1  0.403299  0.524491  \n",
       "sg3-Serpinh1  0.304585  0.423285  "
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# get descriptive stats for all groups\n",
    "df.groupby('groups').describe()['mean_intensity']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}