Day 5 - Univariate Statistics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Univariate Statistics and the SME\n",
    "\n",
    "Yesterday, we learned how to extract measures of spectral power from EEG signals. Today, we're going to extend that work by analyzing **how spectral power correlates with cognitive, behavioral, or physiological variables** of interest. This is the core of almost any analysis you'll carry out in the Computational Memory Lab. \n",
    "\n",
    "The statistical approach we'll use today will be mostly straightforward T-tests to assess differences between normal distributions of power values. In particular, we'll focus on the contrast between power during successfully-encoded versus forgotten words. You'll see this called the **Subsequent Memory Effect** (SME) in the literature. In general, you can take an SME-style approach to any physiological variable.  \n",
    "\n",
    "The basic analysis scheme is as follows: \n",
    "\n",
    "1. Extract spectral power for all encoding events (e.g. wavelet transform)\n",
    "2. Normalize spectral power for each experimental session\n",
    "3. Use a 2-sample T-test to compare distributions of power between remembered and not-remembered items\n",
    "4. *Optional:* Average results into time windows or frequency bands of interest. \n",
    "5. Repeat for each subject in the dataset.\n",
    "\n",
    "Let's get to work!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from CMLLoad import CMLLoad\n",
    "\n",
    "# Remember to update this to where your CMLExamples data is located.\n",
    "load = CMLLoad('./CMLExamples')\n",
    "\n",
    "# We load the dataframe for all sessions.\n",
    "df = load.Index()\n",
    "\n",
    "# Specify which subject and experiment we want\n",
    "sub = 'R1111M'\n",
    "exp = 'FR1'\n",
    "\n",
    "# Select the dataframes for our subject and experiment.\n",
    "df_select = df[(df['subject']==sub) & (df['experiment']==exp)]\n",
    "# Pick the first session\n",
    "df_sess = df_select.iloc[0]\n",
    "\n",
    "# For first session...\n",
    "evs = load.Load(df_sess, 'events')\n",
    "word_evs = evs[evs['type']=='WORD']\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPP total time wavelet loop:  12.540498971939087\n",
      "0 completed successfully!\n",
      "CPP total time wavelet loop:  9.192278861999512\n",
      "1 completed successfully!\n",
      "CPP total time wavelet loop:  13.955910444259644\n",
      "2 completed successfully!\n",
      "CPP total time wavelet loop:  11.25844693183899\n",
      "3 completed successfully!\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "# Instead of just looking at session 1,\n",
    "# we want to analyze data from all sessions.\n",
    "# To do that, we'll use a for loop.\n",
    "\n",
    "recalled = []  #for storing encoding success\n",
    "sess_pows = []  #for storing our computed spectral powers\n",
    "\n",
    "# Iterate over the DataFrame rows for each session\n",
    "for df_sess in df_select.itertuples():\n",
    "    # Load this session's events, and select word events.\n",
    "    evs = load.Load(df_sess, 'events')\n",
    "    word_evs = evs[evs['type']=='WORD']\n",
    "\n",
    "    # Grab the events as a PTSA object.\n",
    "    buf = 1000  #to remove edge effects during wavelet convolution\n",
    "    eeg_ptsa = load.LoadPTSA(df_sess, -buf, buf+1600+buf)\n",
    "    # Select the word events.\n",
    "    eeg_ptsa = eeg_ptsa[word_evs.index]\n",
    "    \n",
    "    from ptsa.data.filters import ButterworthFilter\n",
    "    from ptsa.data.filters import MorletWaveletFilter\n",
    "    \n",
    "    sess = df_sess._asdict()['session']\n",
    "    sr = eeg_ptsa.samplerate\n",
    "    \n",
    "    # Filter out line noise.\n",
    "    freq_range = [58., 62.]\n",
    "    b_filter = ButterworthFilter(timeseries=eeg_ptsa, freq_range=freq_range, filt_type='stop', order=4)\n",
    "    eeg_filtered = b_filter.filter()\n",
    "    \n",
    "    # Get spectral power.\n",
    "    myfreqs = np.logspace(np.log10(3), np.log10(170), 10)  #let's start with 3-170 Hz, 10 frequencies\n",
    "    wf = MorletWaveletFilter(timeseries=eeg_filtered, freqs=myfreqs, width=5, output='power')\n",
    "    # Output is freqs, events, elecs, and time.\n",
    "    pow_wavelet = wf.filter()\n",
    "    pows = pow_wavelet.data[:, :, :, int((buf/1000.)*sr):-1*int((buf/1000.)*sr)]  #remove buffer period\n",
    "    pows = np.log10(pows)\n",
    "    \n",
    "    # Reshape array to average into 100 ms bins.\n",
    "    samps = int(sr/10.)  # samples per window\n",
    "    numwins = int(pows.shape[3]/samps)\n",
    "    reshaped_arr = np.reshape(pows[:, :, :, :samps*numwins], (pows.shape[0], pows.shape[1], pows.shape[2], numwins, samps))\n",
    "    avg_pows = np.mean(reshaped_arr, axis=4)\n",
    "\n",
    "    from scipy.stats import zscore\n",
    "    # z-score power values across events (important to do within-session).\n",
    "    z_pows = zscore(avg_pows, axis=1)\n",
    "    sess_pows.append(z_pows)\n",
    "\n",
    "    # We need to store this data across-sessions.\n",
    "    recalled.extend(list(word_evs['recalled']))\n",
    "\n",
    "    print(str(sess)+' completed successfully!')    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's unpack what just happened up there: \n",
    "\n",
    "```\n",
    "# Load this session's events, and select word events.\n",
    "evs = load.Load(df_sess, 'events')\n",
    "word_evs = evs[evs['type']=='WORD']\n",
    "```\n",
    "\n",
    "You should be pretty used to this syntax by now; this is just how we load the events, and filter the events dataframe for only the **word encoding events**.\n",
    "\n",
    "```\n",
    "# Grab the events as a PTSA object.\n",
    "buf = 1000  #to remove edge effects during wavelet convolution\n",
    "eeg_ptsa = load.LoadPTSA(df_sess, -buf, buf+1600+buf)\n",
    "# Select the word events.\n",
    "eeg_ptsa = eeg_ptsa[word_evs.index]\n",
    "```\n",
    "\n",
    "Also the same as last time! We're loading the 1600 ms that correspond to the word presentation interval, and we're adding **1 second buffers** to account for edge effects of our wavelet filter.  Then we are reassigning the eeg_ptsa data to a selection of just the word encoding events.\n",
    "\n",
    "```\n",
    "# Filter out line noise.\n",
    "freq_range = [58., 62.]\n",
    "b_filter = ButterworthFilter(timeseries=eeg_ptsa, freq_range=freq_range, filt_type='stop', order=4)\n",
    "eeg_filtered = b_filter.filter()\n",
    "```\n",
    "\n",
    "This is our standard stop-band filter to remove the effects of 60 Hz line noise. We could also have applied one to the 120 Hz harmonic, but typically that's less of an issue.\n",
    "\n",
    "```\n",
    "# Get spectral power.\n",
    "myfreqs = np.logspace(np.log10(3), np.log10(170), 20)  #let's start with 3-170 Hz, 10 frequencies\n",
    "wf = MorletWaveletFilter(timeseries=eeg_filtered, freqs=myfreqs, width=5, output='power')\n",
    "# Output is freqs, events, elecs, and time.\n",
    "pow_wavelet = wf.filter()\n",
    "pows = pow_wavelet.data[:, :, :, int((buf/1000.)*sr):-1*int((buf/1000.)*sr)]  # Remove buffer period.\n",
    "pows = np.log10(pows)\n",
    "```\n",
    "\n",
    "The meat of the analysis! We defined our desired frequencies (note that they're **log spaced** since a wavelet is less sensitive to differences at higher frequencies), created a wavelet filter object, ran the convolution, removed the buffer period, and took the log transform. \n",
    "\n",
    "```\n",
    "# Reshape array to average into 100 ms bins.\n",
    "samps = int(sr/10.)  # samples per window\n",
    "numwins = int(pows.shape[3]/samps)\n",
    "reshaped_arr = np.reshape(pows[:, :, :, :samps*numwins], (pows.shape[0], pows.shape[1], pows.shape[2], numwins, samps))\n",
    "avg_pows = np.mean(reshaped_arr, axis=4)\n",
    "```\n",
    "\n",
    "This is new, and technically optional. Here, we're averaging spectral power into 16 100 ms bins that span the entire 1.6 second interval. It's effectively a **data-smoothing technique**. You'll find spectral analyses to often involve some type of temporal windowing, which is important because the precise timing of effects differs somewhat between electrodes, brain regions, and subjects. Averaging into reasonable windows of interest removes some of that variability and let's the true signal come through. There's no \"right\" choice for your windows, but I find that nonoverlapping 100-200 ms works well for these types of data. (And it's really only important for the high frequencies, since the wavelet already smoothed our low-frequency signals). \n",
    "\n",
    "```\n",
    "from scipy.stats import zscore\n",
    "# z-score power values across events (important to do within-session).\n",
    "z_pows = zscore(avg_pows, axis=1)\n",
    "sess_pows.append(z_pows)\n",
    "```\n",
    "\n",
    "Finally, our last step is to **z-score, or normalize, the power values** across events but within each experimental session. Why? Properties of the recording hardware can change from session-to-session, yielding power values that may fall in entirely different ranges from one day to the next. If we averaged these raw values together, we'd get an innaccurate read of the underlying neural signals. We can account for this by z-scoring the power across all the events in each session. Be sure you're z-scoring across events and not some other dimension!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1020,)\n",
      "(10, 1020, 141, 16)\n"
     ]
    }
   ],
   "source": [
    "# Now let's compute some statistics!\n",
    "from scipy.stats import ttest_ind  # this is our 2-sample t-test\n",
    "\n",
    "recalled = np.array(recalled)\n",
    "all_pows = np.concatenate(sess_pows, axis=1)\n",
    "\n",
    "print(recalled.shape)\n",
    "print(all_pows.shape)  # dimensions are freqs, events, elecs, and time bins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compare remembered vs. not-remembered.\n",
    "# Final output in shape [freqs, elecs, time].\n",
    "# equal_var=False applies Welch's t-test,\n",
    "# which assumes unequal population variances\n",
    "tstats, pvals = ttest_ind(all_pows[:, recalled==1, :, :], all_pows[:, recalled==0, :, :], equal_var=False, axis=1)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### T-tests\n",
    "\n",
    "Brief aside on t-tests and how to use them.\n",
    "\n",
    "For more information, see: https://blog.minitab.com/blog/adventures-in-statistics-2/understanding-t-tests-t-values-and-t-distributions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "middletemporal\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVsAAAC5CAYAAACGNJzQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3deZwcVbn/8c93ZhKyQgTCGi5BQVwAWXIRN1xABcTghoLIIiou6AUEidxcULxXFFwRRA27GlFkUfCyRRa5/hAlhE1kkVXDFgKCCIRkZr6/P6o6dHqqZ6prunq6Z57361Wvma6uU+fMZPLMmVPnPEe2CSGEUK6ukW5ACCGMBRFsQwihBSLYhhBCC0SwDSGEFohgG0IILRDBNoQQWiCCbUrSLpLuknSPpC+WXNdGkq6WdIek2yUdUmZ9VfV2S7pJ0m9aUNc0SedJujP9Ol9Xcn2Hpd/LP0s6R9KEJt//DElLJP256tyakhZI+mv68SUl1/eN9Pt5q6QLJU0rs76q946QZElrN6u+sSiCLUkQAr4P7Aq8Cthb0qtKrLIXONz2K4EdgINLrq/iEOCOFtQDcCJwme1XAK8ps15JGwL/AcyyvQXQDezV5GrOAnapOfdF4ErbmwFXpq/LrG8BsIXtrYC7gaNKrg9JGwFvB/7WxLrGpAi2ie2Be2zfZ3s58HNgj7Iqs/2I7UXp58+QBKINy6oPQNIM4F3AaWXWk9a1OrAjcDqA7eW2nyq52h5goqQeYBLwcDNvbvta4Mma03sAZ6efnw28p8z6bF9huzd9eT0wo8z6Ut8BjgRi9dMwRbBNbAj8ver1YkoOfhWSZgLbAH8suarvkvyn6S+5HoCXAo8DZ6bDFqdJmlxWZbYfAr5J0vt6BHja9hVl1VdlXduPpG14BFinBXVWHAhcWmYFkmYDD9m+pcx6xooItgllnCv9N7mkKcD5wKG2/1liPbsDS2zfWFYdNXqAbYEf2N4GeJbm/om9inSsdA9gE2ADYLKkj5RV30iTNJdkKGp+iXVMAuYCx5RVx1gTwTaxGNio6vUMmvxnaC1J40gC7XzbF5RZF/AGYLakB0iGSN4m6acl1rcYWGy70ls/jyT4lmVn4H7bj9teAVwAvL7E+ioek7Q+QPpxSdkVStof2B3Yx+UmNnkZyS+vW9KfmxnAIknrlVjnqBbBNnEDsJmkTSSNJ3m4clFZlUkSyXjmHba/XVY9FbaPsj3D9kySr+0q26X1/Gw/Cvxd0ubpqZ2Av5RVH8nwwQ6SJqXf251ozYPAi4D908/3B35dZmWSdgHmALNtP1dmXbZvs72O7Znpz81iYNv03zYUEMEWSB86fBa4nOQ/6bm2by+xyjcA+5L0MG9Oj91KrG8kfA6YL+lWYGvguLIqSnvQ5wGLgNtIfq7nNbMOSecAfwA2l7RY0seArwNvl/RXkif2Xy+5vpOBqcCC9GfmhyXXF5pIkWIxhBDKFz3bEEJogQi2IYTQAhFsQwihBSLYhhBCC0SwDSGEFohgW0PSQVFf1NdudY2F+ka7CLYDtfoHLOrr3PpG89c2EvWNahFsQwihBcbkooYJ6vZUejLfW0YfE+jOfG/NSeMK1Tfl5ZvWfe/xJ55g+lprZb736C13Fqpv+SD/poN9fdOnFcu3PXF6/ZzZjz/9DNPXmJr5nvv6CtX39/vqpyAY7Otba3Kxf78nnl3RcF3Dsd66UzLPP/ncMtacVP/fyH3FEro9tjR75W8ZX98z9LLMfVmJnwa1kSZ6WZ2EdUtZfrntAbl42012xBnlptLD+1m/4XJ7v7rxMgA7LCi2McLx6722ULmHl/UOfVGGT75186EvyvDKT+1ZqFz/M8VS3H7+AycXKrfPlsWyZs6//qFC5Yqas88bC5Vb/kyxdAnfOXVRoXJFnM8jhcotp5+9uzbIfO+k/gc6YgeJMRlsQwidp7vh/nB7admYbZ09lX5RlYjlAUk3V713VLof2F2S3lnnnqXtARVCaB8SjO9S5tEpWvmA7Cxq9jiy/SHbW9vemiS36wUA6X5cewGvTsucku4TVqvMPaBCCG2iiwi2uQ2yx1Elv+sHgXPSU3sAP7f9gu37gXtI9gmrVdoeUCGE9iHEOGUfnaJdxmzfBDxm+6/p6w1JNrSrqLcn2Cp7QElq5R5QIYQWqQwjdLJ2CbZ782KvFkrYEyxdDXMQwJQSpuuEEMojItgOW7r19PuA7apO590T7DFJ66e92kH3gLI9jzR7/3StNvYmF4fQwUZDsG2HFWQ7A3faXlx17iJgL0mrSdoE2Az4U0bZlu4BFUIYGTEboQGD7HG0F6sOIZDu/3UuySaBlwEH2+5L73OapFnppaXtARVCaB9dZAfaPMFW0kaSrpZ0h6TbJR3SgiYP0LJhBNt71zl/QJ3zXwW+mnH+41WfP0Gyk2oIYZQbxqKGXuBw24skTQVulLTAdpk7Pg8w4mO2IYQwlK5hzEZIZyxVZi09I+kOktlNEWxDCKHaEA/I1pa0sOr1vPSB+MD7SDOBbYA/NrN9eUSwbcA5NxRLorFV98RC5f7zzvMLlbvnK18qVG6d7/6sULkVFxxfqNxVc84rVO64Uz5UqNzkzV9dqNz8nY4pVO4l44o9Ernu7BsLlfu/J54vVK6oz39q1tAX1fj9+VcWqkuCcV11v59LbQ/ZGElTSFaqHmr7n4UaMgwRbEMI7U+ie3zx5/mSxpEE2vm2L2hauxoQwTaE0PYk6BpfbDFSmg7gdOAO299uasMaMGSwTRv6JuDNwExgIvA4sAhYYDtrsUHWfc4AdgeW2N4iPbcm8Iv0vg8AH7T9D0njgR8Bs4B+4BDb12TcM7N8nvaEEDqIRPe4wis/3wDsC9xWlVnwP21f0pS25VS3Xy5pgqQ5wIPAAmA2sHZa5pXAccCDki6qmvc6mLOoyfpF/axdnwCwvSXJ/NlvScpqa2T9CmEMkKB7fFfmMRTbv7ct21tVsgy2OtDC4Isa7gZ2AD4HrG77323vYXsv2++wvSFJ0L0BuLBqkUKmOlm/6mXtehVJ8MT2EuApkl5urcj6FcJYIKHursyjUwzW0nfbfq/tX9t+IesC2/fY/m9gU+D3BepfJWsXUMnadQuwh6SedLnudqyaK2Go8gNIOkjSQkkLl1Fs76sQwsiQoHtcV+bRKeqO2dq+Je9N0mB8V1NalDiDpNe8kGQY4zqSVSCFRSKaEDqYRHfBB2TtItevBUl3S/pROn2i+vxaku4eRv2Ppdm6qM7aZbvX9mHp2MoewDTgr3nLhxBGl2TMtjvz6BR5++CbkjzculJS9b7bPcDLhlF/ZtYuSZMkTU4/fzvQW2cdc2T9CmEMkNTxwwh5W2qSWQHPAjeke4Q1pE7Wr3pZu9YBFqVrmOeQTNuo3CeyfoUw1oyCnm3eRQ0imRHwLuBbwHWS9gJuyltRvaxfZGTtsv0AsHmd+0TWrxDGGHWp8KKGdpE32BrAdj9wWLod+fnA98pqWAghrJTORuhkjfRsV7J9evpgrFimlDGm5/xiiVr+taLYBIyZJ80vVO6ep1cUKjd+18MLlVv3zcVyOP9r6vhC5f6w9esLlXtlwfo+ed6Rhcod+s7/KVSu1TY+4dSGy6z2x7cVqkujYDZC3mC7GbC0+oTt/5P0GpIFCCGEUB7RUQsYsuQKtrbvrXN+ZVLeEEIoSzIbobPzZg3aeklX5LmJ7XcMpxGSpgGnAVuQjA8fCDwH/BCYQpJkZp+sHJSSdgFOBLqB02zHjIQQRhuJrvGjONgCtT3ajwG/BJqdePdE4DLbH0gzfk0iSX5zhO3fSToQ+AJwdHUhSd3A90mmfS0mmZZ2Uav3FgohlEyia9y4oa9rY4MGW9ufrn4t6SPA0bbva1YDJK0O7AgckNa5HFguaXPg2vSyBcDl1ARbYHvgnkp7JP2cJDlNBNsQRpHRMIzQDiPOLyXJj3umpJvSRQuTgT+TpHUE2JPsRDQbAn+ver04PRdCGE3SYYSso1O0Q7DtAbYFfmB7G5JVal8kGbc9WNKNwFRgeUbZrB3gMpPMRNavEDqXJLrG9WQenaIdgu1iYLHtym6X5wHb2r4zzZu7HXAOA8ePK2Wre7wzgMydI2zPsz3L9qwJdPZ8vRDGHEH3+J7Mo1MMNRthds2pLmAnSVtUn7R9UdEG2H5U0t8lbW77LpLlt3+RtI7tJekODf9FMjOh1g3AZmnO24eAvYAPF21LCKE9VXq2nWyo1v8q49yPal4bht1V/BwwP52JcB/wUWA/SQen718AnAkgaQOSKV672e6V9FmSh2fdwBm2bx9mW0II7SbdqaFY0YH7H46EoYJtS+Za2L6ZgdvenJgetdc+DOxW9foSoOX7CYUQWieZjVA4HJ0FnAz8uGkNKmCoqV/xJCmEMPJE4ZkHtq+VNLOp7SlgsN11d8h7E0mTJb26OU0KIYRVjYbZCIO19KeS7iVZRnuJ7WdrL0iTiH+EZIx1DtAR46XrrTOFIz/ceAao80+9oVB9PWtOL1Tu/vnFRkc2GXdCoXL/tXzXQuXO3nurQuVe9sTNhcpdM7tYlrGLFxdb+Lj/WzYuVK7V2btmTir2Z/Y73r5JoXLPslrDZfoyZ2vmMPgDsrUlLax6PS/dc7CtDBZsXwl8BjgO+Jmku0iSziwDXgK8gmT+66+Ad9i+reS2hhDGKgn11E11udR27TOftjPY7rorSB9SpUMKbwQ2BiaS9GBPAq6y/XgrGhpCGMMkNJpzI1TYvh64vsyGpEllFgIP2d696vwRwDeA6baXZpSLrF8hjHYS1O/ZDlFU5wBvIRluWAx8yfbpTWxdLu00unwIcAeweuWEpI1IMnr9LatAZP0KYWwQQl3F5tkOsv9hS7XDcl0kzSDZTPK0mre+AxxJnXwHVGX9SrOFVbJ+hRBGk0rPNuvoEO3Ss/0uSVCdWjmRLhV+yPYtUt0nmFlZv16bdaGkg4CDADaYOqkJTQ4htIyEesbAmG2ZJFWW0d0o6S3puUnAXGCoHSByZ/1Kp4LMA9hy3bXq9ZRDCO1IQuM6pxebJVewlbSm7SdLasMbgNmSdgMmkIzZ/gTYBKj0amcAiyRtb/vRqrK5s36FEDqYBB3es807ZvuwpJ9LenuzG2D7KNszbM8kydp1le33217H9sz0/GKStIuP1hRfmfUrTWKzF1A4A1kIoV0lwwhZR6fIG2zfm157saQHJX1JUrFlNcMkaQNJlwDY7gUqWb/uAM6NrF8hjELqQj3jM49OkXee7aXApZLWBPYl2S/saElXAacDF6azAYbF9jXANRnnZ1Z9Hlm/QhhjNAoWNTQ09cv2k7ZPTLevOZRko8ZzSIYZ/id9sBVCCM011qZ+SVoH2I8k8cxMki1sTgc2INk3bHuGnkEw4nrX3Yilh53ccLn9Njy2UH2/+dgphcqtvsaEQuVO+cXZhcod/bdDC5WbM/VVhcoVNXvjNVpaX++y3kLldl13cqFylz42IOdTLg88t6JQuXm/vrtQuTlfPqjhMn0PP1ioLqDwooZ2kXc2wmySDRh3Be4i2a3hJ7b/UXXNzcBNZTQyhDDGSdA98sMIkg6xfeJQ57Lk/VUxH1gK7Gh7K9vfqw60qfuA43PeL4QQGiDc1ZN5tNj+GecOyFMwb0vXt/2vwS6w/TxwdM77rZTmP/gxsB7QT5KL8kRJW5Ns8jgB6AU+Y/tPGeUjEU0Io50EXSO3K7akvUk2k91EUvX00qnAE3nukTfYvkPSCtsX1zTg3UCP7Qtz3idLL3C47UWSpgI3SloAnAAca/vSdMHDCSSZe6rrj0Q0IYwBBtw9ogteryPJ57028K2q888At+a5Qd7WfwU4IuP8MpJkMYWDre1HSL4IbD8j6Q6SnAfmxQxga5C9MmxlIhoASZVENBFsQxhNJGj9kMFKth8EHgReV/QeeVv/MuDOjPN/BV5atPJa6aZs2wB/JJladrmkb5KMLWftY5M7EU0IoZNpJMZnB7Yi2UjhJJKdbMaTDF8+a3v1QQuS/wHZUyQBt9ZmJN3oYZM0BTgfONT2P4FPA4fZ3gg4jGSK2YBiGecyk8xIOkjSQkkLn3oy1xBLCKFdVMZss47WOhnYm6SjORH4OEnwHVLeYHsR8B1JKwOupE1Jxi6GnYtA0jiSQDvf9gXp6f2Byue/JBkyqJU7EY3tebZn2Z41bc21htvkEEKrqSv7aDHb9wDdtvtsnwm8NU+5vC09EngOuFPS/ZLuJ8lF8DzwhSINrlCS1ut04A7b365662HgzennbyP5TVIrEtGEMCYId/dkHi32XBprbpZ0gqTDgFyrV/LmRnha0uuAXYCtSf58XwRcbnu4uWHfQJJv4bZ0YQTAfwKfINlssofkQdxBkCSiIZnitZvtXkmVRDTdwBmRiCaEUWiYD8iaOEV0X5JO6mdJhjc3At6Xp2Du1qdB9dL0aBrbvyd77BVgu4zrIxFNCGNQ0QdkTZ4i+p50tdgy4Nj0/oeQBPJB5W69pFkkf86vQ83wg+3PN9LaEEJoyPAWNTRziuj+DAysB2ScGyBvboTDSB6GPUAyllo9dBBbzIQQSjasqV/DniI6yAqy1WnyCrLDgM/b/m4jDWxbD9xL/4G5hllWsXSDKYWq226frQuV09wfFSr36U/vWajcDVsVm699+GeyJooM7c4Li609uWlJsaxYRc2//qGW1lfUFquvVqjcn//5QqFyy595ruEy7usvVJclXL9nu7akhVWv56V7DlbkniI6iJatIFuDeMofQhgpNn39dePjUtuzBik97L0KKyvIJO0MPG+7X9LLgVcAt+W5R96pX+fSAXlqQwijk4E+Zx85NHOK6LXABEkbAleS5PY+K0/BvD3be4H/Tpeq3QaskqXY9vdyNzWDpDOAypbmW6TnIutXCGGlorNMmzxFVLafk/Qx4CTbJ0jKlcc7b7A9mGSqw07pUc3AsIItyW+Gk0lSLVZE1q8QAvBiz7Zw+eZNEVW65mAf4GPpuVxxNO+iho2Gvqo429emSWhWOU1k/QohAJjBxmxb6RDgKJJNbm+X9FLg6jwFG55LIWkt4MkmrBwbSmT9CiEAw+/ZNovta0nGbZG0XtrR+488ZXM9IJM0TtJxkp4CHgM2Sc9/TdKnijV7SOVl/VpRbGO8EMLIcNqzzTpGUEPDEnlnIxwNvJ9kjKJ6Ut6NJE/jylBe1q8O338+hLFoGLMRhi3N0TLgdCP3yBts9wE+aft8kn3CKm4DNm+kwgZE1q8QAlAZRnDm0SIDZkIBpzZyg7xjthuQLNWt1d3APeqSdA7JTIO1JS0GvkRk/QohVCm4+KxZBvRibZ/SyA3yBsq/AG9iYMDdE8g1x2wwtveu81Zk/Qoh4Nb2YrNMl1Q34VZNLu5MjWz4eFbaq+wC3idpc2A/4N057xFCCIWN8MyvbmAKDY7TVss7z/bXkvYB5qaVfpWkR/se21cUrTyEEPKojNmOoEdsf2U4N2gkefio+XN90nprsvURH2643OJt9ypU37pTis1+0KXfL1TuHyuKDW5t88YZhcqdfdqNhcrtd+A2hcrtduhhhcp5WbFsYZ/d+qBC5XaaPqlQud1ObPxnE6B7+oaFyt0//8JC5Satt2bDZbrGFXvEY0Z8UUPhHm3FyO8NHEIIQ7BhxciuaqhNU9CwvMnD/8Eg+R9tN/4rLoQQcjJmRf/ITUew/eRw75G3Z3tEzetxwDbAe4CvDbcRdbJ+fYPk4dtykqxjH7X9VEbZyPoVwijXBj3bYcv7gCxrqSxpdvQ3Z73XoLMYmPVrAXBUOpf2eJLkD3Nq6o+sXyGMAQZWtEcimsLyriCr50qSLFvDkiZ3eLLm3BW2e9OX15Msxa21MuuX7eVAJetXCGEUsc2Kvv7Mo1MMN9juSc7NzobpQLK3UM/K+pX5SLY6Ec3jT/+rhCaGEMo0krkRmiHvA7KbWPUBmYD1gOnAZ0toV3Xdc0l2apif9XbGucxvf7oB3DyA7V6+cQf9E4UQkjHbzunFZsn7gOw3Na/7gceBq8vMRSBpf5IHZzvVyZ877I3cQgjtbzSM2eZ9QHZ02Q2plc4ymAO82Xa9PZNXZv0CHiLJ+lVsRngIoW3ZpreEnq2kPYEvA68Etre9cPASxQ13zLYp0qxffwA2l7Q43UztZGAqsEDSzZJ+mF67gaRLINnIjWQY43LgDuDcyPoVwuhT6dlmHcP0Z+B9pLsvlCnvmO0KBlnUUM32+EYbUSfrV+Z0s8j6FcLYU9aYre07AKRhr8YdUt4x28OBY4CLSXqgAK8jWXTwZZLx2xBCKMWYGbMl2Slhru0fVZ2bl+4/tpvt2c1vWnmWdK3OKZMbX+r80UnFUkk8s7zYb+S1X14sUcsW8z5QqNzvniqWOOWYafWG1Ae3dN7xhcodtcW+hcptOqXhP7oA+MTumxYqV9ThHz6tpfWdcEax7+ej7/pCw2X6LyqWtMg2y+v3bNdOF1hVzEtnHwEg6bcks6dqzbX960INKiBv9NiZgUt2AX4LfKt5zQkhhIH6Dct76wbbpbZn1XvT9s7ltKoxeR+QPUEyiFzrvcDS5jUnhBCyteHuug3J27P9MnCapDfz4pjtDsAupHuDlSHdDeIXVadeChxj+7tV14gkEc1uwHPAAbYXldWmEELr2R6sZ1uYpPcCJ5Es0PpfSTfbfmfTKyL/PNszJd0FHAJ8kGTl1l9I5sD+vzIaltZ7F7A1rEw68xBQm+l4V2Cz9Hgt8IP0YwhhlBhiGKEw2xcyMKaUopGdGq4DriuxLUPZCbjX9oM15/cAfpyuMLte0jRJ69t+pPVNDCGUwZTTs22l3IsaJE2XdKik70laKz23g6SNy2veKvYCzsk4nzsZTQihM9nwQm9/5tEpcgVbSdsAdwEfAz4FrJG+tStwXDlNW6X+8cBs4JdZb2ecGzBqXp3169mnhp10PYTQQjYs7+vPPDpF3p7tt4BTbG8JvFB1/jLgjU1v1UC7AotsP5bxXq5kNLbn2Z5le9bkabGLTwidpD99QJZ1dIq8wXY74MyM8w8D6zavOXXtTfYQAsBFwH5K7AA8HeO1IYwuBpb39mUenSLvA7JlwOoZ5zen5KW6kiaRbHvzyapznwKw/UOSvAi7AfeQTP36aJntCSGMAHfWnNoseYPtxcAxkj6YvrakfwO+DlxQSssqFSXpFdeqOffDqs8NHFxmG0III6s/fUDWyfIOIxxOsrZ4CTAR+B0v9iTnltO0EEJIOJ1n28ljtnkXNTwt6fUkf85vSxKkFwGX19lBIYQQmmY0zLPVULFS0jjgGuDAdEVXx+uZNsPT3tj41mlT1p1ZqL6HF11RqNyM7d9VqNzqa04sVO74A7YrVG5FwV33dtx4jaEvyrDaDecXKtf7+EOFyr2wtNhUwTXeViwZ3u+6Xl6o3Jd+9edC5XqXF3vI9PA9WZODBvfoxUfzwtL7Gk4eu/q/vcL/fsSpme9ddciONw6WiKZdDNmztb1C0mYk+46FEELL2aavg+bUZsk7ZvsTkgUNIYTQeoa+3v7Mo1PknY0wHvi4pLcDC4Fnq9+0/flmN6xC0gPAM0Af0Fv750Jk/Qph9DPQX3C4ql3kDbZbA7emn7+q5r1WfAfearte3tzI+hXCaGd3VC82S97ZCG8quyHDEFm/QhjlDLjDFzUMOmYraStJI73duYErJN0oKStRea6sX9WJaLz82dq3QwjtzNDX1595dIqhAulNwNqVF5L+V9L65TZpgDfY3pZkuOBgSTvWvJ8r61d1IhqNn1xGO0MIJXE6jNDJD8iGCra1gWxHkhVkLWP74fTjEpKM6tvXXJIr61cIobP19znz6BQjPUQwKEmTJU2tfA68A6iduR1Zv0IY5TwGpn6ZgX+St/JXybrAhcnsLnqAn9m+LLJ+hTDGlLSoQdI3gHcDy4F7gY/afqrpFTF0sBXwU0mVhOETgFMlPVd9ke1i6xKHYPs+4DUZ5yPrVwhjiA395fRiFwBH2e6VdDxwFDCnjIqGCrZn17z+aRmNCCGEoZTRs7VdnbjkeuADTa8kNWiwtT0q/yQfP3EiM7bcquFys9/y0kL17fe13QqV22jc84XKLXGx2RbjuhrODwLAxJ5i5Z5c1luoHFvuUajY+IJf382PFZsq+PTzKwqVm9BT7Pty1gHFcrHMGL+8ULmeJx9ouMxrF04oVJcHX9SwtqSFVa/n2Z5XoJoDgV8UKJdL7q3MQwhhJLm/bnaypYNl/ZL0W5J83LXm2v51es1coBeYP9x21hPBNoTQ9ux++nuL9cBt7zzY+5L2B3YHdiozP3e7T/2aIOlPkm6RdLukYzOukaTvSbpH0q2Sth2JtoYQSmTTv2J55jEcknYheSA2O92CqzTt3rN9AXib7X+lScx/L+lS29dXXROJaEIY7ezCPdshnAysBixIp5heb/tTZVTU1sE27dL/K305Lj1qu/mRiCaEUc4lBVvbmzb9pnW09TACgKRuSTeTbDa5wPYfay7JlYgmhNDB3E9f7/LMo1O0fbC13Wd7a5KcB9tL2qLmklyJaKqzfvU+93QZTQ0hlMQljdm2UtsH24p0Cd01wC41b+VKRFOd9atnUrGNBkMIIyQdRsg6OkVbB1tJ0yVNSz+fCOwM3FlzWSSiCWHUM+7vyzw6RVs/IAPWB86W1E3yi+Fc27+JRDQhjC3JA7JiK/LaRVsHW9u3AttknI9ENCGMJXZHPQzL0tbBNoQQ4MUHZJ0sgm0IoQOUtqihZVTiUuC2Jelx4ME6b68N1Ns2vQxRX+fWN5q/trLq29j29EYLSbqMqv0Qayy1XTtLqe2MyWA7GEkLB8sgFPVFfSNR11iob7Rr66lfIYQwWkSwDSGEFohgO1CRDO9RX33bSSptq5EMDX19ks6S9JtW1NUEo72+US3GbEMhkob6wTnb9gGS1gP+YfuFIa5vKklvAa6uOvU8ycKXb9r+cdV1a5D8PyhlR9UQKmLqVyhq/arPdwdOrTn3PIDtR1vZqAyvBp4EJgF7AmdJuruSE9l2ZCUKLRHDCKEQ249WDuCp2nOVICbJlWEESTPT13tJ+p2k5yXdJGkrSVtIuk7Ss5J+L2mT6vokvVvSjZKWSbpf0lcljc/R1CVpe+6zfTxJ4F25KrF2GEHSNfGQtsMAAAG7SURBVJJOkXScpKWSlkj6pqT4vxKGJX6Awkg4FjieJOg9BfwMOAmYC2wPTAC+V7lY0jtJNuI7maSneiDJltPH5a1QUpek9wDTgBuGuHwfks3/Xg98FjgU+FDeukLIEsMIYSR82/YlAJK+BVwMvN/21em5k0kCa8Vc4Bu2z0xf3ytpDvBTSV8YYpO+B9LtTlYjyX18pO2Fg1wP8Bfbx6Sf3y3pE8BOwDn5v8QQVhXBNoyEW6s+fyz9eFvNucmSJqWb8G1Hkjh+TtU1XcBEki2qB0up+VbgHyTBdnvgJEnP2D41Z/sgyY+8ziDXhzCkCLZhJFTnyvMg57qqPh4L/DLjXo8PUdf9titLTm+X9FrgGJIHennaV2lPDLmFYYlgGzrBIuAVtu9pwr36SGYmhNBSEWxDJ/gK8BtJDwLnkjy82gLY3vaRQ5RdR1IPLw4j7AucV2ZjQ8gSwTa0PduXS3oXcDRwBEmwvRs4K0fx29OPvSS7MP8I+HLzWxnC4GIFWQghtEAM+ocQQgtEsA0hhBaIYBtCCC0QwTaEEFoggm0IIbRABNsQQmiBCLYhhNACEWxDCKEFItiGEEIL/H8WQOrI22KgWQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 360x180 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Visualize results for one electrode\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Electrode 1 shows a classic tilt response,\n",
    "# electrode 11 shows a theta response.\n",
    "elec = 1\n",
    "plt.figure(figsize=(5, 2.5))\n",
    "ax = plt.subplot(111)\n",
    "plt.matshow(tstats[:, elec, :], fignum=0, cmap='RdBu_r', aspect='auto', vmin=-2.5, vmax=2.5)\n",
    "plt.yticks(np.arange(len(myfreqs)), np.round(myfreqs))\n",
    "ax.invert_yaxis()\n",
    "plt.xlabel('Time Bin', fontsize=14)\n",
    "plt.ylabel('Frequency (Hz)', fontsize=14)\n",
    "cb = plt.colorbar()\n",
    "cb.set_label('T-stat')\n",
    "\n",
    "channels = load.Load(df_select.iloc[0], 'channels')\n",
    "print(channels.iloc[elec]['ind.region'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Ttest_indResult(statistic=-2.5067252687493564, pvalue=0.012342128124738515)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAADTCAYAAACbfYSmAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAaR0lEQVR4nO3de5xU5Z3n8c9XJIAXVKRB5SI4IdkRtVFaI+skOiEoO0GMLriYRDEhYfOKE5HgqoSM6IzsONEVNXEyQWHECaJozHiZjIJEJROjBglEERG8NyIgiooKcvntH+c0Fngaurq76hTU9/169avqXKrOt9GuXz3nOed5FBGYmZntaK+8A5iZWWVygTAzs0wuEGZmlskFwszMMrlAmJlZJhcIMzPLtHfeAVpL586do1evXnnHMDPbrTz99NNvRURN1rY9pkD06tWL+fPn5x3DzGy3IunVxrb5FJOZmWVygTAzs0wuEGZmlmmP6YMwM2uqTZs2UV9fz4YNG/KOUjbt27ene/futG3btsmvcYEws6pTX1/P/vvvT69evZCUd5ySiwjWrl1LfX09vXv3bvLrXCDMWtMj//jJ878en18O26kNGzZUTXEAkMTBBx/MmjVrinqdC4RZa3rs6k+eu0BUtGopDg2a8/u6k9rMLAeSGDdu3Lbla6+9liuuuAKAK664gm7dutGvXz+OPPJIZs6cmUvGsrQgJE0DhgCrI+KoHbZdDFwD1ETEW+m68cAoYAtwYUQ8VI6cZlalbh3Suu93/gO73KVdu3bcc889jB8/ns6dO39q+9ixY7n44otZtmwZ/fv3Z9iwYUV1MLeGcrUgbgUG77hSUg9gEPBawbojgRFA3/Q1/yypTXlimpmVx957783o0aOZPHnyTvfr06cP++yzD++88w4AL774IoMHD6Z///588Ytf5Pnnnwdg1apVnHnmmdTW1lJbW8vjjz/e4oxlKRARMQ94O2PTZOASoHDe0zOAOyJiY0S8DCwHTih9SjOz8rrggguYMWMG7777bqP7LFiwgD59+tClSxcARo8ezU9/+lOefvpprr32Wr7//e8DcOGFF3LyySezaNEiFixYQN++fVucL7dOaklDgRURsWiHzpNuwBMFy/Xpuqz3GA2MBujZs2eJkpqZlUbHjh0577zzuPHGG+nQocN22yZPnszNN9/MSy+9xIMPPgjA+vXrefzxxxk+fPi2/TZu3AjAb3/7W2677TYA2rRpwwEHHNDifLl0UkvaB5gAXJ61OWNdZKwjIqZERF1E1NXUZA5GaGZW0S666CKmTp3KBx98sN36sWPHsnTpUu68807OO+88NmzYwNatWznwwANZuHDhtp8lS5aULFteVzH9BdAbWCTpFaA7sEDSISQthh4F+3YH3ih7QrOdmDznhW0/Zi3RqVMnzj77bKZOnZq5/ayzzqKuro7p06fTsWNHevfuzV133QUkN8AtWrQIgIEDB/Lzn/8cgC1btvDee++1OFsuBSIinomILhHRKyJ6kRSF4yLiTeA+YISkdpJ6A32Ap/LIadaYG+Yu2/Zj1lLjxo3jrbfeanT75ZdfznXXXcfWrVuZMWMGU6dOpba2lr59+3LvvfcCcMMNN/DII49w9NFH079/fxYvXtziXOW6zHUmcArQWVI9MDEiMstlRCyWNAt4DtgMXBARW8qR08yqVBMuS21t69ev3/a8a9eufPjhh9uWG+6HaNC/f3+WLl0KQO/evbf1SRTq2rXrtmLRWspSICLinF1s77XD8iRgUikzmZnZzvlOajMzy+QCYWZmmVwgzMwsk0dzNSuzUl0aO3bQ50ryvla9XCDMMhTzIV6479hmvodZJfIpJjOzHLRp04Z+/fpx1FFHcfrpp7Nu3bq8I32KWxBmVvVau7XXlNN9HTp0YOHChQCMHDmSm266iQkTJrRqjpZyC8LMLGcDBgxgxYoV25avueYajj/+eI455hgmTpy4bf1tt93GMcccQ21tLeeee27Jc7kFYWaWoy1btjB37lxGjRoFwOzZs1m2bBlPPfUUEcHQoUOZN28eBx98MJMmTeL3v/89nTt35u23s2ZQaF0uEGZmOfjoo4/o168fr7zyCv3792fQoEFAUiBmz57NscceCyRDcixbtoxFixYxbNiwbbPPderUqeQZfYrJzCwHDX0Qr776Kh9//DE33XQTkIzQOn78+G3DeS9fvpxRo0YREewwd07JuUCYmeXogAMO4MYbb+Taa69l06ZNnHbaaUybNm3bYH4rVqxg9erVDBw4kFmzZrF27VqAspxiKkuBkDRN0mpJzxasu0bS85L+LOnXkg4s2DZe0nJJSyWdVo6MZmZ5OfbYY6mtreWOO+7g1FNP5etf/zoDBgzg6KOPZtiwYbz//vv07duXCRMmcPLJJ1NbW8sPf/jDkucqVx/ErcDPgNsK1s0BxkfEZkn/BIwHLpV0JDAC6AscBjws6XMe8tvMSiWPu9ALh/sGuP/++7c9HzNmDGPGjPnUa0aOHMnIkSNLnq1BWVoQETEPeHuHdbMjYnO6+ATJzHEAZwB3RMTGiHgZWA6cUI6cZmb2iUrpg/g28J/p827A6wXb6tN1ZmZWRrkXCEkTSGaOm9GwKmO3aOS1oyXNlzR/zZo1pYpoZlaVci0QkkYCQ4BvRERDEagHehTs1h14I+v1ETElIuoioq6mpqa0Yc1sj/LJR051aM7vm1uBkDQYuBQYGhEfFmy6DxghqZ2k3kAf4Kk8MprZnql9+/asXbu2aopERLB27Vrat29f1OvKchWTpJnAKUBnSfXARJKrltoBc9KbP56IiO9FxGJJs4DnSE49XeArmMysNXXv3p36+nqq6dR0+/bt6d69+653LFCWAhER52SsnrqT/ScBk0qXyMyqWdu2bendu3feMSpe7p3UZmZWmVwgzMwskwuEmZllcoEwM7NMLhBmZpbJBcLMzDK5QJiZWSYXCDMzy+QCYWZmmVwgzMwskwuEmZllcoEwM7NMLhBmZpapLAVC0jRJqyU9W7Cuk6Q5kpaljwcVbBsvabmkpZJOK0dGMzPbXrlaELcCg3dYdxkwNyL6AHPTZSQdCYwA+qav+WdJbcqU08zMUmUpEBExD3h7h9VnANPT59OBrxWsvyMiNkbEy8By4IRy5DQzs0/k2QfRNSJWAqSPXdL13YDXC/arT9d9iqTRkuZLml9NM0OZmZVDJXZSK2Nd5sSxETElIuoioq6mpqbEsczMqkueBWKVpEMB0sfV6fp6oEfBft2BN8qczcys6uVZIO4DRqbPRwL3FqwfIamdpN5AH+CpHPKZmVW1vctxEEkzgVOAzpLqgYnA1cAsSaOA14DhABGxWNIs4DlgM3BBRGwpR04zM/tEWQpERJzTyKaBjew/CZhUukRmZrYrldhJbWZmFcAFwszMMrlAmJlZJhcIMzPL1OQCIWl4I+uHtV4cMzOrFMW0IKY2sn5KawQxM7PKssvLXCUdkT7dK71xrXAojCOADaUIZmZm+WrKfRDLScZCEvDiDtveBK5o5UxmZlYBdlkgImIvAEmPRcTJpY9kZmaVoMl9EC4OZmbVpclDbaT9D5OAfsB+hdsiomcr5zIzs5wVMxbT7SR9EOOAD0sTx6w4k+e8kHcEsz1WMQWiL3BSRGxtzQCSxgLfIekIfwb4FrAPcCfQC3gFODsi3mnN45qZ2c4Vcx/EPODY1jy4pG7AhUBdRBwFtAFGAJcBcyOiDzA3XTYzszIqpgXxCvCQpHtILm/dJiIub2GGDpI2kbQc3gDGk8wfATAdeBS4tAXHMDOzIhVTIPYF7gfasv2UoM0WESskXUsyYdBHwOyImC2pa0SsTPdZKalL1usljQZGA/Ts6X5yM7PW1OQCERHfau2DSzoIOAPoDawD7pL0zSIyTSEd6qOuri5aO5+ZWTUr5jLXIxrbFhEvNfP4XwFejog16THuAf47sErSoWnr4VBgdTPf38zMmqmYU0yFQ240aPjW3qaZx38NOFHSPiSnmAYC84EPgJEk81aPBO5t5vubmVkzFXOKabsrniQdAkwEftfcg0fEk5LuBhYAm4E/kZwy2g+YJWkUSRHJHGrczMxKp5gWxHYi4k1JFwEvkNxE19z3mUhSaAptJGlNmJlZTlo6o9znSS5NNTOzPUwxndS/45M+B0gKQ1/g71s7lJkVr1TDjowd9LmSvK9VvmJOMd2yw/IHwKKIWNaKeczMrEIU00k9vZRBzMyssjS5D0JSW0lXSnpJ0ob08UpJnyllQDMzy0cxp5h+ApwAfA94FTgc+DugIzC29aOZmVmeiikQw4HaiFibLi+VtABYhAuEmdkep5jLXFXkejMz240VUyDuAu6XdJqkv5Q0GPj3dL2Zme1hijnFdAnwY+Am4DBgBTATuKoEuczMLGe7bEFIOknSP0XExxFxeUR8NiL2SWd7awccV/qYZmZWbk05xfQjkulGszwCTGi9OGZmVimaUiD6AQ82su1hoH9LAkg6UNLdkp6XtETSAEmdJM2RtCx9PKglxzAzs+I1pUB0BBq7Ga4tsH8LM9wAPBgR/w2oBZYAlwFz09NYc9NlMzMro6YUiOeBUxvZdmq6vVkkdQS+BEwFSPs51pFMQ9owtMd04GvNPYaZmTVPUwrEZOAXks6StBeApL0knQX8C3BdC45/BLAG+FdJf5J0i6R9ga4RsRIgfezSgmOYmVkz7PIy14i4PZ09bjrQTtJbQGdgAzAxIma28PjHAT9IZ5e7gSJOJ0kaDYwG6NmzZwtimJnZjpp0o1xEXAd0A04HLk4fu0fE5BYevx6oj4gn0+W7SQrGKkmHAqSPqxvJNSUi6iKirqampoVRzMysUDHDfb8HPNSaB0+nLX1d0ucjYinJNKPPpT8jgavTx3tb87hmZrZrzZ6TuhX9AJiRDhv+EvAtkpbNLEmjgNdIBgo0M7Myyr1ARMRCoC5j08ByZzEzs0/kXiCsOpRqvuSqsfD2Xe/T7+ulz2FVpZjRXM3MrIq4QJiZWSYXCDMzy+QCYWZmmVwgzMwskwuEmZllcoEwM7NMLhBmZpbJN8qZ5a0pN8GZ5cAtCDMzy+QWhFlLNdYCcMvAdnNuQZiZWaaKKBCS2qRTjj6QLneSNEfSsvTxoLwzmplVm4ooEMAYYEnB8mXA3IjoA8yliGlIzcysdeReICR1B74K3FKw+gySObBJH79W7lxmZtUu9wIBXA9cAmwtWNc1IlYCpI9dsl4oabSk+ZLmr1mzpvRJzcyqSK4FQtIQYHVEPN2c10fElIioi4i6mpqaVk5nZlbd8r7M9SRgqKS/AdoDHSX9Elgl6dCIWCnpUGB1rinNzKpQrgUiIsYD4wEknQJcHBHflHQNMBK4On28N7eQZrsLT0tqrawS+iCyXA0MkrQMGJQum5lZGeV9immbiHgUeDR9vhYYmGceM7NqV6ktCDMzy5kLhJmZZaqYU0xWGSbPeSHvCGZWIdyCMDOzTG5BWPVplWG4j2+F9zCrbC4QZtXE90pYEXyKyczMMrkFYXsWz+Jm1mrcgjAzs0xuQVj+mvqt3+fGzcrKLQgzM8vkAmFmZpnynjCoh6RHJC2RtFjSmHR9J0lzJC1LHw/KM6eZWTXKuw9iMzAuIhZI2h94WtIc4HxgbkRcLeky4DLg0hxzWiXwFUpmZZVrCyIiVkbEgvT5+8ASoBtwBjA93W068LV8EpqZVa+K6YOQ1As4FngS6BoRKyEpIkCXRl4zWtJ8SfPXrFlTrqhmZlWhIgqEpP2AXwEXRcR7TX1dREyJiLqIqKupqSldQDOzKpR7gZDUlqQ4zIiIe9LVqyQdmm4/FFidVz4zs2qVaye1JAFTgSURcV3BpvuAkSRzUY8E7s0hnplRujlCxg76XEne11pP3lcxnQScCzwjaWG67kckhWGWpFHAa8DwnPKZmVWtXAtERPwXoEY2DyxnFjMz217eLQjb0/neBbPdVu6d1GZmVplcIMzMLJNPMZnZ9jwtqaVcIErMlwia2e7KBcKy+VukWdVzgbDm8xVKZns0F4jdVKlOXZmZNXCBqEb+5m8t5VOQVcGXuZqZWSa3IMysNHbVylgxH85/oDxZrFncgjAzs0wV3YKQNBi4AWgD3BIRV5fqWO70NTPbXsW2ICS1AW4C/gdwJHCOpCPzTWVmVj0quQVxArA8Il4CkHQHcAbwXK6p8uQrR8ysjCq5QHQDXi9Yrge+kFOW0mutS099CavtTm4d0jrv09TO7qYcr9I6znPMrIgoyRu3lKThwGkR8Z10+VzghIj4QcE+o4HR6eLngaWNvF1n4K0Sxm0u5ypepWZzruI4V3FKmevwiKjJ2lDJLYh6oEfBcnfgjcIdImIKMGVXbyRpfkTUtW68lnOu4lVqNucqjnMVJ69cFdtJDfwR6COpt6TPACOA+3LOZGZWNSq2BRERmyX9LfAQyWWu0yJicc6xzMyqRsUWCICI+A3wm1Z4q12ehsqJcxWvUrM5V3Gcqzi55KrYTmozM8tXJfdBmJlZjqquQEi6WFJI6px3FgBJ/yDpz5IWSpot6bC8MwFIukbS82m2X0s6MO9MkFz+LGmxpK2Scr/aRNJgSUslLZd0Wd55GkiaJmm1pGfzzlJIUg9Jj0hakv53HJN3JgBJ7SU9JWlRmuvKvDMVktRG0p8klfUmjaoqEJJ6AIOA1/LOUuCaiDgmIvoBDwCX5x0oNQc4KiKOAV4Axuecp8GzwFnAvLyDVPhwMLcCg/MOkWEzMC4i/hI4EbigQv7NNgJfjohaoB8wWNKJOWcqNAZYUu6DVlWBACYDlwAV0/ESEe8VLO5LhWSLiNkRsTldfILkPpTcRcSSiGjshshy2zYcTER8DDQMB5O7iJgHvJ13jh1FxMqIWJA+f5/kQ69bvqkgEuvTxbbpT0X8LUrqDnwVuKXcx66aAiFpKLAiIhblnWVHkiZJeh34BpXTgij0beA/8w5RgbKGg8n9w253IakXcCzwZL5JEulpnIXAamBORFRELuB6ki+2W8t94Iq+zLVYkh4GDsnYNAH4EXBqeRMldpYrIu6NiAnABEnjgb8FJlZCrnSfCSSnBWaUI1NTc1UIZayriG+dlU7SfsCvgIt2aEXnJiK2AP3S/rZfSzoqInLtw5E0BFgdEU9LOqXcx9+jCkREfCVrvaSjgd7AIkmQnC5ZIOmEiHgzr1wZbgf+gzIViF3lkjQSGAIMjDJeD13Ev1fedjkcjH2apLYkxWFGRNyTd54dRcQ6SY+S9OHk3cl/EjBU0t8A7YGOkn4ZEd8sx8Gr4hRTRDwTEV0ioldE9CL5wz6uHMVhVyT1KVgcCjyfV5ZC6WRNlwJDI+LDvPNUKA8HUyQl39CmAksi4rq88zSQVNNwpZ6kDsBXqIC/xYgYHxHd08+tEcBvy1UcoEoKRIW7WtKzkv5McgqsIi77A34G7A/MSS/B/Ze8AwFIOlNSPTAA+A9JD+WVJe3EbxgOZgkwq1KGg5E0E/gD8HlJ9ZJG5Z0pdRJwLvDl9P+rhem347wdCjyS/h3+kaQPosLG/S4/30ltZmaZ3IIwM7NMLhBmZpbJBcLMzDK5QJiZWSYXCDMzy+QCYVZmkm6VdFXeOcx2xQXCdhuSviFpfcZPSKrEMayKJul8SVvS3+u99D6BIXnnsurkAmG7jYiYERH7Ff4AFwGrgJvzyCSpFMPV/CH93Q4kuet4lqROJTjOTinhz4gq5v/4ttuSdCzJEO4jImJlI/t8VtJjkt6V9JakOwu29ZU0R9LbklZJ+lG6vp2k6yW9kf5cL6lduu2U9M7kSyW9Cfxrun5I+m1/naTHJR1TmFPSAknvp8dv35TfLyK2AtOADsAR6Xt9N52c6G1J9ymdYErSlZJ+mj5vK+kDST9JlztI2iDpoHT5xDTjOiUT5JxSkPXRdHTh3wMfNhzXqpMLhO2W0nFz7gauiohHd7LrPwCzgYNIBtNr+BDdH3gYeBA4DPgsMDd9zQSSyWz6AbUk8z78uOA9DwE6AYcDoyUdR/JB/r+Bg4FfAPelheYzwL8D/5a+5i7gfzbxd9wb+A6wHlgm6cvAPwJnkwwN8SrJHBQAjwGnpM+PB94ETk6XBwBLI+IdSd1IBoS8Ks1zMfArSTUFhz4XGE0y1MqrTclqeyYXCNvtpAO+TScZafMnu9h9E8kH+WERsSEi/itdPwR4MyL+X7r+/YLx/78B/H1ErI6INcCVJB+aDbYCEyNiY0R8BHwX+EVEPBkRWyJiOskMZSemP22B6yNiU0TcTTLWz86cKGkdyYf8OcCZEfFummtaRCyIiI0ks/wNUDKvwh9IBg48GPgSyampbumw2ieTFBCAbwK/iYjfRMTWiJgDzAcKx0O6NSIWR8TmiNi0i6y2B3OBsN3RpcBRwMjCYcglfbGg47ph0LxLSOZteErJXMPfTtf3AF5s5P0PY/tvzq+m6xqsiYgNBcuHA+PSUzbr0g/3HulrDiOZqKpw0LNdfSt/IiIOjIjOEXFiRDyclSudAW0t0C0tVPNJisGXSArC4ySD4xUWiMOB4Ttk/SuSFkmDwkmQrIrtUfNB2J4vPV8+AfhSRKwr3BYRvwP222HdmyTf8JH0V8DDkuaRfAie08hh3iD5IG0oMj3Zfp6HHUe4fB2YFBGTMvKeTPJNXgVFoieNF6edacjV8N77kpzSWpGuegz4MsksbX9Ml08jOUXWMIf368C/RcR3d3Icj+BpgFsQthuRdCjJOfeLIuJPTXzNcCVz+gK8Q/LhtwV4ADhE0kVpX8H+kr6Q7jcT+LGSOQI6k0wD+8udHOZm4HuSvpBe+bOvpK+m/Rx/IJmR70JJe0s6i+QDuzluB74lqV/aaf5/gScj4pV0+2PAecBz6RzZj5L0Ybycnioj/T1Ol3Sakik226cd7xUx57hVFhcI2518F+gK3KBP3wvR2HwVxwNPSlpPMpnPmIh4OSLeBwYBp5Oc618G/HX6mqtITtf8GXgGWJCuyxQR89NsPyMpQsuB89NtHwNnpcvvAP8LaNYsahExF/g7ktnYVgJ/QTKJTIPHSa54amgtPAdsKFgmIl4HziCZgncNSYvi/+DPAsvg+SDMzCyTvzWYmVkmFwgzM8vkAmFmZplcIMzMLJMLhJmZZXKBMDOzTC4QZmaWyQXCzMwyuUCYmVmm/w8HMDgTxb0cRQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Let's take a look at exactly what we were comparing\n",
    "# by drilling down into one of the pixels from the plot above\n",
    "freq_idx = 1  # 5Hz\n",
    "time_idx = 2  # 200ms-300ms after word event\n",
    "\n",
    "rec_pows = all_pows[freq_idx, recalled==1, elec, time_idx]\n",
    "nrec_pows = all_pows[freq_idx, recalled==0, elec, time_idx]\n",
    "\n",
    "plt.figure(figsize=(6, 3))\n",
    "ax = plt.subplot(111)\n",
    "plt.hist(nrec_pows, 30, alpha=0.75, label='NRec', color='C1')\n",
    "plt.hist(rec_pows, alpha=0.5, label='Rec', color='C0');\n",
    "plt.legend()\n",
    "plt.vlines(np.mean(rec_pows), ymin=ax.get_ylim()[0], ymax=ax.get_ylim()[1], color='C0', linewidth=3)\n",
    "plt.vlines(np.mean(nrec_pows), ymin=ax.get_ylim()[0], ymax=ax.get_ylim()[1], color='C1', linewidth=3)\n",
    "plt.xlabel('Z-scored Power', fontsize=12)\n",
    "plt.ylabel('Count', fontsize=12)\n",
    "\n",
    "print(ttest_ind(rec_pows, nrec_pows, equal_var=False))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "supramarginal\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVsAAAC5CAYAAACGNJzQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3deZxcVZn/8c+3u7MvBJKgbBIEBBGRbRDBnVWGxWV0YFBBYBgVHEAQZBhQ9CcuiCPCOBrZVBBFAQURISLIOAxoCLsBDNsYBEIIYAhk6fT398e9BZXqW123btet7up+3q/XfaXq1j33nE46T58+95znyDYhhBDK1TXUDQghhNEggm0IIbRBBNsQQmiDCLYhhNAGEWxDCKENItiGEEIbRLBNSdpb0gOSFkj6bMl1bSTpRknzJd0n6Zgy66uqt1vSHZJ+2Ya6pkn6maT706/zLSXXd1z6d3mvpEsljW/x/S+QtEjSvVXn1pE0R9Kf0z/XLrm+M9O/z7slXSlpWpn1VX12giRLmtGq+kajCLYkQQj4T+A9wFbAQZK2KrHKXuB4268HdgaOKrm+imOA+W2oB+Bs4Ne2twTeVGa9kjYA/hXY0fbWQDdwYIuruQjYu+bcZ4EbbG8O3JC+L7O+OcDWtrcBHgROLrk+JG0E7AH8XwvrGpUi2CZ2AhbYftj2SuDHwAFlVWb7Cdvz0tdLSQLRBmXVByBpQ+DvgfPKrCetayrwduB8ANsrbT9XcrU9wARJPcBE4K+tvLntm4ElNacPAL6fvv4+8N4y67N9ve3e9O2twIZl1pf6D+BEIFY/DVIE28QGwF+q3i+k5OBXIWkWsB1wW8lVfZPkP01fyfUAvBZ4GrgwHbY4T9Kksiqz/TjwdZLe1xPA87avL6u+Kq+y/UTahieAddtQZ8VhwLVlViBpf+Bx23eVWc9oEcE2oYxzpf8klzQZuBw41vbfSqxnX2CR7dvLqqNGD7A98F+2twOW0dpfsdeQjpUeAGwCrA9MkvThsuobapJOIRmKuqTEOiYCpwCnlVXHaBPBNrEQ2Kjq/Ya0+NfQWpLGkATaS2xfUWZdwK7A/pIeJRkiebeki0usbyGw0Halt/4zkuBblt2BR2w/bXsVcAWwS4n1VTwlaT2A9M9FZVco6RBgX+Bgl5vYZFOSH153pd83GwLzJL26xDpHtAi2iT8Cm0vaRNJYkocrV5VVmSSRjGfOt/2NsuqpsH2y7Q1tzyL52n5ru7Sen+0ngb9I2iI9tRvwp7LqIxk+2FnSxPTvdjfa8yDwKuCQ9PUhwC/KrEzS3sBJwP62XyyzLtv32F7X9qz0+2YhsH36bxsKiGALpA8djgauI/lPepnt+0qsclfgIyQ9zDvTY58S6xsKnwIukXQ3sC1wRlkVpT3onwHzgHtIvq9nt7IOSZcC/wtsIWmhpMOBrwB7SPozyRP7r5Rc37nAFGBO+j3znZLrCy2kSLEYQgjli55tCCG0QQTbEEJogwi2IYTQBhFsQwihDSLYhhBCG0SwrSHpyKgv6htudY2G+ka6CLb9tfsbLOrr3PpG8tc2FPWNaBFsQwihDUbloobJXd2e3jUm87MX+lYzuas787OsbDV5jBubfT+A53t7WaunJ/Ozya/brFB9i5bXT+y17LklTJq2TuZn644v9hVqdW/dz55e8iwz18nOqe3elYXqY/Wq+vU9/wIz15qc/WFX9t9zw+qWvZB5/pllLzF90oS65XqXF/v63Jf9f/LZFStZe9zYuuXGr1sst/eKZ57Jru+lFaw9YVzdcl3d9b+v61n4/AsseXF5099oG2mCl9dJWLeYldfZ7peLd7gp9t3X4aZ3jeGzkzduutwYFQtGm85aq1C5XW8ottT+3PuWFSp39JbZP4Aa6X7h6ULl+hY/Xqjc6meKLc/vmpr9Q6aR52/7n0Lllsx/rFC5okH6dZ8q9lv/Ixf8oFC5cdPq/FAbwH4XFdskZCV9HNS1fuZn5/Q92hE7SIzKYBtC6DzdRX+1HCbaNmZbZ0+ln1QlYnlU0p1Vn52c7gf2gKS96tyztD2gQgjDhwRju5R5dIp2PiC7iJo9jmz/o+1tbW9Lktv1CoB0P64DgTekZb6d7hNWq8w9oEIIw0QXEWxzG2CPo0p+1w8Bl6anDgB+bHuF7UeABST7hNUqbQ+oEMLwIcQYZR+dYriM2b4NeMr2n9P3G5BsaFdRb0+wNfaAktTOPaBCCG1SGUboZMMl2B7EK71aKGFPsHQ1zJEA62i4fNkhhDxEBNtBS7eefj+wQ9XpvHuCPSVpvbRXO+AeULZnk2bv37hn/OibXBxCBxsJwXY4rCDbHbjf9sKqc1cBB0oaJ2kTYHPgDxll27oHVAhhaMRshCYMsMfRgaw5hEC6/9dlJJsE/ho4yvbq9D7nSdoxvbS0PaBCCMNHF9mBNk+wlbSRpBslzZd0n6Rj2tDkfto2jGD7oDrnD61z/kvAlzLOH1H1+hmSnVRDCCPcIBY19ALH254naQpwu6Q5tsvc8bmfIR+zDSGERroGMRshnbFUmbW0VNJ8ktlNEWxDCKFagwdkMyTNrXo/O30g3v8+0ixgO+C2VrYvj1EZbNeaOZl9Dn1L0+XW+/QXCtXXteihQuVW3nhxoXKfeu0bCpXrWjypULlVjxf7+p79w9zGF2VYcPWdjS/KcMmtxRLfvH5K/UxbAzn4i/sVKrf2vh8qVO5zW/1TsXIPX1Oo3IMnn9B0mdUr62eIG4gEY7rqPmJabHvHeh++cg9NJlmpeqztvxVqyCCMymAbQugwEt1jiz/PlzSGJNBeYvuKlrWrCRFsQwjDngRdA+SFHrisBJwPzLf9jZY2rAkNg23a0LcB7wBmAROAp4F5wBzbWYsNsu5zAbAvsMj21um5dYCfpPd9FPiQ7WcljQW+C+wI9AHH2L4p456Z5fO0J4TQQSS6xxQLtsCuwEeAe6oyC/6b7V+1pG051e2XSxov6STgMWAOsD8wIy3zeuAM4DFJV1XNex3IRdRk/aJ+1q5/BrD9RpL5s2dJymprZP0KYRSQoHtsV+bRiO3f25btbSpZBtsdaGHgRQ0PAjsDnwKm2v472wfYPtD2nrY3IAm6fwSurFqkkKlO1q96Wbu2Igme2F4EPEfSy60VWb9CGA0k1N2VeXSKgVq6n+332f6F7RVZF9heYPuLwGbA7wvUv0bWLqCStesu4ABJPely3R1YM1dCo/L9SDpS0lxJc5e8uLxAU0MIQ0WC7jFdmUenqDtma/uuvDdJg/EDLWlR4gKSXvNckmGMW0hWgRRWnYhmm/VmRCKaEDqJRHfBB2TDRa4fC5IelPTddPpE9fnpkh4cRP1Ppdm6qM7aZbvX9nHp2MoBwDTgz3nLhxBGlmTMtjvz6BR5++CbkTzcukHS9KrzPcCmg6g/M2uXpImSJqWv9wB666xjjqxfIYwCkjp+GCFvS00yK2AZ8Md0j7Cm1Mn6VS9r17rAvHQN80kk0zYq94msXyGMNiOgZ5t3UYNIZgT8PXAWcIukA4E78lZUL+sXGVm7bD8KbFHnPpH1K4RRRl0qvKhhuMgbbA1guw84Lt2O/HLgW2U1LIQQXpbORuhkzfRsX2b7/PTB2OWtb1L5nnhyKV/+yk1Nlzt9+lmF6vvB54vNn541cUzjizLse/4nC5V77Oe/KVTukRsfLVRu5YurCpXb6+LmE6AAXLLH6YXKzV+6slC5ia/ZsFC5ZTf/slC5w49uPrkSwKKJxdp59eXNT0B6fnmxaZcaAbMR8gbbzYHF1Sds/7ekN5EsQAghhPKIjlrAkCVXsLWdmUOvOilvCCGUJZmN0Nl5swZsvaTr89zE9p6DaYSkacB5wNYk48OHAS8C3wEmkySZOTgrB6WkvYGzgW7gPNsxIyGEkUaia+wIDrZAbY/2cOCnQKsT754N/Nr2P6QZvyaSJL85wfbvJB0GfAY4tbqQpG7gP0mmfS0kmZZ2Vbv3FgohlEyia0yxZxjDxYDB1vYnqt9L+jBwqu2HW9UASVOBtwOHpnWuBFZK2gK4Ob1sDnAdNcEW2AlYUGmPpB+TJKeJYBvCCDIShhGGw4jza0ny414o6Y500cIk4F6StI4AHyQ7Ec0GwF+q3i9Mz4UQRpJ0GCHr6BTDIdj2ANsD/2V7O5JVap8lGbc9StLtwBQga/5N1g5wmUlmqrN+LWd1a1oeQmgLSXSN6ck8OsVwCLYLgYW2K7td/gzY3vb9ad7cHYBL6T9+XClb3ePdEMjcOcL2bNs72t5xPJ09Xy+EUUfQPbYn8+gUjWYj7F9zqgvYTdLW1SdtX1W0AbaflPQXSVvYfoBk+e2fJK1re1G6Q8O/k8xMqPVHYPM05+3jwIFAsS1GQwjDVqVn28katf7nGee+W/PeMOiu4qeAS9KZCA8DHwM+Kumo9PMrgAsBJK1PMsVrH9u9ko4meXjWDVxg+75BtiWEMNykOzUUK9p//8Oh0CjYtmWuhe076b/tzdnpUXvtX4F9qt7/Cmj7fkIhhPZJZiMUDkcXAecCP2hZgwpoNPUrniSFEIaeKDzzwPbNkma1tD0FDLS77s55byJpkqQ3tKZJIYSwppEwG2Ggll4s6SGSZbS/sr2s9oI0ifiHScZYTwI6Yrx0/fXX4rRPvKfpcmedfHWh+s648/xC5bz8xULllv7h5sYXZZiw7rRC5Tbda7NC5TY+6AOFyq14qNialS9/p15K5XKM2/bthcq9dM2PC5Wbusl6hcpNHF/skctx13yu6TI//ZfTCtXFwA/IZkiaW/V+drrn4LAyULB9PfBJ4AzgR5IeIEk6sxxYG9iSZP7rz4E9bd9TcltDCKOVhHrG1vt0se3aZz7DzkC7664ifUiVDim8FdgYmEDSgz0H+K3tp9vR0BDCKCahkZwbocL2rcCtZTYkTSozF3jc9r5V508AzgRm2l6cUS6yfoUw0klQv2fboKguBd5JMtywEPic7WJje4MwnEaXjwHmA1MrJyRtRJLR6/+yCkTWrxBGByHUVWye7QD7H7bVcFiui6QNSTaTPK/mo/8ATqROvgOqsn6l2cIqWb9CCCNJpWebdXSI4dKz/SZJUJ1SOZEuFX7c9l1SVr4ZIDvr15uzLpR0JHAkwAZrTW5Bk0MIbSOhns4esx3ynq2kyjK626vOTQROARrNE8md9as6Ec30SRMKtzeEMAQkNGZs5tEpcvVsJa1je0lJbdgV2F/SPsB4kjHbHwKbAJVe7YbAPEk72X6yqmzurF8hhA4mwSjp2f5V0o8l7dHqBtg+2faGtmeRZO36re0P2F7X9qz0/EKStItP1hR/OetXmsTmQKBwBrIQwnCVDCNkHZ0ib7B9X3rt1ZIek/Q5SRuX2K66JK0v6VcAtnuBStav+cBlkfUrhBFIXahnbObRKfLOs70WuFbSOsBHSPYLO1XSb4HzgSvT2QCDYvsm4KaM87OqXkfWrxBGGY2ARQ1NPSCzvcT22en2NceSbNR4Kckww/9LH2yFEEJrjbapX5LWBT5KknhmFskWNucD65PsG7YTsGdrm9h6PZMnMX2XtzRd7ogjiq1MnvfpLxYqt81nPlKo3G1nXleo3B43/bBQud4Zry1UjlsvL1Rs3KZbFSp3zDtOLlSuqLOWFUu6ctdeJxYqpw/t2/iiDI+8udj39awN9mq6zAtjv16oLqDwoobhIu9shP1JNmB8D/AAyW4NP7T9bNU1dwJ3lNHIEMIoJ0H30A8jSDrG9tmNzmXJ+6PiEmAx8Hbb29j+VnWgTT0MfDXn/UIIoQnCXT2ZR5sdknHu0DwF87Z0PdsvDHSB7ZeAU3Pe72Vp/oMfAK8G+khyUZ4taVuSTR7HA73AJ23/IaN8JKIJYaSToGvodsWWdBDJZrKbSKqeXjoFeCbPPfIG2z0lrbK9RvZsSfsBPbavzHmfLL3A8bbnSZoC3C5pDvA14HTb16YLHr5Gkrmnuv5IRBPCKGDA3UOaXeAWknzeM4Czqs4vBe7Oc4O8rf8CcELG+eUkyWIKB1vbT5B8EdheKmk+Sc4D80oGsLXIXhn2ciIaAEmVRDQRbEMYSSRo/5DBy2w/BjwGNP9kPZW39ZsC92ec/zNQ8FF0f+mmbNsBt5FMLbtO0tdJxpZ3ySiSOxFNCKGTaSjGZ/u3ItlI4RySnWzGkgxfLrM9dcCC5H9A9hxJwK21OUk3etAkTQYuB461/TfgE8BxtjcCjiOZYtavWMa5zEQ0ko6UNFfS3Kefa0mTQwjtUhmzzTra61zgIJKO5gTgCJLg21DeYHsV8B+SXg64kjYjGbsYdC4CSWNIAu0ltq9ITx8CVF7/lGTIoFbuRDTVWb9mTpuSdUkIYThTV/bRZrYXAN22V9u+EHhXnnJ5W3oi8CJwv6RHJD1CkovgJeAzRRpcoSSt1/nAfNvfqPror8A70tfvJvlJUisS0YQwKgh392QebfZiGmvulPQ1SccBk/IUzJsb4XlJbwH2BrYl+fV9HnCd7Xq7KOS1K0m+hXvShREA/wb8M8lmkz0kD+KOhCQRDckUr31s90qqJKLpBi6IRDQhjECDfEDWwimiHyHppB5NMry5EfD+PAVztz4NqtemR8vY/j3ZY68AO2RcH4loQhiFij4ga/EU0femq8WWA6en9z+GJJAPKHfrJe1I8uv8utQMP9j+dDOtDSGEpgxuUUMrp4geQv/AemjGuX7y5kY4juRh2KMkY6nVQweDHUYIIYQGBjX1a9BTRAdYQTaVFq8gOw74tO1vNtPA4WpJ9xQunfrOpsudO/NVher76v/8e6FyNx99bqFy/3fpzwuV6y349fUW/HH730d8o/FFGd513/8WKvflpcU2Xu65vFjKj7H3zSlUbptbfleoXO+BWRN2GnvjLQ07ZZn+c6sjmy6zdEVvobos4fo92xmS5la9n217dtX73FNEB9C2FWRrEU/5QwhDxWZ1X934uNj2jgOUHvRehZUVZJJ2B16y3SfpdcCWwD157pF36tdldECe2hDCyGRgtbOPHFo5RfRmYLykDYAbSHJ7X5SnYN6e7UPAF9OlavcAq6o/tP2t3E3NIOkCoLKl+dbpucj6FUJ4WdFZpi2eIirbL0o6HDjH9tck5crjnTfYHkUy1WG39KhmYFDBluQnw7kkqRYrIutXCAF4pWdbuHzrpogqXXNwMHB4ei5XHM27qGGjxlcVZ/vmNAnNGqeJrF8hBAAz0JhtOx0DnEyyye19kl4L3JinYNNzKSRNB5a0YOVYI5H1K4QADL5n2yq2byYZt0XSq9OO3r/mKZvrAZmkMZLOkPQc8BSwSXr+y5I+XqzZDZWW9Wvps0ta2MwQQtmc9myzjiHU1LBE3tkIpwIfIBmjWFF1/naSp3FlKC3r15S112lpQ0MI5RvEbIRBS3O09DvdzD3yBtuDgX+xfTnJPmEV9wBbNFNhEyLrVwgBqAwjOPNok34zoYDvNXODvGO265Ms1a3V3cQ96pJ0KclMgxmSFgKfI7J+hRCqrO5rfE2J+vVibX+7mRvkDZR/At5G/4D7QSDXHLOB2D6ozkeR9SuEgNvbi80yU1LdhFs1ubgzNbPh40Vpr7ILeL+kLYCPAvvlvEcIIRQ2xDO/uoHJNDlOWy3vPNtfSDoYOCWt9EskPdr32r6+aOUhhJBHZcx2CD1h+wuDuUEzycNHzK/r6/Qt459WZo13D2yrI95drL6jbyhUbtZLjxQqt/P5DxQqV3RMbPPpEwuVO/PYXPvk9bPRcysLlVvRW+wLfOM+hze+KMPTPcVmvczYeJtC5WYvKBaM9txseqFy0xfkyiy4hp6uYh1DM+SLGgr3aCuGfm/gEEJowIZVQ7uqoTZNQdPyJg9/lgHyP9qOiashhNIYs6pv6KYj2B70Sqi8PdsTat6PAbYD3gt8ebCNqJP160ySh28rSbKOfcz2cxllI+tXCCPcMOjZDlreB2RZS2VJs6O/I+uzJl1E/6xfc4CT07m0XyVJ/nBSTf2R9SuEUcDAquGRiKawvCvI6rmBJMvWoKTJHZbUnLvedmUPjVtJluLWejnrl+2VQCXrVwhhBLHNqtV9mUenGGyw/SA5NzsbpMPI3kI9K+vXBlk3qE5Es/i5v5XQxBBCmYYyN0Ir5H1AdgdrPiAT8GpgJnB0Ce2qrvsUkp0aLsn6OONc5l9/ugHcbIAdtty0g/6JQgjJmG3n9GKz5H1A9sua933A08CNZeYikHQIyYOz3erkzx30Rm4hhOFvJIzZ5n1AdmrZDamVzjI4CXiH7RfrXPZy1i/gcZKsX//UpiaGENrENr0l9GwlfRD4PPB6YCfbcwcuUdxgx2xbIs369b/AFpIWppupnQtMAeZIulPSd9Jr15f0K0g2ciMZxrgOmA9cFlm/Qhh5Kj3brGOQ7gXeT7r7QpnyjtmuYoBFDdVsj222EXWyfmVON4usXyGMPmWN2dqeDyANejVuQ3nHbI8HTgOuJumBAryFZNHB50nGb0MIoRSjZsyWZKeEU2x/t+rc7HT/sX1s79/6ppXnjr8sZa3jftt0uWt/9M5C9T2/fHmhcvf3vqpQuWuO36xQuR/f+1Shcm/daEqhclfP/F2hcuOW9Da+KMOi9f+uULn7P1Fsm71HbnqsULmvH39uoXK3/eiHhcoV9YVvnNx0maLx0jYr6/dsZ6QLrCpmp7OPAJD0G5LZU7VOsf2LYi1qXt5guzv9l+wC/AY4q3XNCSGE/voMK+tnbVtse8d6H9revZxWNSfvA7JnSAaRa70PWNy65oQQQrZhuLtuU/L2bD8PnCfpHbwyZrszsDfp3mBlSHeD+EnVqdcCp9n+ZtU1IklEsw/wInCo7XlltSmE0H62B+rZFibpfcA5JAu0rpF0p+29Wl4R+efZXijpAeAY4EMkK7f+RDIH9n/KaFha7wPAtvBy0pnHgStrLnsPsHl6vBn4r/TPEMII0WAYoTDbV9I/ppSimZ0abgFuKbEtjewGPGS79qnDAcAP0hVmt0qaJmk920+0v4khhDKYcnq27ZR7UYOkmZKOlfQtSdPTcztL2ri85q3hQODSjPO5k9GEEDqTnWxrlHV0iryLGrYjSaf4OLAF8E2Sh2bvATYDDi6rgWn9Y4H9SXLa9vs441y/UXNJR1IZXx47uZXNCyGUzGagqV8dIW/P9izg27bfCKyoOv9r4K0tb1V/7wHm2c6aCJorGY3t2bZ3tL2jeiaU1MwQQhn60gdkWUenyBtsdwAuzDj/V6DYzPvmHET2EALAVcBHldgZeD7Ga0MYWQys7F2deXSKvA/IlgNTM85vQclLdSVNJNn25l+qzn0cwPZ3SPIi7AMsIJn69bEy2xNCGALurDm1WfIG26uB0yR9KH1vSa8BvgJcUUrLKhUl6RWn15z7TtVrA0eV2YYQwtDqSx+QdbK8wwjHk6wtXgRMAH7HKz3JU8ppWgghJJzOs+3kMdu8ixqel7QLya/z25ME6XnAdXV2UAghhJYZCfNsGwZbSWOAm4DDbF8PXF92o0rnPlavbD4T157/cGIJjWk9dRXLCT9+7azESI2dMXVmoXInHZeVxrixj25UrJ0zX3i8ULmXzvlJ44syvLDw+ULlThjbXajco7v8e6Fya43LvbZpDa9Zq/lZPZMKfm029I70YGt7laTNSfYdCyGEtrPN6lEyz/aHwOFlNiSEEOoyrO7tyzw6Rd7fH8YCR0jaA5gLLKv+0PanW92wCkmPAkuB1UBvbd7KyPoVwshnoG91Zz8eyhtstwXuTl9vVfNZO/4G3mW7Xt7cyPoVwkhnd1QvNkve2QhvK7shgxBZv0IY4Qy4wxc1DDhmK2kbSUO93bmB6yXdniaTqZUr65ekIyXNlTTXvcX2BAshDBHD6tV9mUenaBRI7wBmVN5IukbSeuU2qZ9dbW9PMlxwlKS313yeK+vXmoloxpfRzhBCSZwOI3TyA7JGwbY2kL2dZAVZ29j+a/rnIpKM6jvVXJIr61cIobP1rXbm0SmGeohgQJImSZpSeQ3sCdxbc1lk/QphhPMomPpl+v9K3s4fJa8Crkxmd9ED/Mj2ryPrVwijTEmLGiSdCewHrAQeAj5m+7mWV0TjYCvgYkmVhOHjge9JerH6Itv7l9E42w8Db8o4H1m/QhhFbOgrpxc7BzjZdq+kr5LsBnNSGRU1Crbfr3l/cRmNCCGERsro2ab5XipuBf6h5ZWkBgy2tkfkr+Tq7mHc1BmNL6wxc8tiayXcV+yb5JkFxRbCTXvN6wuV2/6trytU7r3brV+o3C6vmVao3JyHny1U7m8rxhYqt7y3WH7839ybtYtTY/fNK/Z8d+EfrilUrp16H1pYqJwHXtQwQ9Lcqvezbc8uUM1hQLGsQzkUS/cTQght5r66W+Asrl3GX03Sb0jycdc6xfYv0mtOAXqBSwbbznoi2IYQhj27j77elQXLeveBPpd0CLAvsFuZ+bmH+9Sv8ZL+IOkuSfdJOj3jGkn6lqQFku6WtP1QtDWEUCKbvlUrM4/BkLQ3yQOx/dMtuEoz3Hu2K4B3234hTWL+e0nX2r616ppIRBPCSGcX7tk2cC4wDpiTTjG91fbHy6hoWAfbtEv/Qvp2THrUdvMjEU0II5xLCra2N2v5TesY1sMIAJK6Jd1JstnkHNu31VySKxFNCKGDuY/VvSszj04xrHu2ALZXA9tKmkaymmxr29VLdnMlokkzhh0JoHFTS2lrCKEcTsdsO9mw79lWpEvobgL2rvkoVyKaNbJ+jWlrLp0QwmClwwhZR6cY1sFW0sy0R4ukCcDuwP01l0UimhBGPOO+1ZlHpxjuwwjrAd+X1E3yg+Ey27+MRDQhjC7JA7JVQ92MQRnWwdb23cB2GecjEU0Io4ndUQ/DsgzrYBtCCDAyHpBFsA0hdIDSFjW0jUpcCjxsSXoaeKzOxzOAetumlyHq69z6RvLXVlZ9G9ue2WwhSb+maj/EGott185SGnZGZbAdiKS5A2UQivqivqGoazTUN9IN66lfIYQwUkSwDSGENohg21+RDO9RX307SCptq5EMTX19ki6S9Mt21NUCI72+ES3GbEMhkhp943zf9qGSXg08a3tFg+tbStI7gRurTr1EsvDl67Z/UHXdWiT/D0rZUTWEipj6FYpar+r1vsD3as69BIl1a6IAAAIqSURBVGD7yXY2KsMbgCXAROCDwEWSHqzkRLb9/FA2LoweMYwQCrH9ZOUAnqs9VwliklwZRpA0K31/oKTfSXpJ0h2StpG0taRbJC2T9HtJm1TXJ2k/SbdLWi7pEUlfkpRnB8dFaXsetv1VksD78qrE2mEESTdJ+rakMyQtlrRI0tclxf+VMCjxDRSGwunAV0mC3nPAj4BzgFOAnYDxwLcqF0vai2QjvnNJeqqHkWw5fUbeCiV1SXovMA34Y4PLDybZ/G8X4GjgWOAf89YVQpYYRghD4Ru2fwUg6SzgauADtm9Mz51LElgrTgHOtH1h+v4hSScBF0v6TINN+h5NtzsZR5L7+ETbcwe4HuBPtk9LXz8o6Z+B3YBL83+JIawpgm0YCndXvX4q/fOemnOTJE1MN+HbAdgpDbAVXcAEki2qB0qp+S7gWZJguxNwjqSltr+Xs32Q5Eded4DrQ2gogm0YCtW58jzAua6qP08Hfppxr6cb1PWI7cqS0/skvRk4jeSBXp72VdoTQ25hUCLYhk4wD9jS9oIW3Gs1ycyEENoqgm3oBF8AfinpMeAykodXWwM72T6xQdl1JfXwyjDCR4CfldnYELJEsA3Dnu3rJP09cCpwAkmwfRC4KEfx+9I/e0l2Yf4u8PnWtzKEgcUKshBCaIMY9A8hhDaIYBtCCG0QwTaEENoggm0IIbRBBNsQQmiDCLYhhNAGEWxDCKENItiGEEIbRLANIYQ2+P+B7+3GoiN8SwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 360x180 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Electrode 1 shows a classic tilt response,\n",
    "# electrode 11 shows a theta response.\n",
    "elec = 11\n",
    "plt.figure(figsize=(5, 2.5))\n",
    "ax = plt.subplot(111)\n",
    "plt.matshow(tstats[:, elec, :], fignum=0, cmap='RdBu_r', aspect='auto', vmin=-2.5, vmax=2.5)\n",
    "plt.yticks(np.arange(len(myfreqs)), np.round(myfreqs))\n",
    "ax.invert_yaxis()\n",
    "plt.xlabel('Time Bin', fontsize=14)\n",
    "plt.ylabel('Frequency (Hz)', fontsize=14)\n",
    "cb = plt.colorbar(); cb.set_label('T-stat')\n",
    "\n",
    "print(channels.iloc[elec]['ind.region'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise: Compute the SME for all electrodes in R1002P, and plot time-frequency spectrograms averaged across hippocampal electrodes**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Exercise: Compute the SME for all electrodes in R1001P using MNE**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "environmentname",
   "language": "python",
   "name": "environmentname"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}