CSV_to_stats.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "0ec00416",
   "metadata": {},
   "outputs": [],
   "source": [
    "import time\n",
    "import os\n",
    "import math\n",
    "import numpy as np\n",
    "import cv2\n",
    "import pandas as pd\n",
    "import pyautogui as pa\n",
    "from scipy.stats import binom_test, kstest\n",
    "from matplotlib.markers import TICKDOWN\n",
    "import datetime\n",
    "\n",
    "from lib.auto_GUI.auto_GUI_base import AutoGUIBase\n",
    "from lib.auto_GUI.auto_PhotoZ import AutoPhotoZ\n",
    "from lib.analysis.laminar_dist import *\n",
    "####################################################################################\n",
    "#    Aggregates SNR (.dat), latency (.dat), and ROI photos (png)\n",
    "#            scraped from PhotoZ\n",
    "#    Given ROI files in format slice#_rois.dat \n",
    "# 8/23/23 We want to add:\n",
    "#     Take into account the RELATIVE STIM LOCATION within every barrel. Does velocity depend on distance from stimulating location?\n",
    "#            New column: distance (um) from stim loc to relevant barrel pair boundary (1:2)\n",
    "#            New column: distance (um) from stim loc to relevant barrel pair boundary (2:3)\n",
    "#\n",
    "#     METHOD\n",
    "#         To calculate stim location to boundary, take these files produced by MovieMaker2.ipynb:\n",
    "#             electrode.dat\n",
    "#             corners_layer.dat\n",
    "#             corners_barrel.dat\n",
    "#         Calculate: Distance to lower barrel, distance to higher barrel,\n",
    "#                    Distance to lower layer, distance to higher layer\n",
    "#                       (if applicable for each)\n",
    "####################################################################################"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "5b23cb4f",
   "metadata": {
    "scrolled": true
   },
   "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>Date</th>\n",
       "      <th>Slice/Loc/Rec</th>\n",
       "      <th>Slice Edge</th>\n",
       "      <th>Most Medial Barrel</th>\n",
       "      <th>Medial Dir</th>\n",
       "      <th>Slice Location</th>\n",
       "      <th>Stim Location</th>\n",
       "      <th>Actual Stim Time (ms)</th>\n",
       "      <th>Stim Barrel</th>\n",
       "      <th>L4 barrel 1 Latency</th>\n",
       "      <th>...</th>\n",
       "      <th>abs(delta 1)</th>\n",
       "      <th>abs(delta 2)</th>\n",
       "      <th>Distance to Barrel 1</th>\n",
       "      <th>Distance to Barrel 2</th>\n",
       "      <th>Distance to Layer n1</th>\n",
       "      <th>Distance to Layer n2</th>\n",
       "      <th>Distance to Barrel A</th>\n",
       "      <th>Distance to Barrel B</th>\n",
       "      <th>Distance to Layer nA</th>\n",
       "      <th>Distance to Layer nB</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2/8/2023</td>\n",
       "      <td>02_01_02</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>L4</td>\n",
       "      <td>48.4226</td>\n",
       "      <td>NaN</td>\n",
       "      <td>47.7369</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>215.228351</td>\n",
       "      <td>174.000000</td>\n",
       "      <td>180.093736</td>\n",
       "      <td>191.535134</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2/10/2023</td>\n",
       "      <td>04_01_07</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>L4</td>\n",
       "      <td>48.4429</td>\n",
       "      <td>NaN</td>\n",
       "      <td>52.7236</td>\n",
       "      <td>...</td>\n",
       "      <td>2.7731</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>151.086327</td>\n",
       "      <td>213.337695</td>\n",
       "      <td>149.891961</td>\n",
       "      <td>202.354632</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2/14/2023</td>\n",
       "      <td>02_01_01</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>L4</td>\n",
       "      <td>48.4645</td>\n",
       "      <td>NaN</td>\n",
       "      <td>47.5447</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>198.974687</td>\n",
       "      <td>232.730513</td>\n",
       "      <td>132.022414</td>\n",
       "      <td>229.523778</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2/14/2023</td>\n",
       "      <td>03_01_03</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>L4</td>\n",
       "      <td>48.4195</td>\n",
       "      <td>NaN</td>\n",
       "      <td>51.1559</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>0.9262</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>228.000000</td>\n",
       "      <td>162.000000</td>\n",
       "      <td>159.011987</td>\n",
       "      <td>199.155462</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2/28/2023</td>\n",
       "      <td>02_01_01</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>L5</td>\n",
       "      <td>48.4151</td>\n",
       "      <td>NaN</td>\n",
       "      <td>50.3009</td>\n",
       "      <td>...</td>\n",
       "      <td>NaN</td>\n",
       "      <td>1.4949</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>44.346286</td>\n",
       "      <td>328.107622</td>\n",
       "      <td>171.794234</td>\n",
       "      <td>154.795349</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>169</th>\n",
       "      <td>10/19/2023</td>\n",
       "      <td>05_06_02</td>\n",
       "      <td>down/right</td>\n",
       "      <td>1.0</td>\n",
       "      <td>down/left</td>\n",
       "      <td>Anterior</td>\n",
       "      <td>L2/3</td>\n",
       "      <td>48.4228</td>\n",
       "      <td>1.0</td>\n",
       "      <td>52.5402</td>\n",
       "      <td>...</td>\n",
       "      <td>1.2641</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>24.997803</td>\n",
       "      <td>181.431965</td>\n",
       "      <td>314.744870</td>\n",
       "      <td>103.797327</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>170</th>\n",
       "      <td>10/19/2023</td>\n",
       "      <td>06_01_01</td>\n",
       "      <td>down</td>\n",
       "      <td>3.0</td>\n",
       "      <td>right</td>\n",
       "      <td>Posterior</td>\n",
       "      <td>L4</td>\n",
       "      <td>48.5290</td>\n",
       "      <td>2.0</td>\n",
       "      <td>50.6833</td>\n",
       "      <td>...</td>\n",
       "      <td>0.3503</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>118.898651</td>\n",
       "      <td>202.376373</td>\n",
       "      <td>288.667605</td>\n",
       "      <td>69.296479</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>171</th>\n",
       "      <td>10/19/2023</td>\n",
       "      <td>06_02_02</td>\n",
       "      <td>down</td>\n",
       "      <td>3.0</td>\n",
       "      <td>right</td>\n",
       "      <td>Posterior</td>\n",
       "      <td>L4</td>\n",
       "      <td>48.4103</td>\n",
       "      <td>3.0</td>\n",
       "      <td>50.3022</td>\n",
       "      <td>...</td>\n",
       "      <td>0.7656</td>\n",
       "      <td>1.0368</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>151.005549</td>\n",
       "      <td>215.305522</td>\n",
       "      <td>75.027164</td>\n",
       "      <td>306.897033</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>172</th>\n",
       "      <td>10/19/2023</td>\n",
       "      <td>06_03_01</td>\n",
       "      <td>down</td>\n",
       "      <td>3.0</td>\n",
       "      <td>right</td>\n",
       "      <td>Posterior</td>\n",
       "      <td>L2/3</td>\n",
       "      <td>48.5911</td>\n",
       "      <td>3.0</td>\n",
       "      <td>52.3050</td>\n",
       "      <td>...</td>\n",
       "      <td>0.5673</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>322.923386</td>\n",
       "      <td>22.555042</td>\n",
       "      <td>70.935246</td>\n",
       "      <td>350.093066</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>173</th>\n",
       "      <td>10/19/2023</td>\n",
       "      <td>06_04_01</td>\n",
       "      <td>down</td>\n",
       "      <td>3.0</td>\n",
       "      <td>right</td>\n",
       "      <td>Posterior</td>\n",
       "      <td>L2/3</td>\n",
       "      <td>48.4633</td>\n",
       "      <td>2.0</td>\n",
       "      <td>53.4492</td>\n",
       "      <td>...</td>\n",
       "      <td>2.4054</td>\n",
       "      <td>0.5346</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>222.331315</td>\n",
       "      <td>55.298944</td>\n",
       "      <td>239.536169</td>\n",
       "      <td>171.112560</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>174 rows × 37 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "           Date Slice/Loc/Rec  Slice Edge  Most Medial Barrel Medial Dir  \\\n",
       "0      2/8/2023      02_01_02         NaN                 NaN        NaN   \n",
       "1     2/10/2023      04_01_07         NaN                 NaN        NaN   \n",
       "2     2/14/2023      02_01_01         NaN                 NaN        NaN   \n",
       "3     2/14/2023      03_01_03         NaN                 NaN        NaN   \n",
       "4     2/28/2023      02_01_01         NaN                 NaN        NaN   \n",
       "..          ...           ...         ...                 ...        ...   \n",
       "169  10/19/2023      05_06_02  down/right                 1.0  down/left   \n",
       "170  10/19/2023      06_01_01        down                 3.0      right   \n",
       "171  10/19/2023      06_02_02        down                 3.0      right   \n",
       "172  10/19/2023      06_03_01        down                 3.0      right   \n",
       "173  10/19/2023      06_04_01        down                 3.0      right   \n",
       "\n",
       "    Slice Location Stim Location  Actual Stim Time (ms)  Stim Barrel  \\\n",
       "0              NaN            L4                48.4226          NaN   \n",
       "1              NaN            L4                48.4429          NaN   \n",
       "2              NaN            L4                48.4645          NaN   \n",
       "3              NaN            L4                48.4195          NaN   \n",
       "4              NaN            L5                48.4151          NaN   \n",
       "..             ...           ...                    ...          ...   \n",
       "169       Anterior          L2/3                48.4228          1.0   \n",
       "170      Posterior            L4                48.5290          2.0   \n",
       "171      Posterior            L4                48.4103          3.0   \n",
       "172      Posterior          L2/3                48.5911          3.0   \n",
       "173      Posterior          L2/3                48.4633          2.0   \n",
       "\n",
       "     L4 barrel 1 Latency  ...  abs(delta 1)  abs(delta 2)  \\\n",
       "0                47.7369  ...           NaN           NaN   \n",
       "1                52.7236  ...        2.7731           NaN   \n",
       "2                47.5447  ...           NaN           NaN   \n",
       "3                51.1559  ...           NaN        0.9262   \n",
       "4                50.3009  ...           NaN        1.4949   \n",
       "..                   ...  ...           ...           ...   \n",
       "169              52.5402  ...        1.2641           NaN   \n",
       "170              50.6833  ...        0.3503           NaN   \n",
       "171              50.3022  ...        0.7656        1.0368   \n",
       "172              52.3050  ...        0.5673           NaN   \n",
       "173              53.4492  ...        2.4054        0.5346   \n",
       "\n",
       "     Distance to Barrel 1  Distance to Barrel 2  Distance to Layer n1  \\\n",
       "0                     NaN                   NaN                   NaN   \n",
       "1                     NaN                   NaN                   NaN   \n",
       "2                     NaN                   NaN                   NaN   \n",
       "3                     NaN                   NaN                   NaN   \n",
       "4                     NaN                   NaN                   NaN   \n",
       "..                    ...                   ...                   ...   \n",
       "169                   NaN                   NaN                   NaN   \n",
       "170                   NaN                   NaN                   NaN   \n",
       "171                   NaN                   NaN                   NaN   \n",
       "172                   NaN                   NaN                   NaN   \n",
       "173                   NaN                   NaN                   NaN   \n",
       "\n",
       "    Distance to Layer n2  Distance to Barrel A  Distance to Barrel B  \\\n",
       "0                    NaN            215.228351            174.000000   \n",
       "1                    NaN            151.086327            213.337695   \n",
       "2                    NaN            198.974687            232.730513   \n",
       "3                    NaN            228.000000            162.000000   \n",
       "4                    NaN             44.346286            328.107622   \n",
       "..                   ...                   ...                   ...   \n",
       "169                  NaN             24.997803            181.431965   \n",
       "170                  NaN            118.898651            202.376373   \n",
       "171                  NaN            151.005549            215.305522   \n",
       "172                  NaN            322.923386             22.555042   \n",
       "173                  NaN            222.331315             55.298944   \n",
       "\n",
       "     Distance to Layer nA  Distance to Layer nB  \n",
       "0              180.093736            191.535134  \n",
       "1              149.891961            202.354632  \n",
       "2              132.022414            229.523778  \n",
       "3              159.011987            199.155462  \n",
       "4              171.794234            154.795349  \n",
       "..                    ...                   ...  \n",
       "169            314.744870            103.797327  \n",
       "170            288.667605             69.296479  \n",
       "171             75.027164            306.897033  \n",
       "172             70.935246            350.093066  \n",
       "173            239.536169            171.112560  \n",
       "\n",
       "[174 rows x 37 columns]"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Load data, not from Drive for Desktop since PhotoZ complains about loading zda from Drive\n",
    "data_dir = \"C:/Users/jjudge3/Desktop/Data/mm_full_pipeline_targets/\"\n",
    "# data_dir = 'C:/Users/jjudg/Downloads/'\n",
    "cell_type_prefix = 'scnn1a'\n",
    "csv_prefix = cell_type_prefix + '_all_slices-final.csv'\n",
    "\n",
    "# load the all-data CSV dir\n",
    "df = pd.read_csv(data_dir + csv_prefix)\n",
    "df = pd.DataFrame(df)\n",
    "df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "13bf9866",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Barrel pairs in L2/3\n",
      "\t 19 no crossing out of 120 in layer L2/3\n",
      "\t 48 crossing out of 120 in layer L2/3\n",
      "\t 53 no signal out of 120 in layer L2/3\n",
      "Barrel pairs in L4\n",
      "\t 69 no crossing out of 196 in layer L4\n",
      "\t 35 crossing out of 196 in layer L4\n",
      "\t 92 no signal out of 196 in layer L4\n",
      "Barrel pairs in L5\n",
      "\t 10 no crossing out of 28 in layer L5\n",
      "\t 7 crossing out of 28 in layer L5\n",
      "\t 11 no signal out of 28 in layer L5\n",
      "\n",
      " Anterior\n",
      "Barrel pairs in L2/3\n",
      "\t 13 no crossing out of 58 in layer L2/3\n",
      "\t 18 crossing out of 58 in layer L2/3\n",
      "\t 27 no signal out of 58 in layer L2/3\n",
      "\n",
      " Posterior\n",
      "Barrel pairs in L2/3\n",
      "\t 4 no crossing out of 42 in layer L2/3\n",
      "\t 15 crossing out of 42 in layer L2/3\n",
      "\t 23 no signal out of 42 in layer L2/3\n",
      "\n",
      " Anterior\n",
      "Barrel pairs in L4\n",
      "\t 26 no crossing out of 84 in layer L4\n",
      "\t 12 crossing out of 84 in layer L4\n",
      "\t 46 no signal out of 84 in layer L4\n",
      "\n",
      " Posterior\n",
      "Barrel pairs in L4\n",
      "\t 29 no crossing out of 82 in layer L4\n",
      "\t 13 crossing out of 82 in layer L4\n",
      "\t 40 no signal out of 82 in layer L4\n",
      "\n",
      " Anterior\n",
      "Barrel pairs in L5\n",
      "\t 5 no crossing out of 10 in layer L5\n",
      "\t 1 crossing out of 10 in layer L5\n",
      "\t 4 no signal out of 10 in layer L5\n",
      "\n",
      " Posterior\n",
      "Barrel pairs in L5\n",
      "\t 1 no crossing out of 6 in layer L5\n",
      "\t 1 crossing out of 6 in layer L5\n",
      "\t 4 no signal out of 6 in layer L5\n"
     ]
    }
   ],
   "source": [
    "# filter out rows marked in column Exclude \n",
    "df2 = df[df[\"Exclude Reason\"].isnull()]\n",
    "data = {}\n",
    "\n",
    "for stim_loc in ['L2/3', 'L4', 'L5']:\n",
    "    df2 = df[df['Stim Location'].str.contains(stim_loc)]\n",
    "    print(\"Barrel pairs in\", stim_loc)\n",
    "    # count number of barrel pairs with no signal crossing\n",
    "    n_barrel_no_cross = len(df2[df2['# no crossings 2/3'] > 0])\n",
    "    n_barrel_no_cross += len(df2[df2['# no crossings 1/2'] > 0])\n",
    "    print(\"\\t\", n_barrel_no_cross, \"no crossing out of\", len(df2) * 2, \"in layer\", stim_loc)\n",
    "\n",
    "    # count number of barrel pairs with signal crossing\n",
    "    n_barrel_cross = len(df2[df2['# crossings 2/3'] > 0])\n",
    "    n_barrel_cross += len(df2[df2['# crossings 1/2'] > 0])\n",
    "    print(\"\\t\", n_barrel_cross, \"crossing out of\", len(df2) * 2, \"in layer\", stim_loc)\n",
    "\n",
    "    # count number of barrel pairs with no signal\n",
    "    n_barrel_none = len(df2) * 2 - n_barrel_cross - n_barrel_no_cross\n",
    "    print(\"\\t\", n_barrel_none, \"no signal out of\", len(df2) * 2, \"in layer\", stim_loc)\n",
    "    \n",
    "    data[stim_loc] = {}\n",
    "    data[stim_loc]['n_barrel_no_cross'] = n_barrel_no_cross\n",
    "    data[stim_loc]['n_barrel_cross'] = n_barrel_cross\n",
    "    data[stim_loc]['n_barrel_none'] = n_barrel_none\n",
    "    data[stim_loc]['total'] = len(df2) * 2 - n_barrel_none\n",
    "    \n",
    "# Now rostral vs caudal\n",
    "for stim_loc in ['L2/3', 'L4', 'L5']:\n",
    "    df2 = df[df['Stim Location'].str.contains(stim_loc)]\n",
    "    for slice_loc in ['Anterior', 'Posterior']:\n",
    "        print(\"\\n\", slice_loc)\n",
    "        df2_sl = df2[df2['Slice Location'] == slice_loc]\n",
    "        print(\"Barrel pairs in\", stim_loc)\n",
    "        # count number of barrel pairs with no signal crossing\n",
    "        n_barrel_no_cross = len(df2_sl[df2_sl['# no crossings 2/3'] > 0])\n",
    "        n_barrel_no_cross += len(df2_sl[df2_sl['# no crossings 1/2'] > 0])\n",
    "        print(\"\\t\", n_barrel_no_cross, \"no crossing out of\", len(df2_sl) * 2, \"in layer\", stim_loc)\n",
    "\n",
    "        # count number of barrel pairs with signal crossing\n",
    "        n_barrel_cross = len(df2_sl[df2_sl['# crossings 2/3'] > 0])\n",
    "        n_barrel_cross += len(df2_sl[df2_sl['# crossings 1/2'] > 0])\n",
    "        print(\"\\t\", n_barrel_cross, \"crossing out of\", len(df2_sl) * 2, \"in layer\", stim_loc)\n",
    "\n",
    "        # count number of barrel pairs with no signal\n",
    "        n_barrel_none = len(df2_sl) * 2 - n_barrel_cross - n_barrel_no_cross\n",
    "        print(\"\\t\", n_barrel_none, \"no signal out of\", len(df2_sl) * 2, \"in layer\", stim_loc)\n",
    "\n",
    "        data[stim_loc][slice_loc] = {}\n",
    "        data[stim_loc][slice_loc]['n_barrel_no_cross'] = n_barrel_no_cross\n",
    "        data[stim_loc][slice_loc]['n_barrel_cross'] = n_barrel_cross\n",
    "        data[stim_loc][slice_loc]['n_barrel_none'] = n_barrel_none\n",
    "        data[stim_loc][slice_loc]['total'] = len(df2) * 2 - n_barrel_none"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "0d283eee",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'# no crossing': [6.505819318717424, 7.987138329728918, 18.601633295108112], '# crossing': [10.340601836232102, 5.6885382529804, 15.563243006262297]}\n",
      "{'# no crossing': [28.35820895522388, 66.34615384615384, 58.8235294117647], '# crossing': [71.64179104477611, 33.65384615384615, 41.17647058823529]}\n",
      "Binomial test for layer L2/3 : 0.00026080637840575794\n",
      "Binomial test for layer L4 : 0.00026714209499012487\n",
      "Binomial test for layer L5 : 0.0717315673828125\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEiCAYAAAD+lrlwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAtaElEQVR4nO3de7xWc97/8de7HRVCKXRACCNUUuE2o41xGm40k9MMI5kxY36O4xQGGXOI3AzGcFMqzIy6jduYMeJ2SJiGUhuRY5JOVJQpRPX5/bHW3q728Wq3r33t9no/H4/92Nc6f9b6Xutzfa/vWtd3KSIwM7PsaFHsAMzMrHE58ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME3/GSdpN0nRJ/5Z07jou+y1JbxYqtpztdJMUkloWeluNTdIdkq4sdhz1UahykfSopNMaYD3N9n2zvpp14pc0W9Lnkpbn/HUudlxNzCXAxIhoGxG3VJ4oaQ9Jj0v6RNJSSS9J+g5ARDwbEbs1esRVY8wt5w8ljZa0WR7LTZT0o0rjQlL3wkW7toj4aURcW59l0/i/SPd7maRJkvZq6Bjrq77lEhFHRsTYxogxX5KGSbqv2HE0lGad+FP/GRGb5fzNz53o2gA7AK/VMv1vwP8B2wBbA+cCnzZCXOvqPyNiM6AP0A/4RWMHUKT30tnpfm8FTATurc9Kqou9gfanQcslq+drQ+93FhJ/FWmt7v9Jeht4Ox13tKSytFb7T0k9c+bfW9K0tDlknKT7Jf0qnTZY0nPVrL97+rqVpBskzUlrPXdIapNOK5U0V9KFkj6StEDS6TnraSPpvyS9n9bonkvHPSLpnErbfEXScTXs7zGSXkv3baKk3dPxTwEHAb9Pa2W7VlquA7AjcFdEfJn+PR8Rz+XGnzN/n5xmo/9Jj9WvcuetZV+PSpf9VNIHkoblVZiVRMQ84FFgT0ntJP1d0qL0G8vfJXVNt/dr4Fs5+/57SZPS1bycjjsxnbe298ZsSZdKegVYIal7Wv6npWW+WNIVNcUraUy+x6iO/V4F3A/0yFl3f0mT07gXpPu4cc70tc6DnO1fKmkhMFpSC0lDJb0raYmk8ZLa5xNTpfjyKpc0ropvYun59bykmyR9DAxLj/Ez6TmxWNK4OjY/RNL89BhcmK53W0mfSdoqZ7v7pDFttC77lnN8/i3pdUkD0/GtJH2snG9hkrZW8i2oYzq8Lu+thkv+EdFs/4DZwLerGR8ktdj2QBuS2shHwL5ACXBaumwrYGPgfeACYCNgEPAV8Kt0XYOB56pZf/f09e+Ah9NttSWpQf82nVYKrAJ+ma77O8BnQLt0+m0ktbguaVz/kcZ0AvBCzvZ6AUuAjavZ112BFcCh6TYuAd4pnzdd/49qOH4i+WD8O3AcsE2l6aXA3PR1+XE6L93Od4Evc45TXftaCuxFUhnpCXwIHJdO65Ye05Z1lTOwHck3mGtJasHfAzZJj/3/AA/lLFdl33PLLh2u8b2Rs+2ydLttcmK9Kx3uBawEdq8h9jH5HqNqlq2IPz3+vwYm5UzfB9gPaJnGNRM4v5bzoHz715G8z9oA5wP/Arqm4/4b+HNjlQvJ+bUKOCfdjzbAn4ErSN4rrYFv1rD98vj+DGxK8v5alBPTP4Czcua/Cbi1hnUNA+6rYdrxQOc0nhNJzrdO6bQ/ANflzHse8Lf6vLcaNDc25Mqa2l964JYDS9O/h3Le8AfnzHc7cG2lZd8EBgAHAvMB5Uz7J3kkfpLEuQLYOWfa/sB76etS4PPcEyd9I+yXvok+B3pVs1+tgI+BXdLhG4A/1HAMrgTG5wy3AOYBpVHpJKth+a7A74F3gTXApJztlvJ14j8wXW/ucXqOtZNatftaw3Z/B9yUvu5G3QmmvJzfJznZqpwoQG/gk5zhKvtO1cRf43sjZ9tDcqaVx9o1Z9yLwEk1xD5mPY7RRJIPhqUkH7LLgENqKcvzgf+ttK+550Fpup7WOeNm5q4T6ERS8Sn/MClouZCcX3MqzX8PcGfuMa5h++XxfSNn3PXAqPT1icDz6esSYCHQv4Z1DaOGxF/NvGXAsenrfYEPgBbp8FTghPq8txryLwtNPcdFxJbp33E54z/Ieb0DcGH6dWuppKUkn7Kd0795kZZE6v08t92RpFbzUs56J6Tjyy2J5Gt6uc+AzYAOJLWZdyuvNCJWAuOBUyS1AE6m5rbdzrnxRsQakn3vks8ORMTciDg7InYmOU4rSE686rZT+Th9UGmemvYVSftKejr9qr0M+CnJMchXeTnvEBE/i4jPJW0i6b+VNJV9SvKhtaWkknVYb23vjZr2E5IkUmU/81DjMarBuRGxJcl75WjggfLmAkm7ps0oC9P9/w1Vj2nl2BdFxBc5wzsA/5uz7zOB1STXfPLREOVSOcZLSCpVLyppwhxSRwy5y7/P12X3V6CHpJ1IvhEvi4gX89yvCpJ+mNNcsxTYk/Q4R8QLJOfMAEnfIKkQPpwuWt/31nrLQuKvSeUE9eucD4gtI2KTiPgzsADoIkk582+f83oFSXIHkrbDnGmLSWpwe+Ssd4tILnbVZTHwBbBzDdPHAj8ADgE+i4jJNcw3n+QNVh6fSN5c8/KIYS0R8QFJ89Oe1Uyu7jhttw6r/xPJCbFdRGwB3EFycq+PC4HdgH0jYnOSbyXkrDeqXWpttb03yuWznoKKiDUR8SxJM95h6ejbgTdIvqFtDlxO1WNaOfbKwx8AR1ba/9aRtNnXV13lUtlaMUXEwoj4cUR0Bn4C/EG134mV+z7cnuScIP2AG09yHp1KPS6MS9qBpFnvbGCr9EN4RqV9GQuckm7jgZwP1qK9t7Kc+HPdBfw0rXVK0qZKLja2BSaTtDGeK6mlpO8C/XOWfRnYQ1JvSa1JvhICFbXru4CbJG0NIKmLpMPrCihd9m7gRkmdJZVI2l9Sq3T6ZJKml/+i9jfseOAoSYekF60uJGlz/mddMaQX4a5JL6a1UHKxdwhJm29lk0lqgmenx+lY1j5OdWkLfBwRX0jqD3x/HZatbZ2fA0uVXJC8utL0D4Gd6hhX23ujSZG0P8nF3fK7tNqS3IG1PK1tnlWP1d4B/DpNcEjqmJbt+qirXGol6Xh9fTH4E5LkuLqWRa5Mv2XsAZwO5F4MvoekOekYoK7bNVtIap3z14rk2kGQXDtAycX4yhWje4GBJMk/99ty0d5bTvxAREwFfkzSlv0JSa1pcDrtS5ILlYPTaScCD+Ys+xbJxbgnSC6ErnWHD3Bpur5/pV9rnyCp7eTjIuBVYApJm/51rF1m95BcsKrxDRsRb5K84W4l+RbxnyS32H2Zx/a/JGknfYIkgcwg+dAYXM12yo/TGSRtuqeQXBRemcd2AH4G/FLSv4GrSD6w1tfvSC4GLib5sJpQafrNwCAld5aU/4ZhGDA2/ep9Qm3vjSai/K6k5SQJ5hcR8Wg67SKSD9B/kySZuu5+qc7NJN/EHk/L5l8k7dbr43fUXi516Qe8kO7zw8B5EfFeLfM/Q1JuTwI3RMTj5RMi4nmSCtS0iJhdx3ZPJvnAKv97NyJeJ6l8TSapNOwFPJ+7UETMBaaRfEA8mzO+aO8trd0ka/mQNIbkomaj3yteKY4fAmdGxDeLGUdNJL0A3BERo4sdi1lNlNzW/KeIGFnAbdwNzC92ziiXyR9DNAeSNiGpJf+h2LGUkzSA5K6ExSTtpj1Z99qcWaOR1I/ktsr1bb6qbRvdSL4N712obawrN/VsgNJrBItIvlr+qcjh5NqN5JrHMpJrCYMiYkFxQzKrnqSxJM2Y50fEvwu0jWtJmkhH1NEc1ajc1GNmljGu8ZuZZYwTv5lZxmwQF3c7dOgQ3bp1K3YYZmYblJdeemlxRHSsPH6DSPzdunVj6tSpxQ7DzGyDIqna7mXc1GNmljFO/GZmGePEb2aWMRtEG7+ZNb6vvvqKuXPn8sUXX9Q9szW61q1b07VrVzbaaJ0eGJYoRCf/Df23zz77RLG1aNEievXqFXvssUcMGjQoVqxYkddys2bNiv79+0f37t3jhBNOiJUrV1Y735gxY6J79+7RvXv3GDNmTJ3Lz5w5M/bbb7/YeOONY8SIEeu/gxsgl0lhzZo1KxYtWhRr1qwpdihWyZo1a2LRokUxa9asWucDpkZGH8TSINq0aUNZWRkzZsxg44035o477qhx3i+//JIVK1YAcOmll3LBBRfw9ttv065dO0aNGlVl/o8//phrrrmGF154gRdffJFrrrmGTz75pNbl27dvzy233MJFF11UgL3dMLhMCuuLL75gq622Yu1HLFhTIImtttqq3t/GnPjr4Vvf+hbvvPNOlfEzZ87kwgsvZLfdduOtt94iInjqqacYNGgQAKeddhoPPfRQleUee+wxDj30UNq3b0+7du049NBDmTBhQq3Lb7311vTr169+X/OaIZdJYTjpN13rUzbNuo1/9erVfPjhh/VatnPnztWOX7VqFY8++ihHHHEEACtWrGD8+PGMGjWKiOD000/nlVdeoW3btixevJgtt9ySli2Tw9y1a1fmzav64KJ58+ax3XZfPySofL4lS5bktXxN5s+fn/e8ubbZZhtKStbl6YTrpr7lkuUygcKXS1N22WWXcfjhh7N06VLeeOMNhg4dWuyQ6mX+/Pmce+65PPDAA0WNo1kn/nfffZfddsv3mSdri0qd133++ef07t0bSGqXZ5xxBgCdOnWiZ8+ejBw5km984xu1rgOq/5Suab58l69Jly55PVa3ijfffJNdd921Xsvmo77lkuUygcKXS126DX2kQdc3e/hRec/7wgsvcNVVV3H55ZdXfNtqbKtWrar4wK+vzp07Fz3pQzNP/G3atAHgpZdeYtttt61j7rrXVVZWVmX8Aw88wKhRoxg4cCAnn3wyp512GjvskDzitkOHDixdurTiDTN37txqa61du3Zl4sSJFcNz586ltLQ07+Vrsi41UYCFCxeyzz77VBy3QmmocslCmUDjlUtTdPHFF/PYY4/x3nvvsf/++/Puu+/y5JNPMmjQIK666qq15h08eDCbb745U6dOZeHChVx//fUMGjSIiOCSSy7h0UcfRRK/+MUvOPHEE6ts65577uGGG25AEj179uTee+9l8ODBtG/fnunTp9OnTx9OPfVUfvrTn/LZZ5+x8847c/fdd9OuXTtuueUW7rjjDlq2bEmPHj24//77eeaZZzjvvPOApHIwadIklixZwtFHH82MGTMYM2YMDz/8MJ999hnvvvsuAwcO5Prrrwdg1KhRXHfddXTu3JlddtmFVq1a8fvf/77BjmuzTvzlX4u33XbbdTo518Vhhx3GYYcdxpIlS7jvvvs49thj6dChAyNHjqRbt24cdNBBPPDAA5x00kmMHTuWY4+t+ryHww8/nMsvv7zi4uHjjz/Ob3/7WyTltXxN6rvPhW5OKHS5NMcygcKXS1M0YsQIjj/+eO69915uvPFGSktLef7552ucf8GCBTz33HO88cYbHHPMMQwaNIgHH3yQsrIyXn75ZRYvXky/fv048MAD6dSpU8Vyr732Gr/+9a95/vnn6dChAx9//HHFtLfeeosnnniCkpISevbsya233sqAAQO46qqruOaaa/jd737H8OHDee+992jVqhVLly4F4IYbbuC2227jgAMOYPny5bRu3bpKvGVlZUyfPp1WrVqx2267cc4551BSUsK1117LtGnTaNu2LQcffDC9evVquIOKL+42mK222orzzjuPsrIyfvOb31ScpNdddx033ngj3bt3Z8mSJRXNEVOnTuVHP/oRkNwNcuWVV9KvXz/69evHVVddRfv27WtdfuHChXTt2pUbb7yRX/3qV3Tt2pVPP/20CHvedLlMmofp06fTu3dv3njjDXr06FHrvMcddxwtWrSgR48eFdeRnnvuOU4++WRKSkrYZpttGDBgAFOmTFlrufIL9h06dACoKGuA448/npKSEpYtW8bSpUsZMGAAkFzYnzRpEgA9e/bkBz/4Affdd19Fc9ABBxzAz3/+c2655RaWLl1abTPRIYccwhZbbEHr1q3p0aMH77//Pi+++CIDBgygffv2bLTRRhx//PH1PHI1a9Y1/oa0fPnyvOft379/xeuddtqJF198sco8ffv2ZeTIrx/xOWTIEIYMGVJlvpqW33bbbZk7d27eMTVHLpPmraysjMGDBzN37lw6dOjAZ599RkTQu3dvJk+eXG3TV6tWrSpel1+Pqe66TGURUeO1mk033bTO5R955BEmTZrEww8/zLXXXstrr73G0KFDOeqoo/jHP/7BfvvtxxNPPFGl1p8bb0lJCatWrcor3vXlGr+ZNUm9e/emrKyMXXfdlddff52DDz6Yxx57jLKysnW63nHggQcybtw4Vq9ezaJFi5g0adJaFQFIat7jx49nyZIlAGs19ZTbYostaNeuHc8++ywA9957LwMGDGDNmjV88MEHHHTQQVx//fUsXbqU5cuX8+6777LXXntx6aWX0rdvX95444284u3fvz/PPPMMn3zyCatWreIvf/lL3vuaL9f4zazJWrRoEe3ataNFixZ5NfVUZ+DAgUyePJlevXohieuvv77KTQV77LEHV1xxBQMGDKCkpIS9996bMWPGVFnX2LFjKy7u7rTTTowePZrVq1dzyimnsGzZMiKCCy64gC233JIrr7ySp59+mpKSEnr06MGRRx7JggV1P4K6S5cuXH755ey777507tyZHj16sMUWW6zzftdmg3jmbt++faM+/fHPnz+fLl26MG/evIJd3G1OGut4uVzWTbGO18yZM9l9990bbXv2teXLl7PZZpuxatUqBg4cyJAhQxg4cGCV+eoqI0kvRUTfyuObdVNP27Ztufrqq2nbtm2xQ9kgNNbxcrmsGx+v7Bk2bBi9e/dmzz33ZMcdd+S4445r0PU36xq/mdWfa/xNn2v8ZmaWFyd+M7OMceI3M8sYJ34zs4xx4jezJu+yyy5j4sSJPPTQQwwfPryosVx11VU88cQTRY1hffkHXGaWn2EN+yMihi3Le9aG6Ja54rGDLdavvvvLX/5yvZZvClzjN7Mm6+KLL6Znz55MmTKF/fffn5EjR3LWWWdVm3w//PBDBg4cSK9evejVqxf//Oc/mT17Nrvvvjs/+9nP6NOnDx988AEXX3wxe+65J3vttRfjxo0Dkl49DzzwwIp755999llWr17N4MGDK+a96aabgKT75/I+9bt168bVV19Nnz592GuvvSq6ZVi0aBGHHnooffr04Sc/+Qk77LADixcvbqSjVjcnfjNrskaMGMHIkSMZPHgwU6ZMoWfPnrzyyitV+uIHOPfccxkwYAAvv/wy06ZNY4899gCSB9j88Ic/ZPr06UydOrWii+YnnniCiy++mAULFvCnP/2Jww8/vGJaeT9B8+bNY8aMGbz66qucfvrp1cbYoUMHpk2bxllnncUNN9wAwDXXXMPBBx/MtGnTGDhwIHPmzCncQaoHJ34za9Ly7Zb5qaee4qyzzgKSni7L+7fZYYcd2G+//YCau2ju168fo0ePZtiwYbz66qu0bduWnXbaiVmzZnHOOecwYcIENt9882q3+93vfheAffbZh9mzZ1ds56STTgLgiCOOoF27dg1yLBqKE7+ZNUllZWX07t2bK664ghEjRnDUUUcxYcIEevfuzeeff573enK7Va6pp4IDDzyQSZMm0aVLF0499VTuuece2rVrx8svv0xpaSm33XZbxbMaKivvWrm8W+XattNUOPGbWZO0rt0yH3LIIdx+++0ArF69utqH4NTURfP777/P1ltvzY9//GPOOOMMpk2bxuLFi1mzZg3f+973Kp6Ila9vfvObjB8/Hkie3lb+JLemwnf1mFmTtS7dMt98882ceeaZjBo1ipKSEm6//fa1Hq8INXfRPHbsWEaMGMFGG23EZpttxj333MO8efM4/fTTWbNmDQC//e1v84776quv5uSTT2bcuHEMGDCATp06NalO9txJm5lVy5201d/KlSspKSmhZcuWTJ48mbPOOouysrIG3059O2lzjd/MrIHNmTOHE044gTVr1rDxxhtz1113FTuktTjxm5k1sF122YXp06cXO4wa+eKuNYjS0lJKS0uLHYaZ5aGgiV/SBZJekzRD0p8ltZbUXtL/SXo7/d+0bnA1swobwjXArFqfsilY4pfUBTgX6BsRewIlwEnAUODJiNgFeDIdNrMmpnXr1ixZssTJvwmKCJYsWULr1q3rtXyh2/hbAm0kfQVsAswHLgNK0+ljgYnApQWOw8zWUdeuXZk7dy6LFi0qdihWjdatW9O1a9d6LVuwxB8R8yTdAMwBPgcej4jHJW0TEQvSeRZI2rq65SWdCZwJsP322xcqTDOrwUYbbcSOO+5Y7DCsAArZ1NMOOBbYEegMbCrplHyXj4g7I6JvRPTt2LFjocI0M8ucQl7c/TbwXkQsioivgAeB/wA+lNQJIP3/UQFjMDOzSgqZ+OcA+0naRJKAQ4CZwMPAaek8pwF/LWAMZmZWSSHb+F+Q9AAwDVgFTAfuBDYDxks6g+TD4fhCxWBmZlUV9K6eiLgauLrS6JUktX8zMysC/3LXzKwWzfFX6U78ZmYZ48RvZpYxTvxmZhnjxG9mljFO/GZmGePEb2aWMU78Zk1Ic7x10JoeJ34zs4xx4jczyxgnfjOzjHHiNzPLGCd+M7OMKfQzd605GLZF3fPMXpH/vADDltU/HjNbL67xm5lljBO/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME7+ZWcbU2WWDpAOAsohYIekUoA9wc0S8X/DozJqJbkMfyWu+hbOW5D3/7OFHrVdMll351PhvBz6T1Au4BHgfuKegUZmZWcHkk/hXRUQAx5LU9G8G2hY2LDMzK5R8euf8t6TLgFOAAyWVABsVNiwzMyuUfGr8JwIrgTMiYiHQBRhR0KjMzKxgaq3xp7X7+yLi2+XjImIObuM3M9tg1Vrjj4jVJBd283y6hpmZNXX5tPF/Abwq6f+AFeUjI+LcgkVlZmYFk0/ifyT9MzOzZqDOxB8RYxsjEDMzaxw1Jn5J4yPiBEmvAlF5ekT0LGhkZmZWELXV+M9L/x9d35VL2hIYCexJ8uExBHgTGAd0A2YDJ0TEJ/XdhpmZrZsaE39ELEj/r0+fPDcDEyJikKSNgU2Ay4EnI2K4pKHAUODS9diGNQETB29a7BDMLE91/oBL0n6SpkhaLulLSaslfZrHcpsDBwKjACLiy4hYStL1Q/l1g7HAcfUN3szM1l0+v9z9PXAy8DbQBvgRcGsey+0ELAJGS5ouaaSkTYFtcr5NLAC2rlfkZmZWL3n1xx8R7wAlEbE6IkYDB+WxWEuSLpxvj4i9SX4DMDTfwCSdKWmqpKmLFi3KdzEzM6tDPon/s7R9vkzS9ZIuAPJp0J0LzI2IF9LhB0g+CD6U1Akg/f9RdQtHxJ0R0Tci+nbs2DGPzZmZWT7y+QHXKSQfEGcDFwDbAd+ra6GIWCjpA0m7RcSbwCHA6+nfacDw9P9f6xm7mdl6yeeBN83x4Ti13ce/C3ADsDPwKnBRRFyzjus/B/hj+o1hFnA6yYfIeElnAHOA4+sTuFlztO33hxc7BMuA2mr8d5P0wjkJOIbkgu5312XlEVEG9K1m0iHrsh4zM2s4tSX+thFxV/p6hKRpjRGQmZkVVm2Jv7WkvQGlw21yhyPCHwRmZhug2hL/AuDGnOGFOcMBHFyooMzMmqvS0lIAJk6cWLQYauuyIZ979c3MbAOT1w+4zMys+XDiNzPLGCd+M7OMqe0HXH1qW9B39ZiZbZhqu6vnv2qZ5rt6zMw2UL6rx8wsY/J5EMsmkn4h6c50eBdJ9X4co5mZFVc+F3dHA18C/5EOzwV+VbCIzMysoPJJ/DtHxPXAVwAR8Tlfd+NgZmYbmHwS/5eS2pBc0EXSzsDKgkZlZmYFk8+DWK4GJgDbSfojcAAwuJBBmZlZ4dSa+CW1ANqR9MO/H0kTz3kRsbgRYjMzswKoNfFHxBpJZ0fEeKDu546ZmVmTl08b//9JukjSdpLal/8VPDIzMyuIfNr4h6T//1/OuAB2avhwzMys0PJp4x8aEeMaKR4zMyuwWpt6ImINa9f0zcxsA+c2fjOzjHEbv5lZxtSZ+CNix8YIxMzMGkc+NX4k7Qn0AFqXj4uIewoVlJmZFU6diV/S1UApSeL/B3Ak8BzgxG9mtgHK5+LuIOAQYGFEnA70AloVNCozMyuYfBL/5+ltnaskbQ58hC/smpltsPJp458qaUvgLuAlYDnwYiGDMjOzwsnnrp6fpS/vkDQB2DwiXilsWGZmVii1NvVIailJ6evtgL5ASWMEZmZmhVFj4pf0Y5L2/PfT10+SXOi9X9KljRSfmZk1sNqaes4HdgbaAjOBHSJisaRNgCnAdYUPz8zMGlptif/LiPgE+ETSO+VP3YqIzyR92TjhmZlZQ6st8beRtDdJc9DG6Wulf61rWc7MzJqw2hL/AuDG9PXCnNflw2ZmtgGqMfFHxEENsQFJJcBUYF5EHJ126TwO6AbMBk5Im5TMzKwR5NVJ23o6j+Ti8Obp8FDgyYgYLmloOuy7hMxswzdsi7rnmb0i/3kBhi2rfzw1yKfLhnqT1BU4ChiZM/pYYGz6eixwXCFjMDOztdVY45d0QEQ8L6lVRKys5/p/B1xCcktouW0iYgFARCyQtHUN2z8TOBNg++23r+fmzczWz7bfH17sEBpcbTX+W9L/k+uzYklHAx9FxEv1WT4i7oyIvhHRt2PHjvVZhZmZVaO2Nv6vJI0Guki6pfLEiDi3jnUfABwj6Tskt39uLuk+4ENJndLafieSXwebWQGUlpYCMHHixKLGYU1LbTX+o4HHgC9IeuWs/FeriLgsIrpGRDfgJOCpiDgFeBg4LZ3tNOCv9Y7ezMzWWW23cy4m6ZdnZkS83IDbHA6Ml3QGMAc4vgHXbWZmdcjnds4lkv6XpOkmSB67eF5EzM13IxExEZiYvl5C8kQvMzMrgnxu5xxN0jzTGegC/C0dZ2ZmG6B8Ev/WETE6Ilalf2MA32ZjZraByifxL5J0iqSS9O8UYEmhAzMzs8LIJ/EPAU4g6ZhtAcnDWIYUMigzMyucfJ65Owc4phFiMTOzRlDQvnrMzKzpceLPqNLS0opfdZpZtjjxm5llTN798UvaD/gN0AoYEREPFSooM8tDQ/f9XoB+361pqq1b5m0jIvcRiz8nucgr4J/AQ4UNzczMCqG2Gv8dkl4iqd1/ASwFvg+sAT5thNjMzKwAamzjj4jjgDLg75JOBc4nSfqb4KdmmZltsGq9uBsRfwMOB7YEHgTejIhbImJRI8RmZmYFUGPil3SMpOeAp4AZJH3qD5T0Z0k7N1aAZmbWsGpr4/8VsD/QBvhHRPQHfi5pF+DXJB8EZma2gakt8S8jSe5tyHk8YkS8jZO+mdkGq7Y2/oEkF3JXkdzNY2ZmzUBdj168tRFjMTOzRuAuG8zMMsaJ38wsY/Luq8c2HN2GPlLnPAtnLcl73tmt1zskM2tCnPjNmrGJgzctdgjWBLmpx8wsY5z4zcwyxonfzCxjnPjNzDLGF3fNzBpRU7jg7hq/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWVMwX7AJWk74B5gW2ANcGdE3CypPTAO6AbMBk6IiE8KFYdVb9vvDy92CGZWJIWs8a8CLoyI3YH9gP8nqQcwFHgyInYBnkyHzcyskRQs8UfEgoiYlr7+NzAT6AIcC4xNZxsLHFeoGMzMrKpGaeOX1A3YG3gB2CYiFkDy4QBs3RgxmJlZouCJX9JmwF+A8yPi03VY7kxJUyVNXbRoUeECNDPLmIImfkkbkST9P0bEg+noDyV1Sqd3Aj6qbtmIuDMi+kZE344dOxYyTDOzTClY4pckYBQwMyJuzJn0MHBa+vo04K+FisHMzKoqZH/8BwCnAq9KKkvHXQ4MB8ZLOgOYAxxfwBjMzKySgiX+iHgOUA2TDynUds3MrHb+5a6ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWWME7+ZWcY48ZuZZYwTv5lZxjjxm5lljBO/mVnGOPGbmWVMURK/pCMkvSnpHUlDixGDmVlWNXril1QC3AYcCfQATpbUo7HjMDPLqmLU+PsD70TErIj4ErgfOLYIcZiZZVLLImyzC/BBzvBcYN/KM0k6EzgzHVwu6c1GiK2p6wAsbuyNqhArvaYgay2GopQJFKBcmk+ZgM+VcjtUN7IYib+6vYgqIyLuBO4sfDgbDklTI6JvseOwr7lMmiaXS+2K0dQzF9guZ7grML8IcZiZZVIxEv8UYBdJO0raGDgJeLgIcZiZZVKjN/VExCpJZwOPASXA3RHxWmPHsYFy01fT4zJpmlwutVBEleZ1MzNrxvzLXTOzjHHiNzPLGCf+JkDS8mrG/VzS65JekfSkpB0qTZ8gqYukUZJeTud7QNJmjRd5tlRXTjnTBkkKSb6FsJHVcP4MlrRIUln696NixNZUOfE3XdOBvhHRE3gAuL58gqQ2QPuImAdcEBG90vnmAGcXJdoMk9QWOBd4odix2FrGRUTv9G9ksYNpSpz4m6iIeDoiPksH/0Xye4dypcDEdL5PASQJaEM1P4azgruW5IP5i2IHYpYPJ/4NwxnAoznDRwITygckjQYWAt8Abm3c0LJN0t7AdhHx92LHYlV8L6cJdLu6Z88OJ/4mTtIpQF9gRM7oA4Dnygci4nSgMzATOLFRA8wwSS2Am4ALix2LVfE3oFvaBPoEMLbI8TQpTvxNmKRvA1cAx0TEynTcTsAHac+mFSJiNTAO+F6jB5pdbYE9gYmSZgP7AQ/7Am/xRcSS8nMGuAvYp5jxNDVO/E1U2oTw3yRJ/6OcSRXNPEp0L38N/CfwRmPHmlURsSwiOkREt4joRnIt5piImFrk0DJPUqecwWNIvg1bqhi9c1pVm0iamzN8I/AdYDPgf5KczpyIOAY4AjgnnU/AWEmbp69fBs5qtKizp0o5RcSNRYvGylV3/nSUdAywCvgYGFyMwJoqd9mwAZHUCnje3c2a2fpw4jczyxi38ZuZZYwTv5lZxjjxNyJJbSQ9I6lkPddTmvY/8pqkZ9Jxu+X0S1Im6VNJ56fTbpB0cAPsQrPQEOUg6RuSJktaKemiStOOkPSmpHckDa1m2YvSfn06pMN7SRpT31iao0KWkc8V39XT2IYAD6b33NeLpC2BPwBHRMQcSVsDRMSbQO90nhJgHvC/6WK3ktzL/FS9I29e1rscSO4UORc4LndkeuxvAw4leczoFEkPR8Tr6fTt0mlzypeJiFcldZW0fUTMwaCAZeRzxTX+xvYD4K9pjX1i+lPyNyT9Mb0PPx/fJzkh5gBUuse/3CHAuxHxfjrP+8BWkrZtiJ1oBta7HCLio4iYAnxVaVJ/4J2ImJX+yO5+4Nic6TcBl1C1T6W/kTyG1BKFLKNcmTxXnPgbiZLnC+8UEbPTUXsD5wM9gJ1IumFA0k2VvoaW/5U3GewKtEtPhpck/bCazZ0E/LnSuGnl28iyBiyHmnQBPsgZnpuOI72vfF5EvFzNclOBb9Vzt5qVRiijXJk8V9zU03g6AEtzhl+MiLkAksqAbsBzEXFBHetpSfLz80NIeuOcLOlfEfFWuq6NSX6peFml5T4i6c8n6xqqHGpSXW00JG1C0v3GYTUs5/L5WqHLiHRdmT1XnPgbz+dA65zhlTmvV5OWhaSbgIOqWf7+iBhOUoNcHBErgBWSJgG9gLfS+Y4EpkXEh5WWb53GkHUNVQ41mQvk9gTZFZgP7AzsCLyctlR0BaZJ6h8RC3H55Cp0GZXL7LnixN9IIuITSSWSWtcxX121mL8Cv5fUEtgY2Jek3bjcyVT96gpJE9H/rEPIzVIDlkNNpgC7SNqR5KLhScD3I+I1YOvymZR06tY3Ihano3YFZtRzm81KI5RRucyeK078jetx4Jsk/YfUS0TMlDQBeAVYA4yMiBkAaXPCocBPcpeRtBHQnaQd2RqgHNKLf1OBzYE16e2APSLiU0lnA48BJcDdadKvy0HAI/WNpxkqdBll+lxxlw2NSEmPmz+PiFMbebsDgT4RcWVjbrepKlY51ERJH0zPAN+MiHonuubE50ph+a6eRhQR04Gn1+dHKfXUEvivRt5mk1XEcqjJ9sBQJ/2v+VwpLNf4zcwyxjV+M7OMceI3M8sYJ35rciRdoaQDulfSX2Lum44/P70bo3y+f6R9FzXENpfXMX1LST/LGe4s6YEG2vZE+Tm91oic+K1JkbQ/cDTJnRU9gW/zdRcI5wMViT8ivhMRSxsptC2BisQfEfMjYlAjbbtBNKGL2VZkTvzW1HQi+WXySoCIWBwR8yWdS/Iz+qclPQ3Jj6AkdZDULe3Aa6SkGWlHXt+W9LyktyX1T+cfprW7550hqVvuxiVtJulJSdMkvSqpvIO14cDO6TeQEek2y38/0VrS6HT+6ZIOSscPlvSgpAlpHNfnexDS9T+bxjFN0n+k4+/NiYl0X49Jf/A0QtKU9JvST9LppZKelvQn4NV1KglrviLCf/5rMn8kD5gvI+mC4g/AgJxps4EOlYdJ+m5ZBexFUpl5CbibpN+cY4GH0vmHARflLD8D6Ja+Xp7+bwlsnr7uALyTrqcbMCNn2Yph4EJgdPr6GyRdLrcmecD3LGCLdPh9YLtq9nkiya94c8dtArROX+8CTE1fD8jZny2A99KYzwR+kY5vRfIDpB2BUmAFsGOxy9Z/TefPNX5rUiJiOUkndGcCi4Bxkgbnseh7EfFqRKwBXgOejIggqeV2W4cQBPxG0ivAEyQ9a25TxzLfBO5N43+DJMHvmk57MiKWRcQXwOvADnnGsRFwl6RXSboP6JGu/xmgu5LnMJwM/CWS+/8PA36YdmL2ArAVyQcGJJ2cvZfndi0D3GWDNTmRPHxjIjAxTXynAWPqWCy3I681OcNr+Pp9voq1mzer6wvmB0BHYJ+I+CrtU6fWPmOovkfO6uKq6GAsDxcAH5J0wNcC+CJn2r1pnCeRPLCkPIZzIuKxtQKTSklq/GYVXOO3JkXJY/F2yRnVm6QGDfBvoO16rH420CfdTh+SppDKtgA+SpP+QXxdQ69t25NIEjGSdiX5Je6b6xFneRwL0m8wp5L0+1NuDMmFbuLrfoAeA85K+5pB0q6SNl3PGKyZco3fmprNgFvT2zRXkbSxn5lOuxN4VNKCiKiuO966/IWvm0Om8HVX1rn+CPxN0lSSaw1vAETEkvRi8QzgUZLHK5b7A3BH+u1kFTA4IlYq74eqAfCIpPInRU0GLgf+Iul44Glyau0R8aGkmcBDOcuPJGnSmqZkw4uo9MhBs3LussFsA5P+luFVkltelxU7HtvwuKnHbAMi6dsk30JuddK3+nKN38wsY1zjNzPLGCd+M7OMceI3M8sYJ34zs4xx4jczyxgnfjOzjPn/dVoATH2cUAkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def significance_bar(start,end,height,displaystring,linewidth = 1.2,markersize = 8,boxpad  =0.3,fontsize = 15,color = 'k'):\n",
    "    # draw a line with downticks at the ends\n",
    "    plt.plot([start,end],[height]*2,'-',color = color,lw=linewidth,marker = TICKDOWN,markeredgewidth=linewidth,markersize = markersize)\n",
    "    # draw the text with a bounding box covering up the line\n",
    "    plt.text(0.5*(start+end),height,displaystring,ha = 'center',va='center',bbox=dict(facecolor='1.', edgecolor='none',boxstyle='Square,pad='+str(boxpad)),size = fontsize)\n",
    "\n",
    "n_L23 = data['L2/3']['total']\n",
    "n_L4 = data['L4']['total']\n",
    "n_L5 = data['L5']['total']\n",
    "\n",
    "plt.clf()\n",
    "i = 0\n",
    "tick_labels = []\n",
    "colors = {\"# no crossing\": 'tab:blue', \"# crossing\": 'tab:orange'}\n",
    "# in format # no crossing, # crossing, # no signal; l2/3, L4, L5\n",
    "\n",
    "length_of_previous_sl = i\n",
    "\n",
    "all_layers = []\n",
    "totals = []\n",
    "if n_L23 > 0:\n",
    "    all_layers.append('L2/3')\n",
    "    totals.append(n_L23)\n",
    "if n_L4 > 0:\n",
    "    all_layers.append('L4')\n",
    "    totals.append(n_L4)\n",
    "if n_L5 > 0:\n",
    "    all_layers.append('L5')\n",
    "    totals.append(n_L5)\n",
    "\n",
    "\n",
    "# in columns L23, L4, L5\n",
    "crossing_data = {\n",
    "    \"# no crossing\": [data[x]['n_barrel_no_cross'] for x in all_layers],\n",
    "    '# crossing': [data[x]['n_barrel_cross'] for x in all_layers],\n",
    "    #'# no signal': [data[x]['n_barrel_none'] for x in ['L2/3', 'L4', 'L5']],\n",
    "    #\"L2/3\": [3, 17, 30],\n",
    "    #\"L4\": [38, 20, 30],\n",
    "    #\"L5\": [8, 13, 1]\n",
    "}\n",
    "errbars = {}\n",
    "for k in crossing_data:\n",
    "    errbars[k] = [np.sqrt(x) for x in crossing_data[k]]\n",
    "    for j in range(len(totals)):\n",
    "        crossing_data[k][j] /= totals[j] / 100\n",
    "        errbars[k][j] /= totals[j] / 100\n",
    "\n",
    "labels = [x for x in crossing_data.keys()]\n",
    "tick_labels += [all_layers[j] + \"\\n(n=\"+str(totals[j]) + \")\"\n",
    "                            for j in range(len(all_layers))]\n",
    "\n",
    "print(errbars)\n",
    "print(crossing_data)\n",
    "\n",
    "# This is a binomial experiment of n independent trials.\n",
    "# We use the Binomial case of goodness of fit test to a \n",
    "# chi-squared distribution\n",
    "Zs = []\n",
    "for i in range(3):\n",
    "    layer = ['L2/3', 'L4', 'L5'][i]\n",
    "    n = [n_L23, n_L4, n_L5][i]\n",
    "    no_cr = crossing_data['# no crossing'][i]\n",
    "    cr = crossing_data['# crossing'][i]\n",
    "    k = n * cr / 100 # number of crosses\n",
    "    \n",
    "    null_hyp_p = 0.5\n",
    "    \n",
    "    Z = binom_test(k, n, null_hyp_p, alternative='greater')\n",
    "    if Z > 0.5:\n",
    "        Z = 1 - Z\n",
    "    print(\"Binomial test for layer\", layer, \":\", Z)\n",
    "    Zs.append(Z)\n",
    "    \n",
    "##########################\n",
    "for nm in crossing_data:\n",
    "    x = np.arange(len(all_layers))\n",
    "    # else:\n",
    "    plt.bar(x + (i-1) * 0.2, \n",
    "        crossing_data[nm],\n",
    "        0.2,\n",
    "        label=(nm if length_of_previous_sl == 0 else None),\n",
    "        color=colors[nm])\n",
    "    plt.errorbar(x + (i-1) * 0.2, \n",
    "                crossing_data[nm],\n",
    "                yerr=errbars[nm],\n",
    "                fmt='none',\n",
    "                color='black')\n",
    "    i += 1\n",
    "\n",
    "for i in range(3):\n",
    "    if Zs[i] < 0.001:\n",
    "        significance_bar(i-.5, i+0.25, 93, \"P<0.001\", fontsize=10)\n",
    "    \n",
    "plt.title(\"Frequency of Signal Pattern in Barrel Pairs by Layer\")\n",
    "plt.xlabel(\"Stimulation Layer\")\n",
    "plt.ylabel(\"% of Barrel Pairs\")\n",
    "#plt.xticks(rotation='vertical')\n",
    "plt.xticks(np.arange(len(tick_labels)), tick_labels)\n",
    "plt.legend(bbox_to_anchor=(0.82, 1.0), loc='upper center')\n",
    "plt.savefig(data_dir + \"results_count_oct2023.png\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "d7e8f1e3",
   "metadata": {},
   "outputs": [],
   "source": [
    "# use Slice Edge and Medial Dir metadata\n",
    "# to tell which barrel is more medial\n",
    "# Mark this as Medial Crossing or Medial no crossing or lateral etc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "579383fb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Latency between barrel pairs in L2/3\n",
      "\t 0.690 +/- 0.619 (n = 26)\n",
      "\t [1.3583 1.4885 0.2138 0.0093 0.1229 0.0447 1.1051 0.3143 0.1379 1.2641\n",
      " 0.5673 2.4054 1.6877 0.7967 0.2299 0.7214 0.2692 1.2911 0.5468 0.0567\n",
      " 0.1857 0.145  0.1196 1.252  1.0902 0.5346]\n",
      "Latency between barrel pairs in L4\n",
      "\t 1.028 +/- 1.001 (n = 25)\n",
      "\t [2.7731 0.0717 1.222  0.1611 1.3073 1.6698 0.4578 0.695  1.9031 0.3503\n",
      " 0.7656 0.9262 0.957  1.3192 0.0676 2.1705 0.1656 0.7253 0.853  0.543\n",
      " 0.3128 4.6423 0.1724 0.4323 1.0368]\n",
      "Latency between barrel pairs in L5\n",
      "\t 0.764 +/- 0.546 (n = 6)\n",
      "\t [0.1239 0.5704 1.4949 0.5802 1.5165 0.3024]\n"
     ]
    }
   ],
   "source": [
    "for stim_loc in ['L2/3', 'L4', 'L5']:\n",
    "    df2 = df[df['Stim Location'].str.contains(stim_loc)]\n",
    "    dfd1 = df2['abs(delta 1)'].dropna()\n",
    "    dfd2 = df2['abs(delta 2)'].dropna()\n",
    "\n",
    "    print(\"Latency between barrel pairs in\", stim_loc)\n",
    "    dfds = np.concatenate([dfd1, dfd2])\n",
    "    print(\"\\t\", \n",
    "          str(np.mean(dfds))[:5], \n",
    "          \"+/-\", \n",
    "          str(np.std(dfds))[:5], \n",
    "          \"(n =\", \n",
    "          str(len(dfds)) + \")\")\n",
    "    print(\"\\t\", dfds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "d5226216",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Latency between anterior (rostral) barrel pairs in L2/3\n",
      "\t 0.373 +/- 0.430 (n = 10)\n",
      "\t [0.0093 0.1229 0.0447 0.3143 0.1379 1.2641 0.5468 0.0567 0.145  1.0902]\n",
      "Latency between posterior (caudal) barrel pairs in L2/3\n",
      "\t 0.932 +/- 0.701 (n = 8)\n",
      "\t [1.1051 0.5673 2.4054 1.2911 0.1857 0.1196 1.252  0.5346]\n",
      "Latency between anterior (rostral) barrel pairs in L4\n",
      "\t 1.258 +/- 1.284 (n = 9)\n",
      "\t [1.3073 0.4578 0.695  1.9031 0.7253 0.853  0.3128 4.6423 0.4323]\n",
      "Latency between posterior (caudal) barrel pairs in L4\n",
      "\t 0.781 +/- 0.680 (n = 9)\n",
      "\t [0.1611 1.6698 0.3503 0.7656 2.1705 0.1656 0.543  0.1724 1.0368]\n",
      "KS test for ant layer L2/3 and ant layer L4 is significant with p-value  0.03026694667561547\n",
      "KS test for ant layer L4 and ant layer L2/3 is significant with p-value  0.03026694667561547\n",
      "KS test for ant layer L2/3 and post layer L4 is significant with p-value  0.03026694667561547\n",
      "KS test for post layer L4 and ant layer L2/3 is significant with p-value  0.03026694667561547\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAEuCAYAAACESglMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAvQUlEQVR4nO3deZhcVZ3/8fcne4SkUyFIgBCCYQtkmgBRXECJsokjaESRYRiWODjOjAtj1BkZJFEZdIiKCy4MQXAZdMQWkB8uoxIUB4EEYxtIUGRpEggJ0OkESAJJvr8/7q3kptNLdXdVV9Wtz+t56ulbdz11qvpbp849iyICMzPLnyHVToCZmVWGA7yZWU45wJuZ5ZQDvJlZTjnAm5nllAO8mVlOOcDboJP0PklPSXpO0h7VTk+lSDpP0p0VOO9zkl5R7vNWkqSvS7qk2uloNA7wNULSo5JO6MdxiyS9pxJpqgRJw4HPAydFxO4R8Uyn7VMkhaRhvZznV73tl25/Pg2IT0u6QdK4sryQAcq8zufSx6OS/rWUY9N8e7iP1ztZ0q8lbZC0VtIdkk7rX+r7LiL+ISI+NVjXs4QDvA22vYBRwP39PYGks4EevwAyjoiI3YFXAAVgXj+vucv1evsSKtG4NH1nAZ+QdMpATiZpaBfrzgB+AHwLmETyHnwCeGs35yjH67Ia4ABf4yQVJN2alrra0+VJ6bbLgOOAr6SlwK+k6w+V9L+SnpX0oKR3Zc53naSrJP2/tDR3t6Spme2HZ459StLHJU2U9EK2OkXS0WmahneR5pGSrpT0RPq4Ml13MPBguts6Sb/qR340AZcCH+3LcRGxHrgFOCxzrvMlLU/z4WFJ781sO17SSkkfk7Qa+KakeZJulPQdSeuB8yQ1SVoo6UlJqyR9uqsgW0L67iL50psu6VWS7pK0Lj3vVySNyKQtJB2YLl8n6WuSbpP0PDCrU36J5BfTpyLimojoiIhtEXFHRPx9us95kn4r6QuSngXmpa/rW+l7/Jikf5c0JN3/wPQXQEf6y+j7xWul51iTbmuVND2Tzk93ytsPp/s+Ken8TJr3kPRjSesl3Zvm6Z29XcO6EBF+1MADeBQ4oYv1ewDvAF4GjCEpid2U2b4IeE/m+W7A48D5JKXco4CngcPT7dcBzwKvSrd/F/heum0M8CTwYZJS9hjgmHTbbcD7Mtf5AvDlbl7LJ4HfAS8H9gT+jyTAAEwBAhjWzbG9bb8KuKi3/dJ9AzgwXS4APwc+mdn+FmAqIOANwAvAUem244EtwGeBkcBoktL/S8DbSApHo4GbgG+k+f5y4B7gvek5zgPu7O11ptd/XXr9NwFHA69Ot00BlgMf6uZ1XQd0pMcPAUZ1us6h6f4H9JBP56Wv9f3pNUeTlPZvTj8DU4A/AXPS/W8ALi5eDzg2XX8ysAQYl76macDemXR+ulPefhIYDpyavvZCuv176eNlJF/Ijxfzsadr+NHFe1vtBPiRvhHdBPgu9psBtGeeL2LnAH8m8JtOx3wDuDRdvg64JrPtVGBFunwW8Pturnsm8Nt0eSiwGnhVN/v+BTg18/xk4NF0eXtg6+bYbrcDM4Gl7Ah8pQT49cA6YCuwAti3h/1vAj6YLh8PvEgmYJIE+F9nnu8FbAZGZ9adBdyeLp9H7wF+HdBOEsQ/0M2+HwJ+1Ol1ZQP8t3p4Ta9L9x/Vwz7nAW2Z50PT13VYZt17gUXp8reAq4FJnc7zRpIvglcDQzptu46dA/zG7HsHrEmPG0ryJXpIZtun2RHgu72GH7s+XEVT4yS9TNI30p/J64FfA+N6qAbYHzgm/Xm/TtI64GxgYmaf1ZnlF4Dd0+X9SIJzV24GDlPSeuNEoCMi7ulm332AxzLPH0vX9VtaPfBVkgC8pQ+HHhUR40hKml8DfiNpVHrON0v6XVodtY7ky25C5ti1EbGp0/kezyzvT1ICfTKT198gKcmXakJEFCJiWkR8KU3XwUqq4lan7/l/dEpXZ4/3sK14E3vvXtKRPccEYAS7vof7pssfJSk93yPpfkkXAETEr4CvkPzKekrS1ZLGdpeuTu9j8XO4J8kXeDY925f7eI2G5wBf+z4MHEJSVTIWeH26XunfzsOBPg7cERHjMo/dI+J9JVzrcZIqi12kge5/SL4szgG+3cN5niAJfkWT03UDMZakBP/9tE783nT9SknH9XZwRLwEXAMcQFLPPRL4IbAA2Cv9EriNHfkKu+Zt53WPk5R0J2TyemxEHN63l7aLr5H82jgofc8/3ildPaWpswfTdL6jl2tmz/E0SSm683u4CiAiVkfE30fEPiQl+68W7wlExJci4mjgcOBg4CO9XLeztSTVN5My6/bbKaEDv0bDcICvLcMljco8hpHUgW4kuSk5nuQGY9ZTJC1Eim4FDpZ0jqTh6eOVkqaVcP1bgYmSPqTkpugYScdktn+L5Of8acB3ejjPDcC/S9pT0gSSFhs97d+Vkdm8ADaQ/AqYkT5OTfc7Gri7t5Olv3jOJ8nLh0lKqCNJA4qkNwMn9SWBEfEkSb3+5ySNlTRE0lRJb+jLebowhqRq6TlJhwKlfDl3l8YA/gW4RMlN5WI6j5V0dTfHbCX5Mr8s/Qzsn57jOwCS3qn0Rj9J9VIAW9PP2TFKbrw/D2wiqRrrS3q3Ai0kN3pflr7+vytuL8c1GokDfG25jSQAFR/zgCtJbno9TXLj8qedjvkicIaSFjZfiogNJIHq3SSl5tXsuFHYo/TYE0maz60G/kymVUZE/BbYBtwXEY/2cKpPA4uBVuCPwH3pur54jp3zYlZaclwdEatJAjPAUxHxYg/n+YOk50gC0bnA2yPi2fS1foAkkLUDf0PSyqav/o7ky+KB9Dw30nt1SG/mpunZAPwX8P2BnCwibiS5h3IByWfiKZL34+YeDns/SQB9GLgT+G/g2nTbK4G703y9haTa7BGSX1n/RZIPj5FUDy3oR5L/GWgi+Qx+m6TAsDndVq5rNASlNy7MSqKkaeN/R8Q11U6LNQZJnwUmRsS51U5LvXEJ3kom6ZUkzS4HVKI064mSfhzNaZv3VwFzgB9VO131yD3WrCSSridp//3BtHrDrFLGkFTL7EPSfPJz9FydZN1wFY2ZWU65isbMLKcc4M3Mcqqm6uAnTJgQU6ZMqXYyzMzqxpIlS56OiD272lZTAX7KlCksXry42skwM6sbkh7rbpuraMzMcsoB3swspxzgzcxyygHezCynHODNzHLKAd6q4oYbbuDAA6czZMhQCoXpvOMdN9DaWu1UmeVLTTWTtMZwww038JGPXMxhhy3klFOOZf36O7n55jmsWQNXXXUWzc3VTqFZPrgEb4Pusssu401vWsi0abOYMGE4r3jFLN72toU88MBltLRUO3Vm+eEAb4Nu+fLlSMfS1LRj3eTJx9Levpy2tuqlyyxvHOBt0E2bNo2IO+no2LGure1OCoVpTJ5cvXSZ5Y0DvA26iy++mF/+cg7Ll9/O00+/xMMP385NN83hsMMuZvbsaqfOLD98k9UG3VlnnQXAJZe8n1/8YjlNTdM44YTLuPRS32A1K6eamvBj5syZ4cHGzMxKJ2lJRMzsapuraMzMcsoB3swspxzgzcxyygHezCynHODNzHLKAd7MLKcc4M3McsoB3swspxzgzcxyygHezCynHODNzHLKAd7MLKcc4M3McsoB3swspxzgzcxyygHezCynHODNzHLKAd7MLKcqHuAlDZX0e0m3VvpaZma2w2CU4D8ILB+E65iZWUZFA7ykScBbgGsqeR0zM9tVpUvwVwIfBbZV+DpmZtZJxQK8pL8G1kTEkl72u1DSYkmL165dW6nkmJk1nEqW4F8HnCbpUeB7wBslfafzThFxdUTMjIiZe+65ZwWTY2bWWCoW4CPi3yJiUkRMAd4N/Coi/rZS1zMzs525HbyZWU4NG4yLRMQiYNFgXMvMzBIuwZuZ5ZQDvJlZTjnAm5nllAO8mVlOOcCbmeWUA7yZWU45wJuZ5ZQDvJlZTjnAm5nllAO8mVlOOcCbmeWUA7yZWU45wJuZ5ZQDvJlZTjnAm5nllAO8mVlOOcCbmeWUA7yZWU45wNegoUOHMmPGDKZPn8473/lOXnjhhZKOe+SRRzjmmGM46KCDOPPMM3nxxRe73O/666/noIMO4qCDDuL666/fvv68887jgAMOYMaMGcyYMYOlS5eW4+WYWZU4wNeg0aNHs3TpUpYtW8aIESP4+te/3u2+L774Is8//zwAH/vYx7jooov485//TKFQYOHChbvs/+yzzzJ//nzuvvtu7rnnHubPn097e/v27VdccQVLly5l6dKlzJgxo+yvzcwGjwN8jTvuuON46KGHdlm/fPlyPvzhD3PIIYfwpz/9iYjgV7/6FWeccQYA5557LjfddNMux/3sZz/jxBNPZPz48RQKBU488UR++tOfVvplmFkVOMDXsC1btvCTn/yEv/qrvwLg+eef55vf/CbHHnss73nPe5g2bRqtra0ceeSRPPPMM4wbN45hw4YBMGnSJFatWrXLOVetWsV+++23/Xnn/S6++GKam5u56KKL2Lx5c4VfoZlV0rBqJ6CStm7dylNPPVXtZPRqn3322en5xo0bt1ePHHfcccyZMweAvffem+bmZq655hoOPfTQnY6JiF3OK2mXdT3td/nllzNx4kRefPFFLrzwQj772c/yiU98Yqd9n3jiidJfWBXttddeDB06tNrJMKuqXAf4v/zlLxxyyCHVTkavOgfdYh18ZzfeeCMLFy7k7W9/O2eddRbnnnsu+++/PwATJkxg3bp1bNmyhWHDhrFy5cpdvjggKbEvWrRo+/OVK1dy/PHHA8kXCMDIkSM5//zzWbBgwS7H77vvvv18lYPrwQcf5OCDD652MsyqKtcBfvTo0QAsWbKEiRMnVjk1A3fSSSdx0kkn8cwzz/Cd73yH008/nQkTJnDNNdcwZcoUZs2axY033si73/1urr/+ek4//fRdznHyySfz8Y9/fPuN1Z///OdcfvnlADz55JPsvffeRAQ33XQT06dP3+X4rqp9asnq1as5+uijt7/3Zg0tImrmcfTRR0c5rVq1KoBYtWpVWc9babvttlvJ+959993R1tYWERF/+ctf4pWvfGVMnTo1zjjjjNi0aVNERNx7770xZ86c7ccsXLgwpk6dGlOnTo1rr712+/pZs2bF9OnT4/DDD4+zzz47NmzYUKZXNHjq9T036y9gcXQTUxVd1MlWy8yZM2Px4sVlO98TTzzBvvvuy6pVq7qsrrD88XtujUbSkoiY2dU2t6IxM8spB3gzs5xygDczyykHeDOznMp1gB8zZgyXXnopY8aMqXZSbJD4PTfbIdetaMzM8s6taMzMGlDFArykUZLukfQHSfdLml+pa5mZ2a56HapA0hDgCGAfYCNwf0SUMoLXZuCNEfGcpOHAnZJ+EhG/G1CKzcysJN0GeElTgY8BJwB/BtYCo4CDJb0AfAO4PiK2dXV82oX2ufTp8PRROxX+ZmY511MJ/tPA14D3Rqc7sZJeDvwNcA5wfRfHFvcbCiwBDgSuioi7B5xiMzMrSbcBPiLO6mHbGuDK3k4eEVuBGZLGAT+SND0ilmX3kXQhcCHA5MmTS0u1mZn1qtebrJLeKWlMuvzvklokHdWXi0TEOmARcEoX266OiJkRMXPPPffsy2nNzKwHpbSiuSQiNkg6FjiZpErma70dJGnPtOSOpNEkdfkrBpBWMzPrg1IC/Nb071uAr0XEzcCIEo7bG7hdUitwL/C/EXFr/5JpZmZ9VcqMTqskfYOkBP5ZSSMp4YshIlqBIweYPjMz66dSSvDvAn4GnJLWpY8HPlLJRJmZ2cD1WoKPiBck3Q7sl7m5+nRlk2VmZgNVSk/WTwHnAX9hR0elAN5YuWSZmdlAlVIH/y5gakS8WOnEmJlZ+ZRSB78MGFfhdJiZWZmVUoK/HPi9pGUkA4gBEBGnVSxVZmY2YKUE+OuBzwJ/BLocWMzMzGpPKQH+6Yj4UsVTYmZmZVVKgF8i6XLgFnauormvYqkyM7MBKyXAF3ujvjqzzs0kzcxqXCkdnWYNRkLMzKy8um0mKelv0+n6uts+NR1h0szMalBPJfg9SJpHLiGZlak4Zd+BwBtIhiv414qn0MzM+qWnGZ2+KOkrJHXtrwOaSSbdXg6cExFtg5NEMzPrjx7r4NMp9/43fZiZWR0pZagCMzOrQw7wZmY5VcpwwUPTqhozq3Gtq1tpWdFCW0cbk5smM/vQ2TRPbK52sqxKSinBPyTpCkmHVTw1ZtZvratbWXDXAto3tjNp7CTaN7az4K4FtK5urXbSrEpKCfDNwJ+AayT9TtKFksZWOF1m1kctK1oojCpQGF1giIZQGF2gMKpAy4qWaifNqqSUybM3RMR/RcRrgY8ClwJPSrpe0oEVT6GZlaSto42mUU07rWsa1URbh1s0dyaprI9a1WuAlzRU0mmSfgR8Efgc8Argx8BtFU6fmZVoctNkOjZ17LSuY1MHk5smVylFtSsiSnqUum+tKqWK5s/A6cAVEXFkRHw+Ip6KiBuBn1Y2eWZWqtmHzqZ9UzvtG9vZFtto39hO+6Z2Zh86u9pJsypRb98+knaPiOcGIzEzZ86MxYsXD8alzHLJrWjKS1JNl9ABJC2JiJldbStluOCrJH0wItalJysAn4uIC8qYRjMrg+aJzQ7otl1JrWiKwR0gItrZMUa8mZnVqFIC/JC01A6ApPGUVvI3M7MqKiVQfw74P0k3ps/fCVxWuSSZmVk5lDKj07fSMeFnAQJmR8QDFU+ZmVkVZG9UF5/X632NUgcbWwG0ADcDz0lyw1ozy53Owz0AdT3cQymDjb2fpPfqU8BWklJ8kAxhYGaWG9nhHoqKwz3UYym+lDr4DwKHRMQzlU6MmXWt3N3ha71td7W0dbRtL7kX1fNwD6VU0TwOdPS6l5lVTKN0ra+2vA33UEqAfxhYJOnfJP1L8VHphJmZDbbOwz0AdT3cQykBvo1kTtYRwJjMo0eS9pN0u6Tlku6X9MGBJdXMrLKaJzYz9zVzKYwusHL9SgDmvmZuXda/Qwlj0WzfUdotIp4v+cTS3sDeEXGfpDHAEuBtPTWx9Fg0ZgNTD2On1JN6yM+exqIpZbjg10h6AFiePj9C0ld7Oy4inoyI+9LlDenx+/Yp5WZm1m+lVNFcCZwMPAMQEX8AXt+Xi0iaQjJ+zd19S56ZmfVXSR2dIuLxTqtKnoRb0u7AD4EPRcT6LrZfKGmxpMVr164t9bRmZtaLUtrBPy7ptUBIGgF8gLS6pjeShpME9+9GRJcTQ0bE1cDVkNTBl5Rqq0utrdDSAm1tMHkyzJ4NzfV578qsLpRSgv8H4J9I6s9XAjOAf+ztICU9MxYCyyPi8wNIo+VAayssWADt7TBpUvJ3wYJkvZlVRikB/pCIODsi9oqIl0fE3wLTSjjudcA5wBslLU0fpw4otVa3WlqgUEgeQ4bsWG7p8nedmZVDKVU0XwaOKmHdTiLiTpJxa6xBldq9fv780s5X683VzGpNtwFe0muA1wJ7duq5OhYYWumEWf3LBuR585JqmcKOMZxob4cvfan22xmb1aueqmhGALuTfAlke7CuB86ofNIsT2bPTgJ6ezts27Zj2cwqp9eerJL2j4jHBiMx7smab121ojniCJfgy6keel7Wk3rIz556spZSB/+CpCuAw4FRxZUR8cYypc8aRHOzm0WaDaZSWtF8l2RGpwOA+cCjwL0VTJOZmZVBKQF+j4hYCLwUEXdExAXAqyucLjMzG6BSAvxL6d8nJb1F0pHApJ4OMDOrZ62rW5m3aB4A8xbNq9s5WUsJ8J+W1AR8GJgLXAN8qJKJMjOrluzE2wDtG9vrduLtXgN8RNwaER0RsSwiZkXE0cDUQUibmdmg6zzxdmF0YfvE2/WmpNEku+Ap+8wsl9o62mga1bTTunqdeLu/Ad5DEJhZLuVp4u3+BvjabvlvZtZP2Ym3IamDr9eJt7sN8JI2SFrfxWMDsM8gptHMbNBkJ96GpA6+XifeLnnS7cHgoQoaTz10Ba8nzs/yqof8HNCk22ZmVp8c4M3McsoB3swspxzgzcxyygHezCynHODNzHLKAd7MLKdKmdHJrCJa08H5LrhgxxR+nvHJrHxcgreqaG2FBQuS5UmTkgm4FyzYEfSt74rD2V5w8wV1PYa5lY8DfIluuOEGDjxwOkOGDKVQmM473nGDg9EAtLRAIekJzpAhyXKhkKy3viuOYQ4waeykuh7D3MrHVTQluOGGG/jIRy7msMMWcsopx7J+/Z3cfPMc1qyBq646q+GqFcaPH097e3vZzjd/vjo97995CoUCzz77bBlSVH9aVrSwZesWAH784I9pGtXEPrvvQ8uKlrocQ8XKwwG+BJdddhlvetNCxo2bRaEAEybM4m1vW8itt76flpbGC/Dt7e0DHp9j3rykWqZYik/OmzyfN69/55QadxTrpU8u5eF1DwMwduRYNr60kWVrl/HCSy9UOWVWTa6iKcHy5cuRjqUpMwfA5MnH0t6+nLb6mwOgJsyenQT09nbYtm3H8uz6G5G1JqzbvI4h6b+zJEYPH80QhrBu87rqJqwKxo8fj6SyPICynGf8+PFVyQsH+BJMmzaNiDvpyMwB0NZ2J4XCNCbX3xwANaG5GebOTUrsK1cmf+fOdSua/ho3chzbYhsAEcHGlzayLbYxbuS46iasCoq/MGvpUc4qzb5wFU0JLr74Yj7ykTkcdthCDj54Rx18c/NlLnEOQHOzA3q5zNh7BruN2I0/8kfWb15P06gmDhx/IAftcVC1k2ZV5ABfgrPOOguASy55P7/4xXKamqZxwgmXcemljVf/brVp9qGzt7eieeshb6VjU0fdzkJk5eMJP6zPanUShFpN12BpXd3KEXsfwfk3nc/kpsnbg3vLihbaOtq2r8t7q5pa/BxUMk09TfjhAN+D1takXXZbm3taZtXiPxDUbroGUzYPim3jC6MKNI1q2l6qr9fp50pVi5+DagV432TtRrGnZXu7e1pafWpZ0UJhVIHC6AJDNITC6AKFUQVaVrg3WaOoWB28pGuBvwbWRMT0Sl2nUoo9LYvttIt/W1pcirf60NbRxqSxk3Za1zSqibYOt+3tSevq1txUa1XyJut1wFeAb1XwGmVVakeZUnta1trPRGssk5sm076xncLoHb3JOjZ1MLnJbXu7k63Wyg75UK/VWhUL8BHxa0lTKnX+SsgG5O56Wn7pS7VXv2fWlWzLmmwd/Jwj51Q5ZbUrW60FbP9br0M+uJlkN2bP3jHaYVMTdHQkAd6s3Mo9tk9vv0S/zbdLOk8jju2Tt2qtqt9klXShpMWSFq9du7baydmuu56WZuVWiz0vq9n7spomN02mY1PHTuvquVqr6iX4iLgauBqSZpJVTs5O3NPSrLHkrVqr6iV4M7Na0TyxmbmvmUthdIGV61dSGF2o2xusUNlmkjcAxwMTJK0ELo2IhZW6nplZOTRPbK7bgN5ZJVvRnFWpc1t+ufewWfm4isZqhnsPm5WXA7zVjGzvYc/TajZwVW9FY1bU1paU3LOamvCsWdYncelYmNfU+46DKC4dW5XrOsBbzZg8edfewx0deNYs6xPNX9/n3uaVHnlTEjFvwKfpMwd4qxmzZ8PFF8OaNbB5M4wcCS9/OVx2WbVTZnmXtyEKihzgrSq6ai0D0Lng5WF/bDDkbYiCIgd4G3TF1jKFws6tZXbbDaZOhZmZqQva2z1Es1VeXkfedIC3PhvoTaxm4FuvyKwopI+uFNfPKzFdZv2QtyEKijxlXx/V4nRgg22geXDBBUnJfUimke62bXDLLfCGN+w6RHOhkAzfXOl0VUutprtW09Wb/qa7khN9VGvKPpfgbdB111rm1a/eMSRzdojmOfVdiLI6kachCorc0ckG3ezZSeBub09K7sXlf/zHrododv27Wf+4BG+DrjjWfrYVzZw5OwK5A3rp8jR/qJWfA3w/zJvnwbAGymPtD1zn+UP/9PSfOOdH53DAuAOYsfcMB/syqecvUVfR9EFx0CsPhmW1INs5Z83za7j/6fsRon1T+/bJoltX+8M5EMUv0faN7TtNwl0v+eoA3wfFQa88GJbVgraONppGJc1VVzy9glFDR9E0qon1m9dTGF2gMKpAywp/OAci+yU6REPqLl9dRdMHXQ165cGwbKD626/g2p42rv7CjuVFX+h+vx64X0H993B1gO+Drga98mBYA+MJPvo3OBbsXAe/9KmldGzsIAheu99r2Wv3vbb3zJx3/Lz+patKA2TVknrv4dowHZ3Gjx9fc7PEFwoFnn322Wono88G2mmjGNSXLoVHHoHDD4cDD9zR7r2/TSMbrWMO7LgBuPTJpTzS8QjT95zO1PFTyzIaYiPmZ2flGmWyWh2dGibAlyuDy1nibMR/oOw4NEuXJkE9Al7zGpg4sW89V8uZrmoq22ezzK09Gj0/i8qRrw7w1EeAL6daTFMpJFU7CV1q1F9ElVKr6epNLX4+K/nZ9FAFVlblGodm0SK44w7x+tcH69fD6acPrARvBgP7fBYVS+3zZ83n0tsvrau271luJmmDavLkpFoGYNq05G9HB4wdu2PIguLY8GbVkG37DtRd2/csl+D7yK0+Bmb6dPjUp+Cll2DPPZN1ETv6FGSHLDCrhjzN7uQA3wfdTVThAbFK09qaDAl8+OGwahWsXZusv+QSOOOM6qbNrKje275nuYqmD1padpQ03ZO174r5d/DBMGsWvOtdyfply6qbLrOsyU2T6djUsdO6emr7nuUSfB+0tSUl9yz3ZO1eqa0Z5s8X8+f3vl89tuiw+pOd3QmSOvh6nd3JAb4Pupuowj1Zu9Y5IM+bt2v+udWM1Zrmic3Mfc3c7ePNFEYXmHPknLqrfwe3g++TbB18dsahRut52V/lzr+86OvnYLCGr220z2dX6iEP3NEJ92StFW6FtKu+fA7K1XW+3Omqd9kvzRFDRyDEkxue5Oa/uZnT//v0mh5f3x2d+iEbiEaMAAk2b3ZQGihP9DEweWrCVyuyX5rDhwznjkfvYNNLmxg+bDgAD7c/zG4jdmPBXQsq8kVaSQ7wXchWJQwfDj/7GWzYAHvsAfffD0uWwGWXOVDZ4OupCV89zzxUTdkvzUWPLmLsyLGs37yeF196EYDRw0ezasMqjtjriLr7InWA70K2OeSttybBfdgw2LQp2f7QQ/C1ryUPs3Iox/gp3+SbOz2fTwlNk3pQyN4Nz7Hsl2bHpg7GjhzL1thKbEuqp0YNG0XHpo66bAvfMAG+L5MqzANIP9vz/rq3HQeYJmt4fb3B2lUd/G7Dd2PE0BE7jVs+0PHgG0V2zPemUU1sfGkjQzV0ey+hTVs2bc/remsL75usXcg25/v612HrVnjhhaR7/YgRMHp0Ul2zYsXgpcmsqKuqmCvvvpJJYycxRDv6Lm6Lbaxcv5JrT+9x7qeGVO4RJ6v5f1y1m6ySTgG+CAwFromIz1TyeuUye3ZSBw9JkH/ooSTIjxyZjJvS0ZEE+dZW18Pb4Gue2LxLPXC9zzw02DoH5O5a0azbvI5xI8fVdCuanlQswEsaClwFnAisBO6VdEtEPFCpa5aQpgEdX6yDB3j8cTjiiIGlp1HqOK3ysr0vs1U39dj7shq6+tLMg0qORfMq4KGIeDgiXgS+B5xewev1KCL6/PjDH4KpU4OmpmDs2GD06ORbf9q04IgjgvPP7/s5s496nJzCalOx92VhdIGV61dSGF2ouyZ9Vn6VrKLZF3g883wlcEwFr1d2LS1JXfvGjcnz4cPhsceSKppJkzxEgdWWvJZCrf8qGeC7qg/Z5U6EpAuBCwEmVzlillqF88QT4okn4J576HGQLN9ANbNqqmSAXwnsl3k+CXii804RcTVwNSStaCqYnl51NzjW5s2weDGsXJncbN1/f7j2Wt9gNbPaVsk6+HuBgyQdIGkE8G7glgper+xmz04C/MiRcOqpcOaZ8OY3O7ibWX2oWICPiC3APwM/A5YD/xMR91fqepXQ3JyMdFgoJKX3QsEjH5pZ/ahoO/iIuA24rZLXqDQPjmVm9cpT9pmZ5ZQDvJlZTjnAm5nllAO8mVlOOcCbmeVUTQ0XLGkt8Fi109GLCcDT1U5Ejjg/y8v5WV71kJ/7R8SeXW2oqQBfDyQt7m7sZes752d5OT/Lq97z01U0ZmY55QBvZpZTDvB9d3W1E5Azzs/ycn6WV13np+vgzcxyyiV4M7OccoDPkPRcF+v+RdIDklol/VLS/p22/1TSvpIWSvpDut+NknYfvJTXHudl5XWVx5ltZ0gKSXXbAmSw5TE/HeB793tgZkQ0AzcC/1ncIGk0MD4iVgEXRcQR6X5tJEMl286cl4NA0hjgA8Dd1U5LHtRzfjrA9yIibo+IF9KnvyOZmaroeGBRut96ACXz/o2mi+kJG53zctB8iuTLc1O1E5ITdZufDvB9Mwf4Seb5m4GfFp9I+iawGjgU+PLgJq3uOC8rQNKRwH4RcWu105IH9Z6fDvAlkvS3wEzgiszq1wF3Fp9ExPnAPiQzWJ05qAmsI87LypA0BPgC8OFqpyUP8pCfDvAlkHQCcDFwWkRsTte9Ang8Il7M7hsRW4HvA+8Y9ITWAedlRY0BpgOLJD0KvBq4pd5uDNaQus9PB/hepD/RvkESkNZkNm2vUlDiwOIy8FZgxWCntdY5LysrIjoiYkJETImIKST3OU6LiMVVTlpdykN+VnRO1jr0MkkrM88/D5wK7A78IIk3tEXEacApwPvT/QRcL2lsuvwH4H2Dlura5LysvF3yOCI+X7XU1L/c5ad7svaDpJHAb+t5lLla4bw0qxwHeDOznHIdvJlZTjVcgJe0VdJSScsk/VjSuH6c43hJr+3ncXXZnrY7kiZK+p6kv6TDENwm6eABnnOKpGUD3ace+fNZXo3++Wy4AA9sjIgZETEdeBb4p36c43igy38gSQ1z4zpt5fIjYFFETI2Iw4CPA3tVN2V1zZ/PMvHnszEDfNZdwL4AkmZI+l06wNWPJBXS9R/IDJD1PUlTgH8ALkpLWsdJuk7S5yXdDnxW0qsk/Z+k36d/D6naK6ysWcBLEfH14oqIWAr8XslgYvdJ+qOk02HXUo2kuZLmpctHKxlg7C4yQS095jfpue7rT8m0jvnzOTD+fEZEQz2A59K/Q4EfAKekz1uBN6TLnwSuTJefAEamy+PSv/OAuZlzXgfcCgxNn48FhqXLJwA/TJePB26tdh6UMS8/AHyhi/XDgLHp8gTgIZImj1OAZZn95gLzusj/K4r7AS8DRqXLBwGL0+WdzpWXhz+f/nyW89EwP9cyRktaSvIGLAH+V1ITyT/HHek+15P8c0Hyxn5X0k3ATT2c9weR9LwEaCJpy30QyUBZw8v5AuqAgP+Q9HpgG0kptNufxV3k/7dJOj9BkndfkTQD2AoMqP60DvjzWXkN8/lsxCqajRExA9gfGEHvdZxvAa4CjgaW9FCH+Xxm+VPA7ZHUo74VGDWgFNeu+0nypbOzgT2Bo9O8fookD7aw82eumC+i+xEjL0qPP4Jk/JoRA051bfPns3wa/vPZiAEeSLohk/yEmwu8ALRLOi7dfA5wh5LBhvaLiNuBjwLjSHpibiAZp6I7TcCqdPm8sie+dvwKGCnp74srJL2SJDitiYiXJM1Kn0Pyj/BySXso6eD01wARsQ7okHRsut/ZmWs0AU9GxDaS92VoJV9QrfDnsywa/vPZsAEeICJ+T9IV/t3AucAVklqBGST1nEOB70j6I8lkFV9I3+wfA28v3sTq4tT/CVwu6bfU2BteTpFUNr4dOFFJM7T7Sep/bwNmSlpM8s+wIt3/JZJ8vZukTjg7xsz5wFXpTayNmfVfBc6V9DuSn7/Zkmiu+fM5MP58uiermVluNXQJ3swszxzgzcxyygHezCynGjLASxot6Q5J/b7BJOlQSXdJ2ixpbqdtp0h6UNJDkv41s36BpDcOJO21rkx526RkHJY/SLpf0vnp+hGSft1DU8BcKVNeFtKer62S7pE0PV3fUHkJjZmfDRnggQuAlkzHj/54lqQZ24LsyvTDcxVJR4jDgLMkHZZu/jLwr+RbOfL2n4AHIuIIkt6Vn5M0IpIp/X5J48zRWo68/DiwNCKagb8DvgjQgHkJDZifjRrgzwZuVjJ63iJJN0paIem7UjLVUG8iYk1E3Au81GnTq4CHIuLh9E3/HnB6esxjwB6SJpbxtdSaAectSaeSMen+u5N8mW5Jt93Ezu2Q86wceXkYSeAhIlYAUyQVe23eROPkJTRgftbUz4nBIGkE8IqIeFTJwExHAoeTjOnxW+B1wJ2SvkAyWFFn34uIz/RwiX2BxzPPVwLHZJ7fl17jh/1+ETWqjHn7FeCW9LgxwJlpRxKAZcArK/pCakAZ8/IPwOx031eRdOqZRNKppyHyEho3PxsuwJMMLrQu8/yeiFgJoB1jgNwZERf18/xdlQSynQ3WAPv089y1rlx5ezKwFHgjMJVkPJbfRMT6iNgq6UVJYyJiQ5nTX0vKlZefAb6YHlPsELUFoIHyEho0PxsxwG9k57E3NmeWt5LmyQBK8CuB/TLPJ5GUEopGsXNPuDwpV96eD3wm7Yn4kKRHgEOBe9L9RgKbypz2WlOWvIyI9ST5WRwf/ZH0UdQIeQkNmp8NF+Ajol3SUEk9DrA0gBL8vcBBkg4gGe/j3cDfZLYfzI6RAHOljHnbBrwJ+E1av3kI8DCApD2AtWm38twqV14qmRHqhfR+0HuAX6dBqmHyEho3PxsuwKd+DhzLjht3fZbeKF1MMrb2NkkfAg6LiPWS/hn4Gck4H9dGxP3pMcOBA9Pj8mrAeUsy2uF1SsZYEfCxiHg63TaLZCyRRlCOvJwGfEvSVuABYE5mWyPlJTRgfjbkWDSSjgT+JSLOGeTrvh04KiIuGczrDqZK562kFuDfIuLBSpy/ljgvy6sR87Mhm0mmo/TdPpAOD/00DPjcIF9zUFUyb9OWEDfV0j9QJTkvy6sR87MhS/BmZo2gIUvwZmaNwAHezCynHODNzHLKAd6qRtLFSkaLbFUyvdwx6foPSXpZZr/b0vbH5bjmc71sHyfpHzPP95F0Y5muvUjSzHKcy6wUDvBWFZJeQzKp8VHpyHwnsGMMnw8B2wN8RJyazjU6GMYB2wN8RDwREWcM0rXLogqtw6xGOcBbtewNPB0RmwEi4umIeELSB0jG6rld0u0Akh6VNEHSlHT0v2skLUtHATxB0m8l/Tkd/AlJ85QZoz/dd0r24pJ2l/RLSfdJ+qOk09NNnwGmpr8orkivuSw9ZpSkb6b7/17SrHT9eZJaJP00Tcd/lpoJ6fl/k6bjPkmvTdd/O5Mm0td6Wtob8wpJ96a/fN6bbj9e0u2S/ptkjBQziAg//Bj0B8kwwEuBP5HMTP+GzLZHgQmdn5MMCLUF+CuSwskS4FqS3q6nk7RDBpgHzM0cvwyYki4/l/4dBoxNlycAD6XnmQIsyxy7/TnwYeCb6fKhJEMqjALOIxlKoSl9/hiwXxeveREws9O6lwGj0uWDgMXp8hsyr6eJZLyTYcCFwL+n60eS9Io+gGTc/OeBA6r93vpROw+X4K0qIuI54GiSgLUW+L6k80o49JGI+GMkwwffD/wyIoKk1DqlD0kQ8B+SWoFfkAzzvFfPh3As8O00/StIAvnB6bZfRkRHRGwi6cK+f4npGA78Vzosww9IxhsnIu4ADpT0cuAs4IcRsQU4Cfi7dDTDu4E9SL4YIBkh8RHMUo06Fo3VgEhm1lkELEoD3LnAdb0clh0FcFvm+TZ2fJ63sHP1Y1cDTJ0N7AkcHREvSXq0m/2yepoUosvRCUtwEclY4keQpDk7EuG303S+m2Q2omIa3h8RP9spYdLxJCV4s+1cgreqkHSIpIMyq2aQlIgBNpBM9NFfjwJHpdc5iqQKo7MmYE0a3Gexo8Td07V/TTpjj6SDgcnAQLumNwFPpr9IziEZoK7oOpIbzkQ6YB3JIHbvSweuQ9LBknYbYBosp1yCt2rZHfhy2vxxC0kd+IXptquBn0h6MiK6Gpu7Nz9kRzXGvST1/J19F/ixpMUk9wJWAETEM+lN22XAT0jm1y36KvD19NfGFuC8iNiskmd7A+D/SSoOJ3sXyRyfP5T0TuB2MqXwiHhK0nKSqeCKriGpirpPyYXXAm/rSwKscXgsGrMalfYF+CNJU9KOaqfH6o+raMxqkKQTSH5VfNnB3frLJXgzs5xyCd7MLKcc4M3McsoB3swspxzgzcxyygHezCynHODNzHLq/wMJyAaoPFvC5AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "dists = {\"ant\": {'L2/3': [],\n",
    "                 'L4': [],\n",
    "                 'L5': []},\n",
    "         \"post\": {'L2/3': [],\n",
    "                 'L4': [],\n",
    "                 'L5': []}\n",
    "         }\n",
    "plt.clf()\n",
    "i = 1\n",
    "xtick_labels = []\n",
    "df_datasets_ant = []\n",
    "df_datasets_post = []\n",
    "for stim_loc in ['L2/3', 'L4', 'L5']:\n",
    "    df2 = df[df['Stim Location'].str.contains(stim_loc)]\n",
    "    \n",
    "    df2_ant = df2[df2['Slice Location'] == \"Anterior\"]\n",
    "    dfd1_ant = df2_ant['abs(delta 1)'].dropna()\n",
    "    dfd2_ant = df2_ant['abs(delta 2)'].dropna()\n",
    "    \n",
    "    df2_post = df2[df2['Slice Location'] == \"Posterior\"]\n",
    "    dfd1_post = df2_post['abs(delta 1)'].dropna()\n",
    "    dfd2_post = df2_post['abs(delta 2)'].dropna()\n",
    "    \n",
    "    dfds_ant = np.concatenate([dfd1_ant, dfd2_ant])\n",
    "    dfds_post = np.concatenate([dfd1_post, dfd2_post])\n",
    "    \n",
    "    if len(dfds_ant) > 3:\n",
    "        dists['ant'][stim_loc] = dfds_ant\n",
    "        print(\"Latency between anterior (rostral) barrel pairs in\", stim_loc)\n",
    "        print(\"\\t\", \n",
    "              str(np.mean(dfds_ant))[:5], \n",
    "              \"+/-\", \n",
    "              str(np.std(dfds_ant))[:5], \n",
    "              \"(n =\", \n",
    "              str(len(dfds_ant)) + \")\")\n",
    "        print(\"\\t\", dfds_ant)\n",
    "    if len(dfds_post) > 3:\n",
    "        dists['post'][stim_loc] = dfds_post\n",
    "        print(\"Latency between posterior (caudal) barrel pairs in\", stim_loc)\n",
    "        print(\"\\t\", \n",
    "              str(np.mean(dfds_post))[:5], \n",
    "              \"+/-\", \n",
    "              str(np.std(dfds_post))[:5], \n",
    "              \"(n =\", \n",
    "              str(len(dfds_post)) + \")\")\n",
    "        print(\"\\t\", dfds_post)\n",
    "    if len(dfds_post) < 3 and len(dfds_ant) < 3:\n",
    "        continue\n",
    "    \n",
    "    jittered_xs_1 = np.random.normal(i, 0.04, dfds_ant.shape[0])\n",
    "    jittered_xs_2 = np.random.normal(i+2, 0.04, dfds_post.shape[0])\n",
    "    plt.scatter(jittered_xs_1, dfds_ant, alpha=0.5, color='blue')\n",
    "    plt.scatter(jittered_xs_2, dfds_post, alpha=0.5, color='green')\n",
    "\n",
    "    xtick_labels.append(stim_loc + \"\\nRostral\\n(n=\" + str(len(dfds_ant)) +\")\")\n",
    "    xtick_labels.append(stim_loc + \"\\nCaudal\\n(n=\" + str(len(dfds_post)) +\")\")\n",
    "    df_datasets_ant.append(dfds_ant)\n",
    "    df_datasets_post.append(dfds_post)\n",
    "    i += 1\n",
    "\n",
    "# calculate significance\n",
    "for sl1 in ['ant', 'post']:\n",
    "    for sl2 in ['ant', 'post']:\n",
    "        for la1 in ['L2/3', 'L4', 'L5']:\n",
    "            for la2 in ['L2/3', 'L4', 'L5']:\n",
    "                if sl1 == sl2 and la1 == la2:\n",
    "                    continue\n",
    "                elif len(dists[sl1][la1]) < 3 or len(dists[sl2][la2]) < 3:\n",
    "                    continue\n",
    "                else:\n",
    "                    p_val = kstest(dists[sl1][la1], dists[sl2][la2]).pvalue\n",
    "                    if p_val < 0.05:\n",
    "                        print('KS test for ' + sl1 + \" layer \" + la1 + \" and\", \n",
    "                              sl2 + \" layer \" + la2,\n",
    "                             \"is significant with p-value \", p_val)\n",
    "\n",
    "# significance: L2/3 Rostral and L4 Rostral, and \n",
    "#               L2/3 Rostral and L4 Caudal\n",
    "significance_bar(0.75, 2.25, 3.35, \"P<0.05\", fontsize=10)\n",
    "# significance_bar(0.75, 4.25, 3.75, \"P<0.05\", fontsize=10) -- not signficant when correcting for multiple tests\n",
    "\n",
    "df_datasets_post.reverse()\n",
    "plt.boxplot(df_datasets_ant + df_datasets_post)\n",
    "plt.xticks(range(1, len(xtick_labels)+1), xtick_labels)\n",
    "plt.title(\"Latency of L4 Barrel Pair Crossings\")\n",
    "plt.xlabel(\"Stimulation Layer\")\n",
    "plt.ylabel(\"Latency (ms)\")\n",
    "plt.savefig(data_dir + cell_type_prefix + \"_results_latency_prelim_\" + str(datetime.date.today()) +\".png\", bbox_inches=\"tight\")\n",
    "plt.show()\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}