LFP_Coherence_and_Granger.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "LFP_Coherence_and_Granger.ipynb",
      "provenance": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/kayteeflick/nma_cn_flow_state/blob/cleaned_up_LFP_analysis/LFP_Coherence_and_Granger.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nxiMhxrvMqQ8"
      },
      "source": [
        "# Flow state notebook for coherence and granger analyses"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "hozFcItEts4v"
      },
      "source": [
        "This notebook shows the step by step processs of generating these analyses for the neuromatch academy summer 2021.  \n",
        "Written by: Carlos Daniel Carrasco, UC Davis, 2021 "
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "UeXo3ijTMIxb"
      },
      "source": [
        "## Imports"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9qWJi-siDfYx"
      },
      "source": [
        "#title \n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "from scipy import signal\n",
        "from scipy.signal import butter, lfilter, filtfilt"
      ],
      "execution_count": 18,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "0TRRGzo9iwu0",
        "outputId": "6a0f87b9-370b-4cda-de2e-482089275842"
      },
      "source": [
        "# package needed to run Granger \n",
        "!pip install nitime\n",
        "import nitime\n",
        "import nitime.analysis as nta\n",
        "import nitime.timeseries as ts\n",
        "import nitime.utils as tsu\n",
        "from nitime.viz import drawmatrix_channels"
      ],
      "execution_count": 3,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Collecting nitime\n",
            "  Downloading nitime-0.9.tar.gz (6.2 MB)\n",
            "\u001b[K     |████████████████████████████████| 6.2 MB 7.7 MB/s \n",
            "\u001b[?25hRequirement already satisfied: numpy in /usr/local/lib/python3.7/dist-packages (from nitime) (1.19.5)\n",
            "Requirement already satisfied: cython in /usr/local/lib/python3.7/dist-packages (from nitime) (0.29.23)\n",
            "Requirement already satisfied: scipy in /usr/local/lib/python3.7/dist-packages (from nitime) (1.4.1)\n",
            "Requirement already satisfied: matplotlib in /usr/local/lib/python3.7/dist-packages (from nitime) (3.2.2)\n",
            "Requirement already satisfied: networkx in /usr/local/lib/python3.7/dist-packages (from nitime) (2.5.1)\n",
            "Requirement already satisfied: nibabel in /usr/local/lib/python3.7/dist-packages (from nitime) (3.0.2)\n",
            "Requirement already satisfied: cycler>=0.10 in /usr/local/lib/python3.7/dist-packages (from matplotlib->nitime) (0.10.0)\n",
            "Requirement already satisfied: python-dateutil>=2.1 in /usr/local/lib/python3.7/dist-packages (from matplotlib->nitime) (2.8.1)\n",
            "Requirement already satisfied: kiwisolver>=1.0.1 in /usr/local/lib/python3.7/dist-packages (from matplotlib->nitime) (1.3.1)\n",
            "Requirement already satisfied: pyparsing!=2.0.4,!=2.1.2,!=2.1.6,>=2.0.1 in /usr/local/lib/python3.7/dist-packages (from matplotlib->nitime) (2.4.7)\n",
            "Requirement already satisfied: six in /usr/local/lib/python3.7/dist-packages (from cycler>=0.10->matplotlib->nitime) (1.15.0)\n",
            "Requirement already satisfied: decorator<5,>=4.3 in /usr/local/lib/python3.7/dist-packages (from networkx->nitime) (4.4.2)\n",
            "Building wheels for collected packages: nitime\n",
            "  Building wheel for nitime (setup.py) ... \u001b[?25l\u001b[?25hdone\n",
            "  Created wheel for nitime: filename=nitime-0.9-cp37-cp37m-linux_x86_64.whl size=4032151 sha256=1c42108f410ed4a90282ab03fd948bbca271d463617ad3c9090a9975dfe7501f\n",
            "  Stored in directory: /root/.cache/pip/wheels/8e/1b/71/d7d584087121ac331c444b913abeb116b527dad2a0f1d9157a\n",
            "Successfully built nitime\n",
            "Installing collected packages: nitime\n",
            "Successfully installed nitime-0.9\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "q8F_e8RsMFRv"
      },
      "source": [
        "## Data Retrieval - spikes"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "-VOCPJJFCk_i"
      },
      "source": [
        "#title Data retrieval - spikes\n",
        "import os, requests\n",
        "\n",
        "fname = []\n",
        "for j in range(3):\n",
        "  fname.append('steinmetz_part%d.npz'%j)\n",
        "url = [\"https://osf.io/agvxh/download\"]\n",
        "url.append(\"https://osf.io/uv3mw/download\")\n",
        "url.append(\"https://osf.io/ehmw2/download\")\n",
        "\n",
        "for j in range(len(url)):\n",
        "  if not os.path.isfile(fname[j]):\n",
        "    try:\n",
        "      r = requests.get(url[j])\n",
        "    except requests.ConnectionError:\n",
        "      print(\"!!! Failed to download data !!!\")\n",
        "    else:\n",
        "      if r.status_code != requests.codes.ok:\n",
        "        print(\"!!! Failed to download data !!!\")\n",
        "      else:\n",
        "        with open(fname[j], \"wb\") as fid:\n",
        "          fid.write(r.content)\n",
        "\n",
        "alldat = np.array([])\n",
        "for j in range(len(fname)):\n",
        "  alldat = np.hstack((alldat, np.load('steinmetz_part%d.npz'%j, allow_pickle=True)['dat']))"
      ],
      "execution_count": 4,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Qh_K2aup8jrn"
      },
      "source": [
        "## Data retrieval - LFP"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "CL13fYEYDl6d"
      },
      "source": [
        " \n",
        "import os, requests\n",
        "\n",
        "fname = ['steinmetz_st.npz']\n",
        "fname.append('steinmetz_wav.npz')\n",
        "fname.append('steinmetz_lfp.npz')\n",
        "\n",
        "url = [\"https://osf.io/4bjns/download\"]\n",
        "url.append(\"https://osf.io/ugm9v/download\")\n",
        "url.append(\"https://osf.io/kx3v9/download\")\n",
        "\n",
        "for j in range(len(url)):\n",
        "  if not os.path.isfile(fname[j]):\n",
        "    try:\n",
        "      r = requests.get(url[j])\n",
        "    except requests.ConnectionError:\n",
        "      print(\"!!! Failed to download data !!!\")\n",
        "    else:\n",
        "      if r.status_code != requests.codes.ok:\n",
        "        print(\"!!! Failed to download data !!!\")\n",
        "      else:\n",
        "        with open(fname[j], \"wb\") as fid:\n",
        "          fid.write(r.content)\n",
        "\n",
        "dat_LFP = np.load('steinmetz_lfp.npz', allow_pickle=True)['dat']\n"
      ],
      "execution_count": 5,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "5K7UT7dyj_6R"
      },
      "source": [
        "`dat_LFP`, `dat_WAV`, `dat_ST` contain 39 sessions from 10 mice, data from Steinmetz et al, 2019, supplemental to the main data provided for NMA. Time bins for all measurements are 10ms, starting 500ms before stimulus onset (same as the main data). The followin fields are available across the three supplemental files. \n",
        "\n",
        "* `dat['lfp']`: recording of the local field potential in each brain area from this experiment, binned at `10ms`.\n",
        "* `dat['brain_area_lfp']`: brain area names for the LFP channels. \n",
        "* `dat['trough_to_peak']`: measures the width of the action potential waveform for each neuron. Widths `<=10` samples are \"putative fast spiking neurons\". \n",
        "* `dat['waveform_w']`: temporal components of spike waveforms. `w@u` reconstructs the time by channels action potential shape. \n",
        "* `dat['waveform_u]`: spatial components of spike waveforms.\n",
        "* `dat['ss']`: neurons by trials. Exact spikes times for each neuron and each trial, reference to the stimulus onset. A (neuron,trial) entry can be an empty list if that neuron did not fire at all on that trial. \n",
        "* `dat['%X%_passive']`: same as above for `X` = {`lfp`, `ss`} but for  passive trials at the end of the recording when the mouse was no longer engaged and stopped making responses. \n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "xFCTeF8P8JGs"
      },
      "source": [
        "## Functions"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "fogo6aADEKV6"
      },
      "source": [
        "# commmon functions\n",
        "\n",
        "def select_brois(data, brois, data_type, selection=True):\n",
        "  '''\n",
        "  Args:\n",
        "    data: a numpy array of the Steinmetz 2019 data\n",
        "    brois: a list of brain regions of interest\n",
        "    data_type: \"LFP\" or \"spikes\"\n",
        "    selection: if True, filters for sessions that contain all brois. If False, filters\n",
        "        for sessions that contain any of the brois\n",
        "  Returns:\n",
        "    an embedded 1-d numpy array, with length = number of filtered session\n",
        "  '''\n",
        "  selected_data = np.array([])\n",
        "\n",
        "  if data_type == \"LFP\":\n",
        "    brain_area = \"brain_area_lfp\"\n",
        "  elif data_type == \"spikes\":\n",
        "    brain_area = \"brain_area\"\n",
        "\n",
        "  if selection == True:\n",
        "    for i in range(len(data)):\n",
        "      if all(item in data[i][brain_area] for item in brois):\n",
        "        selected_data = np.hstack((selected_data, data[i]))\n",
        "  else:\n",
        "    for i in range(len(data)):\n",
        "      if any(item in data[i][brain_area] for item in brois):\n",
        "        selected_data = np.hstack((selected_data, data[i]))\n",
        "\n",
        "  return selected_data\n",
        "\n",
        "def sel_neurons():\n",
        "  print('todo')\n",
        "  return\n",
        "\n",
        "def spks_to_rate(spks):\n",
        "  print('todo')\n",
        "  return\n",
        "\n",
        "def butter_bandpass_backend(lowcut, highcut, fs, order=6):\n",
        "  nyq = 0.5 * fs\n",
        "  low = lowcut / nyq\n",
        "  high = highcut / nyq\n",
        "  b, a = butter(order, [low, high], btype='band')\n",
        "  return b, a\n",
        "\n",
        "def butter_bandpass_filter_good_func(data, lowcut, highcut, fs, order=6):\n",
        "  b, a = butter_bandpass_backend(lowcut, highcut, fs, order=order)\n",
        "  y = filtfilt(b, a, data)\n",
        "  return y \n",
        "\n",
        "def get_behavioral_idx(recording_session):\n",
        "  '''\n",
        "  Args:\n",
        "    alldat: spiking data as loaded in originally \n",
        "    recording session: what recording session as an integer you want to pull behavioral indices for \n",
        "  Returns:\n",
        "    series of arrays that you can use to index behavior \n",
        "  '''\n",
        "  response = recording_session['response'] # right - nogo - left (-1, 0, 1)\n",
        "  vis_right = recording_session['contrast_right'] # 0 - low - high\n",
        "  vis_left = recording_session['contrast_left'] # 0 - low - high\n",
        "  gocue_idx = recording_session['gocue']\n",
        "  rt_idx = recording_session['reaction_time']\n",
        "  is_correct = np.sign(response)==np.sign(vis_left-vis_right)\n",
        "  return response, vis_right, vis_left, gocue_idx, rt_idx, is_correct\n",
        "\n",
        "def concatenate_lfp(brain_area_raw_dat,zeropadsize=0): \n",
        "  '''\n",
        "  Args:\n",
        "    data: a 2d array of raw brain area LFP should be 2d after you select a specific session and area to pull. \n",
        "    zeropadsize: if you'd like to add padding to the concatenation. Default is 0 (no padding)\n",
        "  Returns:\n",
        "    1D Array of Concatenated trials for a given brain area recording. \n",
        "  '''\n",
        "  concatenated_data = np.array([])\n",
        "  N = zeropadsize\n",
        "  for i in range(len(brain_area_raw_dat)):\n",
        "    single_trial = brain_area_raw_dat[i,:]\n",
        "    trial_padded = np.pad(single_trial,(N,N),'constant')\n",
        "    concatenated_data = np.concatenate([concatenated_data,trial_padded])\n",
        "  return concatenated_data\n",
        "\n",
        "def spikes_avg(data):\n",
        "  '''\n",
        "  Args:\n",
        "    data: a single session data frame from Steinmetz 2019 data\n",
        "  Returns:\n",
        "    A 2d numpy array of trial averages. Rows are neurons. Column are time step\n",
        "  '''\n",
        "  avg = np.mean(data['spks'], axis = 1)\n",
        "\n",
        "  return avg\n",
        "\n",
        "def add_avg(data):\n",
        "  '''\n",
        "  Adds the trial average to the data structure\n",
        "\n",
        "  Args: \n",
        "    data: a single session dataframe from Steinmetz\n",
        "  '''\n",
        "  data[\"spks_avg\"] = spikes_avg(data)\n",
        "  return data\n",
        "\n",
        "def multi_add_avg(alldat):\n",
        "  \"\"\"\n",
        "  Adds trial averages to all sessions in the dataframe\n",
        "  Args:\n",
        "    alldat: a numpy array of sessions\n",
        "  \"\"\"\n",
        "  for i in range(len(alldat)):\n",
        "    alldat[i] = add_avg(alldat[i])\n",
        "  return alldat\n",
        "\n",
        "def concat_trials(data, bin_start, bin_end):\n",
        "  \"\"\"\n",
        "  Args:\n",
        "    data: a single session data frame from Steinmetz 2019 data\n",
        "    bin_start, bin_end: integers indicating which bin to start and end at\n",
        "  Returns:\n",
        "    A 2d numpy array of with trials concatenated horizontally with shape N×TK \n",
        "    (number of neurons by number of time points times number of trials).\n",
        "  \"\"\"   \n",
        "  NN = len(data['spks'][0])\n",
        "  X = np.reshape(data['spks'][:,:,bin_start:bin_end], (NN,-1))\n",
        "\n",
        "  return X\n",
        "\n",
        "def trial_result(data):\n",
        "  '''\n",
        "  Args:\n",
        "    data: a single session data frame\n",
        "  Returns:\n",
        "    numpy logical arrarys corresponding to right trials, correct right, incorrect\n",
        "    right trials, miss right trials, no_go trials\n",
        "  '''\n",
        "  # response: right - nogo - left (-1, 0, 1)\n",
        "  # vis left/right 0 - low - high\n",
        "  response, vis_right, vis_left, gocue_idx, rt_idx, is_correct = get_behavioral_idx(data)\n",
        "  ind = {}\n",
        "\n",
        "  # trial types\n",
        "  ind[\"right_go\"] = np.array(vis_right > vis_left)\n",
        "  ind[\"left_go\"] = np.array(vis_left > vis_right)\n",
        "  ind[\"nogo\"] = np.logical_and(vis_left == 0, vis_right == 0)\n",
        "  ind[\"equal\"] = np.logical_and(vis_right >0, vis_right == vis_left)\n",
        "\n",
        "  # right trials by response\n",
        "  ind[\"right_cor\"] = np.logical_and(response == -1, vis_right > vis_left)\n",
        "  ind[\"right_error\"] = np.logical_and(response == 1, vis_right > vis_left)\n",
        "  ind[\"right_miss\"] = np.logical_and(response == 0, vis_right > vis_left)\n",
        "\n",
        "  # left trials by response\n",
        "  ind[\"left_cor\"] = np.logical_and(response == 1, vis_right < vis_left)\n",
        "  ind[\"left_error\"] = np.logical_and(response == -1, vis_right < vis_left)\n",
        "  ind[\"left_miss\"] = np.logical_and(response == 0, vis_right < vis_left)\n",
        "\n",
        "  # nogo trials by response\n",
        "  ind[\"nogo_cor\"] = np.logical_and(response == 0, ind[\"nogo\"])\n",
        "  ind[\"nogo_error\"] = np.logical_and(response != 0, ind[\"nogo\"])\n",
        "\n",
        "  return ind\n",
        "\n",
        "def frequencyband_coherence_analysis(raw_dat,correct_idx,region_of_interest,timepoints,freq_band,corr):\n",
        "  '''\n",
        "  Args:\n",
        "    raw_dat = array of LFP data for a specific session you would like to look at should be (brain_area,340,250)\n",
        "    correct_idx = idx of correct or incorrect trials to look at one vector so dimension should match second dimension ie. len = 340 \n",
        "    region_of_interest = (i = 1:12) what brain area you'd like to look at in a digit form. if using loop can go through all brain areas and do pair wise comparison. \n",
        "    timepoints =  (tp = [0,50]) what timepoints you'd like to run the coherence analysis over\n",
        "    freq_band =  (fb = [5,8]) --> Theta, What frequence band you'd like to average coherence analysis over. \n",
        "  Returns:\n",
        "    should return a vector for each comparison you chose to do. \n",
        "  '''\n",
        "\n",
        "  coherence_array = np.empty((raw_dat.shape[0]))\n",
        "  zeropadsize = 0 \n",
        "\n",
        "  for i in range(raw_dat.shape[0]):\n",
        "\n",
        "    if corr == 1:\n",
        "      brain_area_1 = raw_dat[region_of_interest,correct_idx==1,timepoints[0]:timepoints[1]] # 1 is correct response\n",
        "      brain_area_2 = raw_dat [i,correct_idx==1,timepoints[0]:timepoints[1]]\n",
        "    else:\n",
        "      brain_area_1 = raw_dat[region_of_interest,correct_idx==0,timepoints[0]:timepoints[1]] # 1 is correct response\n",
        "      brain_area_2 = raw_dat [i,correct_idx==0,timepoints[0]:timepoints[1]]\n",
        "\n",
        "    brain_area_1 = concatenate_lfp(brain_area_1,zeropadsize)\n",
        "    brain_area_2 = concatenate_lfp(brain_area_2,zeropadsize)\n",
        "    f, Cxy = signal.coherence(brain_area_1, brain_area_2, 100, nperseg=100)\n",
        "    Cxy_avg = np.mean(Cxy[freq_band[0]:freq_band[1]])\n",
        "    coherence_array[i] = Cxy_avg\n",
        "\n",
        "  return coherence_array\n",
        "\n",
        "\n",
        "def granger_causality_analysis(averaged_data,TR,freq_band): \n",
        "  '''\n",
        "  Args:\n",
        "    averaged_data = array of averaged LFP Data over a specific ROI and then the timeseries should look like if taking all brain regions (12,250)\n",
        "    TR = the sampling rate of your signal in seconds, for our case this is .01 = 100 hz = 10 ms \n",
        "    freq_band = the lower and upper frequency bands you'd like to look at [5,8] for theta . \n",
        "  Returns:\n",
        "    G which is a granger class with a bunch of different parameters related to your analysis including causality metric.\n",
        "    g2 which is a form of granger where we are examing the difference between Fx→y and Fy→x to get a clearer understanding of which is leading. \n",
        "  '''\n",
        "  f_ub = freq_band[1]\n",
        "  f_lb = freq_band[0]\n",
        "  pdata = tsu.percent_change(averaged_data)\n",
        "  time_series = ts.TimeSeries(pdata, sampling_interval=TR)\n",
        "  G = nta.GrangerAnalyzer(time_series, order=1)\n",
        "  freq_idx_G = np.where((G.frequencies > f_lb) * (G.frequencies < f_ub))[0]\n",
        "  g2 = np.mean(G.causality_xy[:, :, freq_idx_G] - G.causality_yx[:, :, freq_idx_G], -1)\n",
        "  return G, g2\n",
        "\n",
        "\n",
        "def get_high_firing_spikes(spike_data, min_rate, sr=100):\n",
        "  '''\n",
        "  Takes spike data and returns only the neurons with a high enough firing rate\n",
        "  averaged across all trials/time\n",
        "\n",
        "  Parameters\n",
        "  --------------------------------\n",
        "  spike_data : 3d array of spikes with dimensions (neurons, trials, time)\n",
        "  min_rate : The firing rate to filter above (in spikes/sec)\n",
        "  sr : sampling rate of data\n",
        "\n",
        "  Returns\n",
        "  --------------------------------\n",
        "  spikes_subset : 3d array of spikes with dimensions (neurons (filtered), trials, time)\n",
        "  '''\n",
        "  spike_counts = np.sum(np.reshape(spike_data, (spike_data.shape[0], -1)), axis=1)\n",
        "  total_time = spike_data.shape[1] * spike_data.shape[2] / sr\n",
        "  spike_rate = spike_counts / total_time\n",
        "  spikes_subset = spike_data[spike_rate > min_rate]\n",
        "  return(spikes_subset)\n"
      ],
      "execution_count": 55,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gzh18ak_kLgf"
      },
      "source": [
        "#LFP Coherence Analysis "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "B0TIpOSrkNmG",
        "outputId": "87b3037d-9d16-43ba-b3e8-77cf909e3a96"
      },
      "source": [
        "#@title Pulling relevant indices \n",
        "temp_dat = alldat[11]\n",
        "idx = trial_result(temp_dat) \n",
        "print(idx.keys())\n",
        "# don't have enough right error trials so collapsing across errors \n",
        "no_go_error = idx['nogo_error']\n",
        "right_error = idx['right_error']\n",
        "no_go_error = no_go_error.astype('uint8')\n",
        "right_error = right_error.astype('uint8')\n",
        "#Indices of interest \n",
        "correct_right_idx = idx['right_cor']#only at correct trials and for regions of interest\n",
        "error_idx = right_error + no_go_error\n"
      ],
      "execution_count": 45,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "dict_keys(['right_go', 'left_go', 'nogo', 'equal', 'right_cor', 'right_error', 'right_miss', 'left_cor', 'left_error', 'left_miss', 'nogo_cor', 'nogo_error'])\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "t_NkisQPmlni"
      },
      "source": [
        "#@title Running coherence @ Theta Right and Incorrect Trials \n",
        "#pulling data and idx of interest\n",
        "session11 = dat_LFP[11]\n",
        "session11_lfp = session11['lfp']\n",
        "session11_lfp_roi = np.delete(session11_lfp,[0,8,9],axis=0)\n",
        "#parameters \n",
        "timepoints = [0,50]\n",
        "fb = [5,8]\n",
        "\n",
        "#running analysis \n",
        "theta_array = np.empty([9,9])\n",
        "theta_array_incorr = np.empty([9,9])\n",
        "for i in range(theta_array.shape[0]): \n",
        "\n",
        "  theta_array[i,:] = frequencyband_coherence_analysis(session11_lfp_roi,correct_right_idx,i,timepoints,fb,1)\n",
        "  theta_array_incorr[i,:] = frequencyband_coherence_analysis(session11_lfp_roi,error_idx,i,timepoints,fb,1)"
      ],
      "execution_count": 46,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 321
        },
        "id": "tpNlq2W5nAIv",
        "outputId": "87145d2f-ae9c-4fe6-e596-0db68e515156"
      },
      "source": [
        "#@title Plotting Right Correct Coherence \n",
        "positions = (0,1, 2, 3,4,5,6,7,8)\n",
        "labels = session11['brain_area_lfp']\n",
        "roi_labels = np.delete(labels,[0,8,9])\n",
        "plt.xticks(positions, roi_labels,rotation = 'vertical')\n",
        "plt.yticks(positions, roi_labels)\n",
        "plt.imshow(theta_array)\n",
        "plt.colorbar()\n",
        "plt.title('Theta Coherence Right Correct Trials Pre Stim')\n"
      ],
      "execution_count": 47,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "Text(0.5, 1.0, 'Theta Coherence Right Correct Trials Pre Stim')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 47
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUEAAAEeCAYAAAAU1qcTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deZxcRb338c83CZCExIQlsgUSdoSAECIIooKALIq4k3D1gqLR676LyuPKFb1wUXjcQOVBFkGURdRIECUCCrIGQlgEAU1AkbAoBJJMMt/nj6pJTjrTMx36nOme6d/79Tqv6bN0VXVP96/rVNU5JduEEEKnGtbqAoQQQitFEAwhdLQIgiGEjhZBMITQ0SIIhhA6WgTBEEJHG5AgKOmLks4biLzKJmmOpHe1uhwDSdIzkrZp8FhL2q7qMg12kj4r6QcNHnu2pBOrLlOVJM2XtH+ry9GIUoJg/tL0LN2Snius/0eTaTcdQCVtJumHkv4u6WlJ90j6kqT1m0l3MJK0f/4fPZPfi3slvaN4jO0xth8oIa9jJV3XwHGHSLoml+cxSb+X9Lpm829WXz+Akl5e+Iwvzj8Gxe/BVsXjbX/VdqU/prkMi3P+D0s6VdLwktI+Ln9vnpb0qKRZksbmfWsEbdu72J5TRt5VKyUI5i/NGNtjgL8BRxS2nV9GHs+XpA2B64FRwD62xwIHA+OBbQewHJLULs0Pj+T/1QuAjwLfl7RjKwoi6c3AT4FzgInAJsDngSOeR1ojGtlWBtvXFj7zu+TN4wuf+79VXYY6XpzLdCBwNPDu2gPWtjySXgl8FZiRvz8vAn5SQlnbg+1SF+Ah4KCabV8ELiJ90J8G5gPTCvs3By4GHgMeBD6Utx8KLAO6gGeA2/P2dwB357QeAN7TR3lOBOYBw/o4Zl/gJuBf+e++hX1zgK8Af8j5XQlsXNj/UuCPwFPA7cD+Nc/97/zc54DtgJ2A3wBPAPcCby0cfzbwbeBXOa8/AdsW9u9SeO6jwGfz9mHA8cBfgMfze71hnde6P7CwZts/gbcU1g1slx9vBPwC+Hd+b04Erqs59r3Affk9+DYg0hdlCbAi/++e6qUsIv1ofrKP/80w4ATgr7mc5wDj8r7JOf/jcjrXAMfm9/sb+b04EVgPOCUf8yjwPWBUIY8jgbn5Nf6F9Ln771z2Jbn83+qjjD3lGFH4vP8MOC+n+a687bzCc34K/IP0mbsG2KXmc3Bifrwx8Mv83j4BXEudz3Lx/1bI41u9vU95/ztJ36MngdnApDrpfgK4rM6+maTv57L8Pv2iNg7k1/7T/H48Tfo+7gB8Jv9PFwCvLjsWNRyzSk+wfhBcAhwODAdOAm4ofMhvIf36rwtsQwpshxSee15Neq8h1eIEvBJ4Fphapzw3AF/qo7wb5g/B24ERwIy8vlHePyd/MXYg1SbnAF/L+7YgfdEOz6/j4Lw+ofDcv5GC1whgXP6HvyOv7wEsAnYufPgfB/bK+88HLsz7xgJ/Bz4OjMzre+d9H86vcyLpC38GcEGd17s/OQjmMr8O6Ab26O3LBFyYl9HAzrn8tUHwl6Sa9VakH7JD875ji8f2Upad8vO37uOYdwL358/FGOAS4Nya4HMOsH7+/xwLLAc+mN/DUaSAeHn+X48lBfWTchp7kQLRwfn92ALYqfD/e1cDn/mechSDYBfw+pzmKNYMgu/MZVkP+CYwt7DvbFYFwZNIQXudvLwcUJ1yFP9vO5OC7HF13qcj8/v6ovw+nQD8sU66Lyf9iH8JeBmwXs3+leXtLQ6w6vt/SM7rHFJl53P5Nb0beLDKQNfn/6/0BOsHwasK6zsDz+XHewN/qzn+M8D/Kzz3vH7yvAz4cJ199wHv7eO5bwdurNl2PXBs4YtwQmHf+4Ar8uNPk7+Qhf2zgWMKz/1yYd9RwLU1x58BfKHwYfpBYd/hwD358Qzgtjqv4W7gwML6ZqQv4Yhejt2fFPSeApaSajsfqTnGpFrr8JzOjoV9vdUE9yusXwQcnx8fS99B8GX5+SP7OOa3wPsK6zv2vDZWfbm3Kew/tvh5Iv1QLmb1GvU+PV+6/P5/o07ec3j+QfCaXr4DvX6OST8gZlUN92xWBcEvAz+nUMProxwm1TyfJP1wn0gKwr29T78GjiusDyNVJibVSfsw0o/HU6Qa36nA8NryFo5/iNWD4G8K+47IafQ8f2wu3/j+XmMVy0C2Uf2j8PhZYGRum5gEbC7pqZ4F+CypbahXkg6TdIOkJ/Lxh5NOG3rzOCko1LM56VSr6K+kGkG9so/JjycBb6kp+341+S0oPJ4E7F1z/H8AmzaQ15akD3ZvJgGXFtK8mxTc6r2Hj9geT2oTPB14VZ3jJpCCTfE1LOjluHpl7s/j+e/a/H/+mstUfG21ZSquTyDVYm8pvD9X5O3Q9/vajN7eJwAkDZf0NUl/kfRvUsCA3j/DJ5NqbFdKekDS8f3kO9X2Bra3tX2C7e46ZZoEnFZ4T54g/WAUP/cr2f617SNItekjST82a9PR82jh8XPAItsrCuvQ+OemVO3QUL+A9Ks8vrCMtX143u/iwZLWI7UfngJskr/Ms0j/wN5cBbyhj06JR0gfiKKtgIcbLPu5NWVf3/bXCse45vjf1xw/xvZ/NZhXvWErC4DDatIdabvP12B7Kak2u6uk1/dyyGOkU8uJhW1bNlDWlVn0s//eXPY39XFM7f9nq1ym4peqNp/i+iLSl2yXwnszzqnzgJx/vQ6y/srfl76eezQpkBxEaiKZnLev8Rm2/bTtj9vehtR08TFJB5ZQpgWktvTiZ2aU7T/2mYDdbfu3wO+AKb2kO+i0QxC8EXha0qcljcq/klMkvSTvfxSYXAhi65LaUR4Dlks6DHh1H+mfSqrx/EjSJABJW+ThA7uRAugOko6WNELSUaTT9V82UPbzgCPyEI/hkkbmISgT6xz/y5zX2yWtk5eXSHpRA3n9EthM0kckrSdprKS9877vAf9deH0TJB3ZQJrYXgb8L6lNtnbfClIb3BcljZa0E/CfjaSbPQpMlLRunbwNfAz4P5LeIekFkoZJ2k/SmfmwC4CPStpa0hhSL+VPbC9v8PV1A98HviHphbDy/39IPuSHwDskHZjz3iK/zp7yNzReci2NJTVFPE6qpX613oGSXitpO0kitV2uIDVnNOt7wGck7ZLzGSfpLXXKcKSk6ZI2yKMc9iK1xd+QD6nqfRoQLQ+C+Yv2WmB3UmPpIuAHpF9ISL1KAI9LutX208CHSG1PT5J+VS/vI/0nSL2/XcCfJD1Namf6F3C/7cdz/h8nfSg/BbzW9qIGyr6A9Iv+WVJQXgB8kjrvay77q4HppBrOP4Cvk4J6f3k9TWq8PyI/7z7ggLz7NNJ7cGV+fTeQ2lobdRawlaTehqV8gPS/+AdwLikoLW0w3d+RRgL8Q1Kv76ftn5HaSt9Jek8eJbVl/bxQtnNJPagPkhrYP9hg/j0+TTqlvCGffl5FalvE9o2kjqpvkD4Tv2dVzfM04M2SnpR0+lrm2ZdzSKf1DwN3sSqY9Gb7XN5nSG3V37F9dbMFsH0p6bN3YX5P7iS1+/XmSVLnxX2kNsfzgJO9avjbD4Gd86n1Zc2WbaApN0yG0BBJXwc2tX1Mq8sSQhlaXhMM7U3STpJ2K5wGHQdc2upyhVCWCIKhP2NJ7YKLSVcJ/C+rTlVDGDCSzpL0T0l31tkvSadLul/SHZKmNpRunA6HEAYDSa8gtY2eY3tKL/sPJ7UXH05qEz/Ndr9t41ETDCEMCravIY1nrOdIUoC07RuA8ZL6GoMKRBAMIQwdW7D6gPCF1Bn8XTSQd7co1cYbDvfkLdepLP17Fk7o/6AmDP/3kkrTB/DIXofnlZhBxckPrzf+vRzdI6pNH8AVf8OGdVWX9tLFT9C1dHFTb9IhB6zvx59Y0f+BwC13LJ1PGgLV40zbZ9Y7viyDNghO3nIdbpy9NhcvrJ19P/beytIGGH/VfZWmD9C1U3XvD8Cw5WWM2a1v2bhqg/iSDav/+D83odqTrfX/3liAeT7m/ea0ptNY9MQK/jS73rUDq1tns78ssT2tieweZvUrmibSwJVfcTocQqiQWeHuhpYSXA78Z+4lfinwL9t/7+9Jg7YmGEJofwa6S2o3kXQB6S5IG0taCHyBdCsubH+PdAns4aSrg54lXQnUrwiCIYTKGNPlck7Zbc/oZ7+B969tuhEEQwiVKqsmWJUIgiGEyhhYEUEwhNDJ2r0mWHrvsKRn6mx/W76eb76k2yX9QNL4Xo7bX1Ij9/ILIbQ5AyvshpZWGZCaoKRDSVM7Hmb7YaW5UI8h3SL9qYEoQwhh4BnT1eY1wYE6Hf4c8Ime273nG6me1bMzB8lvkrq1+52sO4QwSBhWtHcMHLDB0rsAt/a2Q9JI0u3PjwD2ZPVJh0IIg1gaJ9jY0ioDfsWIpF0lzVWaaeso0tyzD9q+L4/zOa+P586UdLOkmx97vLrLhUIIZRErGlxaZaCC4HxgKoDtebZ3J817OmptErF9pu1ptqdN2Gh4BcUMIZTJQLcbW1ploILgScApNbOw9QTAe0izyfVMe9jnqPAQwuDS7jXBKjpGRufr+nqcavtUSROAX+ee4adIs1vNtr1E0kzgV5KeBa4l3dI9hDDIGehye9+npfQgaPf+im3/CPhRnX1XkNoGQwhDSLpipHW1vEbEFSMhhMoYsaLN79gXQTCEUKluR00whNCh4nQ4hNDRjOiqeqKVJrV36UIIg17UBEMIHcsWKzptiMxAuWfhhEpnhPvjqd+rLG2Al36q2tnsAMY+VO20nitGVvvxGba82ssI9vzYbZWmDzD7t1MrTX/CrUsrS3tYVzlX9HZHTTCE0KlSx0jUBEMIHStOh0MIHSxdNtfeNzuJIBhCqExcMRJC6HjdcTocQuhU0TESQuhoRqyIa4dDCJ2su81rglXMO/y5PLfwHXkukb0lPSRp48IxK+cWlnSspMfysfMl/UzS6LLLFUIYeLbo8vCGllYpNQhK2gd4LTDV9m7AQcCCBp76E9u7294FWAYcVWa5QgitkSZfH9bQ0iplnw5vBiyyvRTA9iIAqbE2AUkjgPWBJ0suVwihRdq9Y6Ts0l0JbCnpz5K+I+mVDT7vKElzgYeBDYFflFyuEEILGNHtxpZWKTUI2n6GNIH6TOAx4CeSjiXVitc4vPD4J3kazk2BecAne0u/OO9w15LFZRY9hFCRFQxraGmV0nO2vcL2HNtfAD4AvAl4HNigcNiGwKJenmtSLfAVddJeOe/wOiPXL7voIYSSpZuqdlbHyI6Sti9s2h34KzAHeHs+ZjjwNuDqOsnsB/ylzHKFEFojTb4+rKGlVcruGBkD/F9J44HlwP2kU+Mu4LuSbgcEXAGcV3jeUZL2IwXlhcCxJZcrhNAiHXVnadu3APvW2X10neecDZxdZjlCCO3BVqm1PEmHAqcBw4Ef2P5azf6tSPObj8/HHG97Vl9pxhUjIYRKlTUGMDelfRs4mHTGeJOky23fVTjsBOAi29+VtDMwC5jcV7rtPYAnhDComXR7/UaWBuwF3G/7AdvLgAuBI3vJ8gX58Tjgkf4SjZpgCKEyRnR1l9bzuwWrX4G2ENi75pgvAldK+iDpwouD+ks0aoIhhEqtxTjBjXvGAedl5vPIbgZwtu2JwOHAuZL6jHNREwwhVKbnipEGLbI9rY/9DwNbFtYn5m1FxwGHAti+XtJIYGPgn/USjZpgCKFS3QxraGnATcD2kraWtC4wHbi85pi/AQcCSHoRMJJ09Vpdg7YmOPzfSxh/1X2VpV/1vMA3/E+18xoDHLL57pWmr/XWqzR9L1tWafoPHL11pekDbPfIHZWmP2zCRpWlrWUrmk7DprSbqtpeLukDwGzS8JezbM+X9GXgZtuXAx8Hvi/po6ROkmPzlWh1DdogGEIYHMq8OUIe8zerZtvnC4/vAl62NmlGEAwhVKbn2uF2FkEwhFCZdO1wB102F0IIqyv3srkqRBAMIVSqwatBWiaCYAihMmX2DlclgmAIoTJGLC/vsrlKRBAMIVQqTodDCB1rMPQOl9JtI+lqSYfUbPuIpF9LujOvj5Z0vqR5ku6UdJ2kMWXkH0JoX+1+e/2ycr6AdB1f0XTgpML6h4FHbe9qewrpQueukvIPIbSjBqfbHApTbv4MeE2+qBlJk4HNWf3eX5tRuOOD7XttL5U0WdI9uZZ4t6SfSRpdUrlCCC1U8k1VK1FKELT9BHAjcFjeNB24iNXnFj4L+LSk6yWdWDMr3Y7Ad2y/CPg38L7e8inOO7yse0kZRQ8hVMjA8u5hDS2tUmbOxVPi6Xl9JdtzgW2Ak0nzDt+Ub3UDsMD2H/Lj80jTbq6hOO/wusNGllj0EEJVOuV0GODnwIGSpgKj88xzq7H9jO1LbL+PFOwO79lVe2iJ5QohtEjPTVU7IgjafoY0ofpZ1NQCASS9TNIG+fG6wM6kidkBtpK0T358NHBdWeUKIbRWR7QJFlwAvJhegiCwLfB7SfOA24CbgYvzvnuB90u6G9gA+G7J5QohtILb/3S47MnXL4NVId32Q8CU/Pgc4Jza50gCWG77bWWWJYTQeoNhsHRcMRJCqEy6djhupdWnYm0xhDD0OGqCIYROFjdQCCF0LDvaBEMIHS5Oh0MIHUysiI6RanjkunTttGVl6Y99qNprk6ueGB1g9iNzK03/PQv36f+gJhww7qFK0z/jQ3tWmj7A4n03qTT9ja//Z3WJq/kaXAyRCSF0Nqd2wXYWQTCEUKnoHQ4hdCwTHSMhhI7W2uuCGxFBMIRQqe7uCIIhhA5lx+lwCKHDxelwCKGjtfsQmUqGckt6vSRL2qmwbS9J10i6V9Jtkn5QnFVO0mWSbqiiPCGE1rHV0NIqVV3PMoN0i/wZAJI2AX4KfNr2jrb3AK4Axub944E9gXGStqmoTCGEAdZRc4z0kDSGNFvccayafe79wI9sX99znO2f2X40r74R+AVwIWtO4h5CGKzcmTXBI4ErbP8ZeFzSnqSbpq4x+1zBDNK8JBfkx71abd7hrsVlljmEUBU3uDRA0qG5Se1+ScfXOeatku6SNF/Sj/tLs4qOkRnAafnxhfQR1GDlqfL2wHW2LalL0hTbd9Yea/tM4EyAF4zdos2bW0MIUN4QGUnDgW8DBwMLSXOXX277rsIx2wOfAV5m+0lJL+wv3VKDoKQNgVcBu0oyMJwU439EavP7eS9PeytphrkH86RLLyAFzs+VWbYQQmuU2Du8F3C/7QcAJF1IOvO8q3DMu4Fv234y5e1+b7NT9unwm4FzbU+yPdn2lsCDwFXAMZL27jlQ0htzLXAGcGg+fjIpWEa7YAhDQM+1ww22CW7c09yVl5k1yW0BLCisL8zbinYAdpD0B0k3SDq0vzKWfTo8A/h6zbaLSUFtOnBKrp52A9cA9wCTgJVDY2w/KOlfkva2/aeSyxdCGEgGN37Z3CLb05rMcQSpeW1/YCJwjaRdbT/V1xNKY/uAXradXlh9eS9Pq43k2J5aZrlCCC1U3unww0DxTsoT87aihcCfbHeRmtj+TAqKN9VLtL3vex1CGOQaOxVusPPkJmB7SVtLWpd0dnl5zTGXkWqBSNqYdHr8QF+JRhAMIVSrpCEytpcDHwBmA3cDF9meL+nLkl6XD5tNGpp3F3A18Enbj/eVblw7HEKoTsl3kbE9C5hVs+3zhccGPpaXhkQQDCFUq81H9EYQDCFUK26lFULoaFETrIhh2PLuypJfMbLat0brrVdp+lD9vMBnTLy+/4OasOufjq40/Y1GVF9DGbGk2gjQPXZkZWl7eAnvj4maYAihs7X7TVUjCIYQqhVBMITQyRSzzYUQOtZa3CuwVSIIhhAqpOgYCSF0uKgJhhA6WgTBEEJHa/MgWNpdZPI8w+cV1kdIekzSLwvbXi/pDkl3S5on6fVl5R9CaENOvcONLK1SZk1wMTBF0ijbz5EmQ1l5w0NJLwZOAQ7Od4/eGviNpAds31FiOUII7aRTaoLZLOA1+XHPNJo9PgF81faDkG6jD5wEfBJA0ofyNHl35AlUQgihcmUHwQuB6ZJGArsBxTlCdmHNuYdvztsBjgf2sL0b8N7eEi/OO9wV8w6HMCjIjS2tUmoQzKe1k0m1wFl9H72GO4DzJb0NWF4n/TNtT7M9bZ111m+qrCGEAWI1trRIFbfXv5zU9ndBzfa7SNNpFu0JzM+PX0OaWHkqaVLl6LkOYbBr9Nb6LawJVhFozgKesj1P0v6F7acAP5X0O9sPSZoMfBZ4s6RhwJa2r5Z0HWkClTFA3WnyQgiDg6q7410pSg+CthcCp/eyfa6kTwO/kLQO0AV8Km9fBzhP0jhAwOl9zRMaQhhE2rx3uLQgaHtML9vmAHMK65cAl/RyXBewX1llCSG0kU4JgiGEUKvVPb+NiCAYQqhW3EUmhNDJOq5jJIQQVhOnwyGEjhVtgiGEjhdBsBoeLpaNW7ey9Ictr/Y/52XLKk0f4IBxD1WaftXzAs/b+8eVpv/yc95TafoAizet4qKsVcbPXVpZ2uou6TsQQTCE0MnidDiE0NkiCIYQOlZ0jIQQOl4EwRBCR4sgGELoVCJOh0MIna7Ng2C1g5hCCJ3N6drhRpZGSDpU0r2S7pd0fB/HvSlPAzytvzQHPAhKWiFprqQ7Jf1U0ui8/ZmBLksIYQCUdHt9ScNJU3AcBuwMzJC0cy/HjQU+zOoTvdXViprgc7Z3tz0FWEadmeVCCENDibPN7QXcb/sB28tIs1se2ctxXwG+DixpJNFWnw5fC2zX4jKEEKrUeE1w454pdfMysyalLYAFhfWFedtKkqaS5iv6VaPFa1nHSJ5N7jDgirV4zkxgJsB6I8dXVLIQQmnWbia5Rbb7bcOrJ0/Ydipw7No8rxVBcJSkufnxtcAPG32i7TOBMwHGjpvY5n1OIQQo9aaqDwNbFtYn5m09xgJTgDmSADYFLpf0Ots310u0FUHwOdu7tyDfEEILlDhO8CZge0lbk4LfdGDlrYxs/wvYeGW+0hzgE30FQGh9m2AIYagrqXfY9nLgA8Bs4G7gItvzJX1Z0uueb/HaabD0aEkLC+un2j61ZaUJITRv7doE+0/OngXMqtn2+TrH7t9ImgMeBHubnzhvj1ppCEOM8tLO2qkmGEIYitq8CzOCYAihUjHlZgihs0VNMITQseLO0iGEjhdBMITQyaImWJHuEWLJhtUVf8+P3VZZ2gAPHL11pekDnPGhPStNf6MR1Q5+qHpe4Gu/fUal6QO87I43Vpr+sturu4beDw8vKaFykqnKoA2CIYRBwNE7HELodFETDCF0qphoKYQQIgiGEDqZ3N5RMIJgCKE60TESQuh47V0RjCAYQqhWu3eMlH4PP0mbSrpQ0l8k3SJplqQd8r6PSFoiaVzh+I0kXS3pGUnfKrs8IYQWK+nO0lUpNQgqzW5yKTDH9ra29wQ+A2ySD5lBmiegOIx+CfB/gE+UWZYQQhtocM7hVtYWy64JHgB02f5ezwbbt9u+VtK2wBjgBFIw7Nm/2PZ1NDhRcghhkGnzmmDZbYJTgFvq7JtOmjH+WmBHSZvYfrTk/EMIbUSAutu7UXAg5/WYAVxouxu4GHjL2iYgaWbP7PTLlywuvYAhhPK1++lw2TXB+cCbazdK2hXYHvhNnhR5XeBBYK06QoqTr6+/0Zbt/fMSQmj5qW4jyq4J/g5YT9LMng2SdgNOB75oe3JeNgc2lzSp5PxDCG1G3Y0trVJqELRt4A3AQXmIzHzgJGB/Uq9x0aWkdkIkPQScChwraaGkncssVwihhTqsYwTbjwBvbeC4jxUeTy67HCGENuD27xiJK0ZCCJVq9ytGIgiGEKoVQTCE0KnipqohhM5mp6WNRRAMIVQqaoIhhI4WN1UNIXQuAzFEphoeAc9NqO7S59m/nVpZ2gDbPXJHpekDLN53k/4PasKIJdV+uBdvWu2l7VVPjA7wh90uqTT9A0a+u7K0rbISKimdigzaIBhCGByiTTCE0NnavHd4IG+lFULoQGXeSkvSoZLulXS/pON72f8xSXdJukPSbxu5SUsEwRBCZZSvHW5k6TctaTjwbeAwYGdgRi83W7kNmGZ7N+BnwP/0l24EwRBCtbobXPq3F3C/7QdsLyPdqf7I4gG2r7b9bF69AZjYX6LRJhhCqJQabxPcWNLNhfUz842Ue2wBLCisLwT27iO944Bf95dpBMEQQnXW7l6Bi2xPKyNbSW8DpgGv7O/Yyk+HJVnSeYX1EZIek/TLvH5sXr9N0n2SZkvat+pyhRAGglddP9zf0r+HgS0L6xPzttVIOgj4HPA620v7S3Qg2gQXA1MkjcrrB7NmwX9iew/b2wNfAy6R9KIBKFsIoWJldYyQ5izfXtLWktYl3Zn+8tXykvYAziAFwH82kuhAdYzMAl6TH88ALqh3oO2rSZMpzax3TAhhkHB5c4zYXg58AJgN3A1cZHu+pC9Lel0+7GTS/OY/lTRX0uV1kltpoNoELwQ+n0+BdwPOAl7ex/G3Au8ZiIKFECpW4mBp27NIlarits8XHh+0tmkOSBC0fYekyaRa4Ky+jwbSvRjX3JhmsZsJsM7YDcoqXgihSu19wciA9g5fDpxCmnluo36O3YNU3V1Ncd7h0ZvEvMMhDAZrMUSmJQYyCJ4FPGV7nqT96x0k6ZWk2t4BA1WwEEKFIggmtheSJmHvzVGS9gNGAw8Cb7K9Rk0whDC4yEYrOjwI2h7Ty7Y5wJz8+Gzg7KrLEUJokagJhhA6WgTBEELHMo3eHKFlIgiGECoVvcMhhM4WQTCE0LFs6G7v8+EIgiGEarV3DIwgGEKoVrQJVmRYF6z/9xWVpT/h1n5vQ9aUYRP6u3KweRtf39CdhJ637rEjK01//Nxq/wfLbh9fafpQ7bzAAFef9f3K0t7rkEXlJBRBMITQsQw0dq/AlokgGEKoUHSMhBA6XZwOhxA6VpwOhxA6m8FxOhxC6GRxOhxC6FhxOhxC6Hht3jvc55Sbkq6WdEjNto9I+rWkO/P6aEnnS5on6U5J10la40aqIYROVOrk65Xob97hC0gTHBdNB04qrH8YeNT2rranAMcBXeUVMYQwaJlUE2xkaUYkEygAAAvCSURBVJH+guDPgNfk2d7J02ZuDiwoHLMZ8HDPiu17bS/Nx18m6RZJ8/N0meTtz0g6OW+/StJekuZIeqAwiXIIYSgYzDVB208ANwKH5U3TgYtYfSbRs4BPS7pe0omSti/se6ftPYFpwIck9Vwwuz7wO9u7AE8DJwIHA28AvlyvPJJmSrpZ0s1dS59p+EWGEFpoMAfBrHhKPD2vr2R7LrANcDKwIXCTpBfl3R+SdDtwA7Al0BMglwFX5MfzgN/b7sqPJ9criO0zbU+zPW2d9aLZMYT259Q73MjSIo30Dv8c+IakqcBo27fk0+KVbD8DXAJcIqkbOFzSJsBBwD62n5U0B+i57UiXvTL0dwNLczrdkqLHOoShwuAV1d3tqQz91gRzgLuadNp7Qe1+SS+TtEF+vC6wM/BXYBzwZA6AOwEvLbPgIYRBos1PhxutdV0AXMqaPcUA2wLflSRSUP0VcDGwLvBeSXcD95JOiUMInWSo3F7f9mWACusPAVPy43OAc3p52lJWdajUpjem8PiL9faFEIaAuGwuhNDJPBRqgiGE8LzYsCKCYAihk8WttEIIncqA4y4yIYSO5bipagihw7V7TVBu8+7reiQ9RhqU3aiNgZImUh2S6Q9EHpH+4Ep/ku0JzWQo6YqcbyMW2T60mfyej0EbBNeWpJttT4v0W5dHpD+00x+sGrmBQgghDFkRBEMIHa2TguCZkX7L84j0h3b6g1LHtAmGEEJvOqkmGEIIa4ggGELoaBEEQwgdLYJgWI2kl1SY9jBJL6gq/UI+H6k6jzB0DMmOEUm/YPUZ8VZju6lpPSWNBI4CngR+AXwKeDnwF+Artpse9Z/zeB+wH+m1XAd81/aSZtPuJa+dgRl5earMAbWSfgy8F1gB3AS8ADjN9sll5dFLnn+zvVVFae9k+56S0toRmAnslDfdDXzf9r0lpf9a4CvAJNIlsgJsu/IfosFkqAbBV+aHbwQ2Bc7L6zNIE8V/tMn0LyJNML8+sAFwJykY7gfsbvu1zaRfyONpVpX9aGC87bc0m3ZOfzKrAl8X6YsyLd81vDSS5treXdJ/AFOB44FbbO9WZj41eS6wvWVFaZcSYCXtQ5qc7AzgNlKA2gN4N/BG201PRyHpftJ3YJ6H4he9JEPyBgq2fw8g6X9rajW/kHRzCVnsbHtKnhlvoe2eoHtFnmK0DFNs71xYv1rSXWUkLOl6Uo3sQuBNtu+T9GDZATBbR9I6wOuBb9nuklT1F7Kp9CWdXm8XML6ZtAs+D8ywPaew7TJJvwO+QJ2pKdbSAuDOCIB9G5JBsGB9SdvYfgBA0tak2luzlgHYXi7pkZp9Zc0veKukl/bUCCTtDZQRwAEeBbYANgEmAPfRZODowxnAQ8DtwDWSJgH/bjZRSU+Tytwz901P+QWMajL5dwAfJ08FW2NGk2n32LYmAALpB1xSWYOaPwXMkvR7Cq/F9qklpT8kDMnT4R6SDiWNkn+A9OWYBMy0fWWT6f6TVIsSqW3wwp5dwFttb9JM+jmPu4Edgb/lTVuRZu1bTmrXaep0UtI40qnSDGB7Ug3nENs3NpNuA/kKGG57eZX5NCPXxk6w/cde9j1oe+sS8rjF9p519t1qe2oJeVwJPAPMI83vDYDtLzWb9lAypIMggKT1WNXwfI/t3n7d1zbNY/rab/tHJeQxqZ881uY2Yv3l9ULgraSAuFVZ7WmSpgCfBHbJm+YDp9ieV0LaI0kdLtsBdwBnlRVYJW0ILLH9bBnp1cmj54d0jV2U90N6p+0pzaYz1A3JIJiHYWxi+768/hZWnSLNtv1oywrXD0mjgS7bXXl9R+Bw4K+2Lykpj5HAWNuP1Wx/IbCR7btLyONI4BTgJFadxk8DPgN8wvbPm0z/J6QOnWtJ7Wd/tf3hZtJsIM8tgell9GwP0A/p/wBXNXvmM9QN1SB4JvBH22fn9fuAXwOjgeW239tk+rVDcEy6WeXVts/r/VkNp30NcFzurNgOuBE4H9gZuNH2Z5pJP+dxJnBFbVCV9Abg1bb/q4Q8bgeOrO1syb3SP7f94ibTn2d71/x4BOm9afoUspd8JgBvIdWSNwcutf2JsvOpQm43XZ/UHthFDJHp1VDtGHkJ8J7C+jO2PwQg6boS0j+ll20bAm+TNMX28U2kvUFPDRY4BrjA9gclrQvcQqpJNWtP2zNrN9q+VNKJJaQPMKK33mbbD+Xe4mZ1FdJcnpoayyFpLKm99GhgB9JQlq1tTywxj8v72t/sWNacxthm0+gEQzUIjqgZFvD2wuOmhzj0DMGplT/Yt5DGwj3v5AuPXwWcnPNcJqmsGWtG97GvrKuIlkvayvbfihtzW2cZbXcvltTTyyxgVF4vo7bzT1IN/ATgOtvOteQy7UMawnIB8CdW9XKXStIGpI6vkT3bbF9TRV6D1VANgt2SNrX9DwDbdwJI2oJCL1nZbK8ooUZyh6RTgIdJjf5XAkgqa3wawD8l7VXbEyxpL+CxOs9ZW18ArpL0VdIPA6Q2weOBTzebuO3hzabRh88A04HvABfk9seybQocTDrNPhr4FanWP7+sDCS9C/gwMBGYC7wUuJ704xp62B5yC/A20iVarwDG5uWVedt/lpD+hr0s2wJfAs5vMu1RpEBxGvDiwvZ9gbeX9P7sRRq790XgiLx8CXgQ2LvE/8OLgXNIQfAW4Nzia2r3BdgG+CxpiMkS0ri7HSrIZz3gWNIP0AdKTHceqQY4N6/vBFzS6ve13ZYh2TECK8cIfpY0PMOk4Rlfs/3rEtJ+kNUH6nYDjwNzgBNtNz0YuGqSNiFdmzyFVe/Pb4GjbL+/lWVrtdwhtYntPxS27Ur6YXqlS6qF5uFbryHVBicDl5OG+jxcUvo32X6JpLmkH7elkubb3qXfJ3eQoXo6jO0rgCtqt0v6iO1vNpn8UcAC23/PaR4DvIn0q9vUeyppHn3f/KGUa26dhgl9QdJU0pfwGFLN+eIy0h+Ihv8KfZOaDijb8/Ldab5aRgaSziH9AM0CvuTcZFOyhbkZ5TLgN5KeZO2mqe0IQ7YmWE8ZF8BLuhU4yPYTkl5BGvT6QWB34EW239xE2rNIX7SF9BIMXcIgaUk7sOrmCYuAn5DG7vU5QHst83iMPhr+XadzqR301KDq7Fs5NKfJPLqBxXm1+H+uZBhLvqnIONLQqGVlpj3YDdmaYB/K6IUbbvuJ/Pgo4EzbFwMX51OPZswm9QhvBlxEaiy/rck0a91DGmT8Wtv3A0hq6s46vai84b9CfXVCNXtdMgC2K7uXZ50B91NJA8ojANboxJuqllH1HZ4H6AIcCPyusK+pHxbbp9neh9SR8zhwlqR7JH0h1+DK8Ebg76Q703xf0oGUPETD9grbV9g+htQreT8wR9IHysynIjdLenftxtzbeksvx7ebK0htjD3tm9eTOnneL+lrLSxXWxqSp8OFO4yssQsYZbvZdrvPkS5lW0S6scFU284fuB/Zflkz6feS3x7AWcBuZTXK53TXB44k1dZeRerJvdQlXWZVdcN/VXKn0aWkuwUVh/esC7zBeehVu6q5muYrwIa2398z4L6M0/mhZEgGwYEg6aWkU9YrbS/O23YAxti+tYT0R5CuiZ1Oqm3OIZ1ONnXNbR/5bUC6POwo2weWkF6x4f/Cihr+KyXpANJrAJhv+3d9Hd8uJN3R04Em6Q/AybYvy+u3u8lLFoeaCIJtRlJPO9rhpKsWLiRda7u4zye2mYFu+A+rSDoP+AdpwP3xpEv+ns09xb+PILi6CIJtRuledj8GLrb9ZKvLEwYfSaNIV4psRmp+uD1v35d0M9dzW1m+dhNBMITQ0TpxiEwIQ9pADbgfKiIIhjD0LKCPAfdhdREEQxh6BmLA/ZARbYIhDFH53o3T8zKKdAnjBbb/3NKCtZkIgiF0gKoG3A8FnXjZXAgdQdIISUdIOp80x869pEsmQ0HUBEMYYobKgPuBEkEwhCEmBtyvnQiCIYSOFm2CIYSOFkEwhNDRIgiGEDpaBMEQQkeLIBhC6Gj/H3Zh0ZzXZu9IAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 432x288 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 321
        },
        "id": "M4nPx2rLpFMv",
        "outputId": "ef620978-ec04-4b9c-f1eb-aabcaac1af5d"
      },
      "source": [
        "#@title Plotting Error \n",
        "positions = (0,1, 2, 3,4,5,6,7,8)\n",
        "plt.xticks(positions, roi_labels,rotation = 'vertical')\n",
        "plt.yticks(positions, roi_labels)\n",
        "plt.imshow(theta_array)\n",
        "plt.colorbar()\n",
        "plt.title('Theta Coherence Incorrect Trials Pre Stim')"
      ],
      "execution_count": 48,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "Text(0.5, 1.0, 'Theta Coherence Incorrect Trials Pre Stim')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 48
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUEAAAEeCAYAAAAU1qcTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deZxcVZ338c83GwESEyCRHcKOELYQ2dQRRhwWRVwhYXDAcYyOCzi4AOoji444I4PiM4yCDg8gCKIsgoZFlMgiyBoSwiIISAKKhEUgENJJf58/zulwU+nqrlD3dlV3/d6v13113aV+99TSvzr3nnvPkW1CCKFTDWt1AUIIoZUiCYYQOlokwRBCR4skGELoaJEEQwgdLZJgCKGjtW0SlHSipPNbXY7XQ9IsSf/S6nKEaki6StIRDW77mKR9qy5TVSRtIuklScNbXZaqtCwJ5je2Z+qW9Eph/h+bjN10ApW0vqT/lfRnSS9KekDSSZLWbCbuYCRpb0kLWl2OZkiypC3rrPtS4bu3WNKywvy82u1tH2D73ArLunf+n3gpf/celPSRkmKPkvRfkhbk+I9J+k5h/QpJ2/bjtsfYXlbG/ttRy5JgfmPH2B4DPA4cVFh2QavKBSBpbeAWYHVgT9tjgXcC44EtBrAcktS2tfUqSRrRy7JKaiO2v1H4Ln4CuKXwXdy+sP+B/DyezOV5A3As8ANJ29Vu1Nv71I/jganAbsBYYG/gruaKOri1+z/YKEnn5V/DeZKm9qyQtIGkSyQ9LelRSUfl5fsDXwIOzb909+TlH5F0f471iKSP97HfY4AXgcNtPwZge77to23PyfH2knS7pL/lv3vVxNhU0s15f9dKmlAo+x6SfifpeUn3SNq7sG6WpH+XdDPwMrC5pG0l/UrSs7lWcEhh+3MknSHpl3lfv5e0RWH99oXnPiXpS3n5MEnHSfqjpGckXZyTf79yGb/Wx+t7a+H1zZd0ZF4+Ln+eT0v6k6Sv9CQVSUfmeN+W9AxwYn5t35M0U9IiYJ96n3uOMTzX6v6Yy3WnpI0l3ZA3uSd/Jw5t5HX28XksP90haQtJv8nv4UJJF0gaXyfWbpLukPRC/ixO62//Ti4HngO2q/M+rSbpVEmP57jfl7R6nZBvBi6z/WSO/Zjt83L5fgRsAlyZ36cvSpqkVIseUXg/vp4/35ckXSlpnfy6X8j/C5MafX/bgu2WT8BjwL41y04EFgMHAsOBU4Bb87phwJ3AV4FRwObAI8B+heeeXxPvXaRanIC3k77QU+qU51bgpD7KuzbpS/lhYAQwPc+vk9fPAv4IbE2qTc4CvpnXbQg8k1/XMFIN8xlgYuG5jwPb59jjgPnAR/L8LsBCYLu8/Tn5+bvl9RcAF+V1Y4E/A58DRuf53fO6o/Pr3AhYDTgTuLDO690bWFCY7+v1bUr6AZkOjATWAXbO684Dfp7LMQn4A/DRvO5IYCnwmfw6Vs+v7W/AW/J7tUY/n/sXgLnANvlz3qnwmRjYsoHv4pHATTWvtfh5jMzL/iWv3zJ/hqsBE4EbgO/09t0mHV18OD8eA+zR3/udX/f7gK78unp7n74NXEH6Xo4FrgROqRP7K/n1fBLYAVBf/4v5czIwovB+PEz6XxoH3Jc/x31zec4D/l+rc8oq5Z9WF6C3Nz4vOxG4rjC/HfBKfrw78HjN9sf3vPn0kgR72eflwNF11j0EfKKP534YuK1m2S3AkYUvylcK6z4JXJ0fHwv8qOa51wBHFJ57cmHdocCNNdufCZyQH58D/LCw7kDggfx4OnB3nddwP/COwvz6+R9tRC/bLv+nbOD1HU+qadTGGA4sISfvvOzjwKz8+MhePtNzgPMK8/197g8CB9d5vc0kwZNrtplFToK9PP+9xfecFZPgDcBJwIR+yrA30A08DzwLzAam9fY+kZL9ImCLwrI9gUfrxB4OfAq4GXgVeLLnu1db3jw/iZWT4JcL6/8LuKowfxAwu7/3uZ2mVT2fMND+Unj8MjA6V8s3BTaQ9Hxh/XDgxnqBJB0AnECqvfTUKubW2fwZUlKoZwPgTzXL/kSq5dUr+5j8eFPgQ5IOKqwfCVxfmJ9feLwpsHvNax0B/KiBfW1MqrH1ZlPgMkndhWXLgHWBJ+o8p2hV9zmB9DqL71vtezafldW+F3197n293mb0Vi4AJK0LnA68jVQLG0Y6KujNR4GTgQckPUo62vhFnW2ftL1RA+WZSK4hS1peLNL7shKnBo4zgDPyIfM/A2dLus32/XX2V+upwuNXepkfwyDS7ucE65lP+qUbX5jG2j4wr1+haxxJqwGXAKcC69oeD8wkfVl6cx3wPtU/Cf4k6R+yaBMaSx7zSTXBYtnXtP3Nwjau2f63NduPsf2vDe5r8z7WHVATd7TtRl5Df/vsrfFoIammWXzfat+z3ro0qn0v+vrc6+27WX11tfSNvH4H228ADqfO98r2Q7anA28E/gP4mV7f1QbF8iwkJZ7tC+/JOKdGlb6D2K/YPoN8vrGX2B1hsCbB24AXJR0rafV8QnyypDfn9U8BkwpJbBTpnM3TwNJcK/yHPuKfRmqVO1fSpgCSNpR0mqQdSQl0a0mHSRqRT7RvB9T7VS86HzhI0n653KOVLomo96v/i7yvD0samac3S3pTA/v6BbC+pM/mk+djJe2e130f+PfC65so6eAGYvbnAmBfSYfk92YdSTvnGsjFeZ9j836PIb0fjervc/8h8DVJWynZUdI6ed1T1P9BaMZY4CXgb5I2JJ2X7JWkwyVNtN1zqAvpsPd1y7F+AHxb0hvzfjaUtF+dMnw2f99Wz5/PEfk13J03qep9aluDMgnmf6h3AzsDj5J+DX9IOlEL8NP89xlJd9l+ETiK9E/4HHAY6URyvfjPAnuRai6/l/Qi8GvSSfqHbT+T9/850qHzF4F3217YQNnnAweTWrCfJtVevkCdzyKX/R+AaaQa6F9ItYjVGtjXi6ST9gfl5z0E7JNXn056D67Nr+9W0jm3pth+nHRe8nO8dj5rp7z6M6TzV48ANwE/Bs5ehdj9fe6nkT7ja4EXgP8lNRxAOk98rlKL9SGU5yRgCum78Uvg0j623R+YJ+kl0vs/zfYrJZThWFJjxa2SXiAdyWxTZ9uXSefx/kJ6/z4FfMD2I3n9KcBX8vv0+RLK1vaUT2aGEEJHGpQ1wRBCKEskwRDCoCDpbEl/lXRvnfWS9F1JD0uaI2lKI3EjCYYQBotzSOdV6zkA2CpPM4DvNRI0kmAIYVCwfQOpsa2eg0kX19v2rcB4SX1d7wtEEgwhDB0bsuKF5AtY8WL8XrX7HSN1TVh7uCdtPLKy+A8smFhZbIDhLyyuND6AR4+qeAcVhx9e71r2cnSPqDY+gCv+DxvWVV3sVxc9S9eri5p6k/bbZ00/82xjvXDdOefVeaT+AnqcZfusZvbfiEGbBCdtPJLbrtm4svh7HfOJymIDjL/uoUrjA3RtW937AzBsaVPX+fZrybhqk/jitav/+r8ysdqDrTX/XF03f3N/dXrTMRY+u4zfX1PvPoAVjVz/j4ttT+1/y7qeIN062WMjGriLKw6HQwgVMsvc3dBUgiuAf8qtxHsAf7P95/6eNGhrgiGE9megu6TzJpIuJPWwM0Gpp/MTSJ1yYPv7pNtZDyTdPfMyqfu5fkUSDCFUxpiuknrmz51P9LXepNsAV0kkwRBCpcqqCVYlkmAIoTIGlkUSDCF0snavCZbeOpy7Cept+eH5fr55SoML/VC9DEiT+zprpF++EEKbM7DMbmhqlQGpCSqNAPdvpJ6Mn1AaOvEIUlfuz/f55BDCoGVMV5vXBAfqcPjLwOd7um7PnWMu70wzJ8nvkJq1bxqgMoUQqmZY1t45cMAult6eOgM8SxpN6h78IGBXYL0BKlMIoWLpOsHGplYZ8DtGJO0gabbSANmHAtuSBs95KF/nU3fMCUkzlAavvuPpZ6q7XSiEUBaxrMGpVQYqCc4jjcOA7bm2dwau4rXxHxpi+yzbU21PnbhOryMKhhDaiIFuNza1ykAlwVOAU2tGVOtJgA+QRobrGSqxz6vCQwiDS7vXBKtoGFkj39fX4zTbp0maCFyVW4afB+4FrrG9WNIM4JeSXiYNpD22gnKFEAaYgS63dz8tpSdBu/dXbPtc4Nw6664mnRsMIQwh6Y6R1tXyGhF3jIQQKmPEsjbvsS+SYAihUt2OmmAIoUPF4XAIoaMZ0VX1QCtNau/ShRAGvagJhhA6li2WddolMgPlgQUTKx0R7nenfb+y2AB7fLHa0ewAxj5W7bCey0ZX+/UZtrTa2wh2PebuSuMDXPPrKZXGn3jXq5XFHtZVzh293VETDCF0qtQwEjXBEELHisPhEEIHS7fNtXdnJ5EEQwiViTtGQggdrzsOh0MInSoaRkIIHc2IZXHvcAihk3W3eU2winGHv5zHFp6TxxLZXdJjkiYUtlk+trCkIyU9nbedJ+lnktYou1whhIFniy4Pb2hqlVKToKQ9gXcDU2zvCOwLzG/gqT+xvbPt7YElwKFlliuE0Bpp8PVhDU2tUvbh8PrAQtuvAtheCCA1dk5A0ghgTeC5kssVQmiRdm8YKbt01wIbS/qDpP+R9PYGn3eopNnAE8DawJUllyuE0AJGdLuxqVVKTYK2XyINoD4DeBr4iaQjSbXilTYvPP5JHoZzPWAu8IXe4hfHHe5avKjMoocQKrKMYQ1NrVL6nm0vsz3L9gnAp4EPAM8AaxU2WxtY2MtzTaoF/l2d2MvHHR45es2yix5CKFnqVLWzGka2kbRVYdHOwJ+AWcCH8zbDgcOB6+uEeSvwxzLLFUJojTT4+rCGplYpu2FkDPB/JY0HlgIPkw6Nu4DvSboHEHA1cH7heYdKeispKS8Ajiy5XCGEFumonqVt3wnsVWf1YXWecw5wTpnlCCG0B1ul1vIk7Q+cDgwHfmj7mzXrNyGNbz4+b3Oc7Zl9xYw7RkIIlSrrGsB8Ku0M4J2kI8bbJV1h+77CZl8BLrb9PUnbATOBSX3Fbe8LeEIIg5pJ3es3MjVgN+Bh24/YXgJcBBzcyy7fkB+PA57sL2jUBEMIlTGiq7u0lt8NWfEOtAXA7jXbnAhcK+kzpBsv9u0vaNQEQwiVWoXrBCf0XAecpxmvY3fTgXNsbwQcCPxIUp95LmqCIYTK9Nwx0qCFtqf2sf4JYOPC/EZ5WdFHgf0BbN8iaTQwAfhrvaBREwwhVKqbYQ1NDbgd2ErSZpJGAdOAK2q2eRx4B4CkNwGjSXev1TVoa4LDX1jM+Oseqix+1eMC3/qf1Y5rDLDfBjtXGl+rrVZpfC9ZUmn8Rw7brNL4AFs+OafS+MMmrlNZbC1Z1nQMm9I6VbW9VNKngWtIl7+cbXuepJOBO2xfAXwO+IGkfyM1khyZ70Sra9AmwRDC4FBm5wj5mr+ZNcu+Wnh8H/CWVYkZSTCEUJmee4fbWSTBEEJl0r3DHXTbXAghrKjc2+aqEEkwhFCpBu8GaZlIgiGEypTZOlyVSIIhhMoYsbS82+YqEUkwhFCpOBwOIXSswdA6XEqzjaTrJe1Xs+yzkq6SdG+eX0PSBZLmSrpX0k2SxpSx/xBC+2r37vXL2vOFpPv4iqYBpxTmjwaesr2D7cmkG527Stp/CKEdNTjc5lAYcvNnwLvyTc1ImgRswIp9f61PoccH2w/aflXSJEkP5Fri/ZJ+JmmNksoVQmihkjtVrUQpSdD2s8BtwAF50TTgYlYcW/hs4FhJt0j6es2odNsA/2P7TcALwCd7209x3OEl3YvLKHoIoUIGlnYPa2hqlTL3XDwknpbnl7M9G9gc+BZp3OHbc1c3APNt35wfn08adnMlxXGHRw0bXWLRQwhV6ZTDYYCfA++QNAVYI488twLbL9m+1PYnScnuwJ5VtZuWWK4QQov0dKraEUnQ9kukAdXPpqYWCCDpLZLWyo9HAduRBmYH2ETSnvnxYcBNZZUrhNBaHXFOsOBCYCd6SYLAFsBvJc0F7gbuAC7J6x4EPiXpfmAt4HsllyuE0Apu/8PhsgdfvxxeS+m2HwMm58fnAefVPkcSwFLbh5dZlhBC6w2Gi6XjjpEQQmXSvcPRlVafirXFEMLQ46gJhhA6WXSgEELoWHacEwwhdLg4HA4hdDCxLBpGquHRo+jaduPK4o99rNp7k6seGB3gmidnVxr/4wv27H+jJuwz7rFK45951K6VxgdYtNe6lcafcMtfqwuu5mtwcYlMCKGzOZ0XbGeRBEMIlYrW4RBCxzLRMBJC6GitvS+4EZEEQwiV6u6OJBhC6FB2HA6HEDpcHA6HEDpau18iU8ml3JLeK8mSti0s203SDZIelHS3pB8WR5WTdLmkW6soTwihdWw1NLVKVfezTCd1kT8dQNK6wE+BY21vY3sX4GpgbF4/HtgVGCdp84rKFEIYYB01xkgPSWNIo8V9lNdGn/sUcK7tW3q2s/0z20/l2fcDVwIXsfIg7iGEwcqdWRM8GLja9h+AZyTtSuo0daXR5wqmk8YluTA/7tUK4w53LSqzzCGEqrjBqQGS9s+n1B6WdFydbQ6RdJ+keZJ+3F/MKhpGpgOn58cX0UdSg+WHylsBN9m2pC5Jk23fW7ut7bOAswDeMHbDNj/dGkKA8i6RkTQcOAN4J7CANHb5FbbvK2yzFXA88Bbbz0l6Y39xS02CktYG/h7YQZKB4aQcfy7pnN/Pe3naIaQR5h7Ngy69gZQ4v1xm2UIIrVFi6/BuwMO2HwGQdBHpyPO+wjYfA86w/Vzat/vtZqfsw+EPAj+yvantSbY3Bh4FrgOOkLR7z4aS3p9rgdOB/fP2k0jJMs4LhjAE9Nw73OA5wQk9p7vyNKMm3IbA/ML8grysaGtga0k3S7pV0v79lbHsw+HpwH/ULLuElNSmAafm6mk3cAPwALApsPzSGNuPSvqbpN1t/77k8oUQBpLBjd82t9D21Cb3OIJ0em1vYCPgBkk72H6+ryeUxvY+vSz7bmH2bb08rTaTY3tKmeUKIbRQeYfDTwDFnpQ3ysuKFgC/t91FOsX2B1JSvL1e0Pbu9zqEMMg1dijcYOPJ7cBWkjaTNIp0dHlFzTaXk2qBSJpAOjx+pK+gkQRDCNUq6RIZ20uBTwPXAPcDF9ueJ+lkSe/Jm11DujTvPuB64Au2n+krbtw7HEKoTsm9yNieCcysWfbVwmMDx+SpIZEEQwjVavMreiMJhhCqFV1phRA6WtQEK2IYtrS7svDLRlf71mi11SqND9WPC3zmRrf0v1ETdvj9YZXGX2dE9TWUEYurzQDdY0dXFtvDS3h/TNQEQwidrd07VY0kGEKoViTBEEInU4w2F0LoWKvQV2CrRBIMIVRI0TASQuhwURMMIXS0SIIhhI7W5kmwtF5k8jjD5xfmR0h6WtIvCsveK2mOpPslzZX03rL2H0JoQ06tw41MrVJmTXARMFnS6rZfIQ2GsrzDQ0k7AacC78y9R28G/ErSI7bnlFiOEEI76ZSaYDYTeFd+3DOMZo/PA9+w/SikbvSBU4AvAEg6Kg+TNycPoBJCCJUrOwleBEyTNBrYESiOEbI9K489fEdeDnAcsIvtHYFP9Ba8OO5wV4w7HMKgIDc2tUqpSTAf1k4i1QJn9r31SuYAF0g6HFhaJ/5Ztqfanjpy5JpNlTWEMECsxqYWqaJ7/StI5/4urFl+H2k4zaJdgXn58btIAytPIQ2qHC3XIQx2jXat38KaYBWJ5mzgedtzJe1dWH4q8FNJv7H9mKRJwJeAD0oaBmxs+3pJN5EGUBkD1B0mL4QwOKi6Hu9KUXoStL0A+G4vy2dLOha4UtJIoAv4Yl4+Ejhf0jhAwHf7Gic0hDCItHnrcGlJ0PaYXpbNAmYV5i8FLu1luy7grWWVJYTQRjolCYYQQq1Wt/w2IpJgCKFa0YtMCKGTdVzDSAghrCAOh0MIHSvOCYYQOl4kwWp4uFgyblRl8YctrfaT85IllcYH2GfcY5XGr3pc4Lm7/7jS+G877+OVxgdYtF4VN2W9ZvzsVyuLre6S/gciCYYQOlkcDocQOlskwRBCx4qGkRBCx4skGELoaJEEQwidSsThcAih07V5Eqz2IqYQQmdzune4kakRkvaX9KCkhyUd18d2H8jDAE/tL+aAJ0FJyyTNlnSvpJ9KWiMvf2mgyxJCGAAlda8vaThpCI4DgO2A6ZK262W7scDRrDjQW12tqAm+Yntn25OBJdQZWS6EMDSUONrcbsDDth+xvYQ0uuXBvWz3NeA/gMWNBG314fCNwJYtLkMIoUqN1wQn9Aypm6cZNZE2BOYX5hfkZctJmkIar+iXjRavZQ0jeTS5A4CrV+E5M4AZAKuNHl9RyUIIpVm1keQW2u73HF49ecC204AjV+V5rUiCq0uanR/fCPxvo0+0fRZwFsDYcRu1eZtTCAFK7VT1CWDjwvxGeVmPscBkYJYkgPWAKyS9x/Yd9YK2Igm+YnvnFuw3hNACJV4neDuwlaTNSMlvGrC8KyPbfwMmLN+vNAv4fF8JEFp/TjCEMNSV1DpseynwaeAa4H7gYtvzJJ0s6T2vt3jtdLH0GpIWFOZPs31ay0oTQmjeqp0T7D+cPROYWbPsq3W23buRmAOeBHsbnzgvj1ppCEOM8tTO2qkmGEIYitq8CTOSYAihUjHkZgihs0VNMITQsaJn6RBCx4skGELoZFETrEj3CLF47eqKv+sxd1cWG+CRwzarND7AmUftWmn8dUZUe/FD1eMC33jGmZXGB3jLnPdXGn/JPdXdQ+8nhpcUqJwwVRm0STCEMAg4WodDCJ0uaoIhhE4VAy2FEEIkwRBCJ5PbOwtGEgwhVCcaRkIIHa+9K4KRBEMI1Wr3hpHS+/CTtJ6kiyT9UdKdkmZK2jqv+6ykxZLGFbZfR9L1kl6S9N9llyeE0GIl9SxdlVKToNLoJpcBs2xvYXtX4Hhg3bzJdNI4AcXL6BcD/wf4fJllCSG0gQbHHG5lbbHsmuA+QJft7/cssH2P7RslbQGMAb5CSoY96xfZvokGB0oOIQwybV4TLPuc4GTgzjrrppFGjL8R2EbSurafKnn/IYQ2IkDd7X1ScCDH9ZgOXGS7G7gE+NCqBpA0o2d0+qWLF5VewBBC+dr9cLjsmuA84IO1CyXtAGwF/CoPijwKeBRYpYaQ4uDra66zcXv/vIQQWn6o24iya4K/AVaTNKNngaQdge8CJ9qelKcNgA0kbVry/kMIbUbdjU2tUmoStG3gfcC++RKZecApwN6kVuOiy0jnCZH0GHAacKSkBZK2K7NcIYQW6rCGEWw/CRzSwHbHFB5PKrscIYQ24PZvGIk7RkIIlWr3O0YiCYYQqhVJMITQqaJT1RBCZ7PT1MYiCYYQKhU1wRBCR4tOVUMInctAXCJTDY+AVyZWd+vzNb+eUllsgC2fnFNpfIBFe63b/0ZNGLG42i/3ovWqvbW96oHRAW7e8dJK4+8z+mOVxbbKClRSnIoM2iQYQhgc4pxgCKGztXnr8EB2pRVC6EBldqUlaX9JD0p6WNJxvaw/RtJ9kuZI+nUjnbREEgwhVEb53uFGpn5jScOBM4ADgO2A6b10tnI3MNX2jsDPgP/sL24kwRBCtbobnPq3G/Cw7UdsLyH1VH9wcQPb19t+Oc/eCmzUX9A4JxhCqJQaPyc4QdIdhfmzckfKPTYE5hfmFwC79xHvo8BV/e00kmAIoTqr1lfgQttTy9itpMOBqcDb+9u28sNhSZZ0fmF+hKSnJf0izx+Z5++W9JCkayTtVXW5QggDwa/dP9zf1L8ngI0L8xvlZSuQtC/wZeA9tl/tL+hAnBNcBEyWtHqefycrF/wntnexvRXwTeBSSW8agLKFECpWVsMIaczyrSRtJmkUqWf6K1bYl7QLcCYpAf61kaAD1TAyE3hXfjwduLDehravJw2mNKPeNiGEQcLljTFieynwaeAa4H7gYtvzJJ0s6T15s2+Rxjf/qaTZkq6oE265gToneBHw1XwIvCNwNvC2Pra/C/j4QBQshFCxEi+Wtj2TVKkqLvtq4fG+qxpzQJKg7TmSJpFqgTP73hpIfTGuvDCNYjcDYOTYtcoqXgihSu19w8iAtg5fAZxKGnlunX623YVU3V1BcdzhNdaNcYdDGAxW4RKZlhjIJHg28LztuZL2rreRpLeTanv7DFTBQggViiSY2F5AGoS9N4dKeiuwBvAo8AHbK9UEQwiDi2y0rMOToO0xvSybBczKj88Bzqm6HCGEFomaYAiho0USDCF0LNNo5wgtE0kwhFCpaB0OIXS2SIIhhI5lQ3d7Hw9HEgwhVKu9c2AkwRBCteKcYEWGdcGaf15WWfyJd/XbDVlThk3s787B5k24paGehF637rGjK40/fna1n8GSe8ZXGh+qHRcY4Pqzf1BZ7N32W1hOoEiCIYSOZaCxvgJbJpJgCKFC0TASQuh0cTgcQuhYcTgcQuhsBsfhcAihk8XhcAihY8XhcAih47V563CfQ25Kul7SfjXLPivpKkn35vk1JF0gaa6keyXdJGmljlRDCJ2o1MHXK9HfuMMXkgY4LpoGnFKYPxp4yvYOticDHwW6yitiCGHQMqkm2MjUIv0lwZ8B78qjvZOHzdwAmF/YZn3giZ4Z2w/afjVvf7mkOyXNy8Nlkpe/JOlbefl1knaTNEvSI4VBlEMIQ8Fgrgnafha4DTggL5oGXMyKI4meDRwr6RZJX5e0VWHdP9veFZgKHCWp54bZNYHf2N4eeBH4OvBO4H3AyfXKI2mGpDsk3dH16ksNv8gQQgsN5iSYFQ+Jp+X55WzPBjYHvgWsDdwu6U159VGS7gFuBTYGehLkEuDq/Hgu8FvbXfnxpHoFsX2W7am2p45cLU47htD+nFqHG5lapJHW4Z8D35Y0BVjD9p35sHg52y8BlwKXSuoGDpS0LrAvsKftlyXNAnq6Hemyl6f+buDVHKdbUrRYhzBUGLysut6eytBvTTAnuOtJh70X1q6X9BZJa+XHo4DtgD8B44DncgLcFtijzIKHEAaJNj8cbrTWdSFwGSu3FANsAXxPkkhJ9ZfAJcAo4BOS7gceJB0ShxA6yVDpXt/25YAK848Bk/Pj84Dzennaq7zWoFIbb0zh8Yn11k/QufUAAAo+SURBVIUQhoC4bS6E0Mk8FGqCIYTwutiwLJJgCKGTRVdaIYROZcDRi0wIoWM5OlUNIXS4dq8Jym3efF2PpKdJF2U3agJQ0kCqQzL+QOwj4g+u+JvantjMDiVdnffbiIW2929mf6/HoE2Cq0rSHbanRvzW7SPiD+34g1UjHSiEEMKQFUkwhNDROikJnhXxW76PiD+04w9KHXNOMIQQetNJNcEQQlhJJMEQQkeLJBhC6GiRBMMKJL25wtjDJL2hqviF/Xy26n2EoWNINoxIupIVR8Rbge2mhvWUNBo4FHgOuBL4IvA24I/A12w3fdV/3scngbeSXstNwPdsL242di/72g6Ynqfny7ygVtKPgU8Ay4DbgTcAp9v+Vln76GWfj9vepKLY29p+oKRY2wAzgG3zovuBH9h+sKT47wa+BmxKukVWgG1X/kM0mAzVJPj2/PD9wHrA+Xl+Ommg+H9rMv7FpAHm1wTWAu4lJcO3Ajvbfncz8Qv7eJHXyn4YMN72h5qNneNP4rXE10X6R5maew0vjaTZtneW9I/AFOA44E7bO5a5n5p9zre9cUWxS0mwkvYkDU52JnA3KUHtAnwMeL/tpoejkPQw6X9grofiP3pJhmQHCrZ/CyDpv2pqNVdKuqOEXWxne3IeGW+B7Z6ke3UeYrQMk21vV5i/XtJ9ZQSWdAupRnYR8AHbD0l6tOwEmI2UNBJ4L/DftrskVf0P2VR8Sd+ttwoY30zsgq8C023PKiy7XNJvgBOoMzTFKpoP3BsJsG9DMgkWrClpc9uPAEjajFR7a9YSANtLJT1Zs66s8QXvkrRHT41A0u5AGQkc4ClgQ2BdYCLwEE0mjj6cCTwG3APcIGlT4IVmg0p6kVTmnrFvesovYPUmw38E+Bx5KNga05uM3WOLmgQIpB9wSWVd1PxFYKak31J4LbZPKyn+kDAkD4d7SNqfdJX8I6R/jk2BGbavbTLuX0m1KJHODV7Uswo4xPa6zcTP+7gf2AZ4PC/ahDRq31LSeZ2mDicljSMdKk0HtiLVcPazfVszcRvYr4DhtpdWuZ9m5NrYV2z/rpd1j9rerIR93Gl71zrr7rI9pYR9XAu8BMwlje8NgO2Tmo09lAzpJAggaTVeO/H8gO3eft1XNeYRfa23fW4J+9i0n32sSjdi/e3rjcAhpIS4SVnn0yRNBr4AbJ8XzQNOtT23hNijSQ0uWwJzgLPLSqyS1gYW2365jHh19tHzQ7rSKsr7Ib3X9uRm4wx1QzIJ5ssw1rX9UJ7/EK8dIl1j+6mWFa4fktYAumx35fltgAOBP9m+tKR9jAbG2n66ZvkbgXVs31/CPg4GTgVO4bXD+KnA8cDnbf+8yfg/ITXo3Eg6f/Yn20c3E7OBfW4MTCujZXuAfkj/E7iu2SOfoW6oJsGzgN/ZPifPPwRcBawBLLX9iSbj116CY1JnldfbPr/3ZzUc+wbgo7mxYkvgNuACYDvgNtvHNxM/7+Ms4OrapCrpfcA/2P7XEvZxD3BwbWNLbpX+ue2dmow/1/YO+fEI0nvT9CFkL/uZCHyIVEveALjM9ufL3k8V8nnTNUnnA7uIS2R6NVQbRt4MfLww/5LtowAk3VRC/FN7WbY2cLikybaPayL2Wj01WOAI4ELbn5E0CriTVJNq1q62Z9QutH2ZpK+XEB9gRG+tzbYfy63FzeoqxFyaTjWWQ9JY0vnSw4CtSZeybGZ7oxL3cUVf65u9ljXHGNtsjE4wVJPgiJrLAj5ceNz0JQ49l+DUyl/sO0nXwr3u8IXHfw98K+9ziaSyRqxZo491Zd1FtFTSJrYfLy7M5zrLOHe3k6SeVmYBq+f5Mmo7fyXVwL8C3GTbuZZcpj1Jl7BcCPye11q5SyVpLVLD1+ieZbZvqGJfg9VQTYLdktaz/RcA2/cCSNqQQitZ2WwvK6FGMkfSqcATpJP+1wJIKuv6NIC/StqttiVY0m7A03Wes6pOAK6T9A3SDwOkc4LHAcc2G9z28GZj9OF4YBrwP8CF+fxj2dYD3kk6zD4M+CWp1j+vrB1I+hfgaGAjYDawB3AL6cc19LA95CbgcNItWn8HjM3T2/Oyfyoh/tq9TFsAJwEXNBl7dVKiOB3YqbB8L+DDJb0/u5Gu3TsROChPJwGPAruX+DnsBJxHSoJ3Aj8qvqZ2n4DNgS+RLjFZTLrubusK9rMacCTpB+jTJcadS6oBzs7z2wKXtvp9bbdpSDaMwPJrBL9EujzDpMszvmn7qhJiP8qKF+p2A88As4Cv2276YuCqSVqXdG/yZF57f34NHGr7U60sW6vlBql1bd9cWLYD6Yfp7S6pFpov33oXqTY4CbiCdKnPEyXFv932myXNJv24vSppnu3t+31yBxmqh8PYvhq4una5pM/a/k6T4Q8F5tv+c455BPAB0q9uU++ppLn03flDKffcOl0mdIKkKaR/wiNINedLyog/ECf+K/QdahqgbM/NvdN8o4wdSDqP9AM0EzjJ+ZRNyRbk0yiXA7+S9ByrNkxtRxiyNcF6yrgBXtJdwL62n5X0d6SLXj8D7Ay8yfYHm4g9k/SPtoBekqFLuEha0ta81nnCQuAnpGv3+rxAexX38TR9nPh3ncaldtBTg6qzbvmlOU3uoxtYlGeLn3Mll7HkTkXGkS6NWlJm7MFuyNYE+1BGK9xw28/mx4cCZ9m+BLgkH3o04xpSi/D6wMWkk+V3Nxmz1gOki4zfbfthAElN9azTi8pP/Feor0aoZu9LBsB2ZX151rngfgrpgvJIgDU6sVPVMqq+w/MFugDvAH5TWNfUD4vt023vSWrIeQY4W9IDkk7INbgyvB/4M6lnmh9IegclX6Jhe5ntq20fQWqVfBiYJenTZe6nIndI+ljtwtzaemcv27ebq0nnGHvOb95CauT5lKRvtrBcbWlIHg4XehhZaRWwuu1mz9t9mXQr20JSxwZTbDt/4c61/ZZm4veyv12As4Edyzopn+OuCRxMqq39Pakl9zKXdJtV1Sf+q5IbjS4j9RZUvLxnFPA+50uv2lXN3TRfA9a2/ameC+7LOJwfSoZkEhwIkvYgHbJea3tRXrY1MMb2XSXEH0G6J3YaqbY5i3Q42dQ9t33sby3S7WGH2n5HCfGKJ/4vqujEf6Uk7UN6DQDzbP+mr+3bhaQ5PQ1okm4GvmX78jx/j5u8ZXGoiSTYZiT1nEc7kHTXwkWke20X9fnENjPQJ/7DaySdD/yFdMH9caRb/l7OLcW/jSS4okiCbUapL7sfA5fYfq7V5QmDj6TVSXeKrE86/XBPXr4XqTPXH7WyfO0mkmAIoaN14iUyIQxpA3XB/VARSTCEoWc+fVxwH1YUSTCEoWcgLrgfMuKcYAhDVO67cVqeVifdwnih7T+0tGBtJpJgCB2gqgvuh4JOvG0uhI4gaYSkgyRdQBpj50HSLZOhIGqCIQwxQ+WC+4ESSTCEISYuuF81kQRDCB0tzgmGEDpaJMEQQkeLJBhC6GiRBEMIHS2SYAiho/1/OYos6W6IBKcAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 321
        },
        "id": "prRk5HlHpTny",
        "outputId": "76692515-a271-447a-931a-e4643a23cdf1"
      },
      "source": [
        "#@title Difference matrix Correct - Incorrect \n",
        "difference_coherence = np.subtract(theta_array,theta_array_incorr)\n",
        "plt.imshow(difference_coherence)\n",
        "positions = (0,1, 2, 3,4,5,6,7,8)\n",
        "plt.xticks(positions, roi_labels,rotation = 'vertical')\n",
        "plt.yticks(positions, roi_labels)\n",
        "plt.colorbar()\n",
        "plt.clim(0,.25) \n",
        "plt.title('Theta Coherence Corr - Incorr: -500 - 0 ms')"
      ],
      "execution_count": 49,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "Text(0.5, 1.0, 'Theta Coherence Corr - Incorr: -500 - 0 ms')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 49
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUcAAAEeCAYAAAAZyNdUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deZxcVZ3//9ebLGRDCKDs+yqbLAEEF2QA2RcdlMTBQcVhnJGvOuooLj8X9Osyol/1O26Mky8gTALKYtTIohAUUSRBIIRFQliSKGsAgYSk03n//jink9uVqu7q1L3dVZ3P8/G4j9TdzjlVXfnUuffcc45sE0IIobcNhroAIYTQjiI4hhBCHREcQwihjgiOIYRQRwTHEEKoI4JjCCHU0VbBUdLnJF061OVYF5JmSXrvUJcjhFCOQQ2Okl4sLKskLSus/0OLabccWCVtJem/Jf1V0guS7pf0eUnjW0m3U0l6haRvSnos/40eyuubD2GZOvZHSJIlvVT4zv+wsE+Svirpmbx8VZIK+/eXNEfS0vzv/i2W5d8kPS7pb5KmStqwlfSGo0ENjrYn9CzAY8DJhW2XDWZZaknaFPg9MBY4zPZGwDHAJsAug1gOSRryGr2k0cCvgb2B44BXAIcBzwCHrEN6I2rWR5ZQzMrVK2eLZX9N4TtfDPLnAKcBrwH2A04G/jnnNxr4KXApMBG4GPhp3j5gko4FzgOOAnYAdgY+v25vZxizPSQL8AhwdM22zwFXAJcALwDzgEmF/VsDVwJPAQ8DH8jbjwNWAF3Ai8Bdefu7gftyWguAf+6jPF8E5gIb9HHM4cDtwPP538ML+2YBXwB+l/O7Hti8sP+1wK3Ac8BdwJtqzv3f+dxlwK7AnsANwBLgAeDtheMvAr4D/CLndRuwS2H/3oVznwA+mbdvQPpP8RApyF0BbNrgvb43nzuhj8/j1bnsz+W/1Sk1ZfweMBN4CTg6/80/DtwNLAdGrsP3Zhbw3vz6TcAi4CPAk8BfgXcXjh0LfB14NP/NbgHG5n2n5DI/l9N8dc13s1jOXQEDZ5N+1H+zjt95A7s22HcrcE5h/WzgD/n1m4HFgAr7HwOOW8dy/A/wpcL6UcDjDY7dMZf73cBC4FngfcDB+fN5DvjPwvG7Ajfnz/tp4PKqYkjVy9Bl3Dg4vgycAIwAvlz4gmwAzAE+A4wm/dotAI4tnHtpTXonkmp9Ao4AlgIHNijPH4DP91HeTfMX453ASGBKXt8s759FCjq75/+Us4Cv5H3bkILRCfl9HJPXX1k49zFSUBsJbJy/iO/O6wfkL9pe+fiLWFODGwlcBkzP+zYiBYmPAGPy+qF53wfz+9wW2BD4ATCtwfudDlzcx+cxCpgPfDL/Pf6OFKj3KJTxeeB1+T2PyX/zO4HtyEFqHb43s+gdHFcC5+fynJD/xhPz/u/k47fJ36fD8/venRSwj8nnfSy/l9GF7+bqcrImQFwCjM/bticFhu0HUHYDfwEeB64Cdizse77n75TXJwEv5Nf/BvyyJq2fAx9Zx8/wLuCMwvrmuWyb1Tm2571/P/8N30z6P3oN8Kr82T4JHJGPnwZ8qvA3f/1QxZhWlyG/fKvjFtszbXcDPyJdZkD6pXql7fNtr7C9APgvYHKjhGz/wvZDTm4m1ebe0ODwzUhBpZETgQdt/8j2StvTgPtJlz89/p/tP9teRqqV9dwXOhOYmd/XKts3ALNJ/5l7XGR7nu2VpJrwI7b/X87rT6Qa89sKx19t+4/5+MsKeZ1EqgV83fbLtl+wfVve9z7gU7YX2V5O+kE5vcFlYn+fx2uBCaQfgBW2byT9h51SOOantn+X3/PLedu3bS/Mn1EZuoDzbXfZnkm6ctgj35p4D/BB24ttd9u+Nb/vM4Bf2L7BdhdwASngHV5It145P2f7JdvLbD9mexPbjw2grEeQgs2epCD588JnP4EUIHs8D0zI9x1r9/Xs32gAeRfVy4t+0vtC/j5dT/phmWb7SduLgd+SfsAh/T12ALbOx9+yjmUccu0YHB8vvF4KjMlfoB2ArSU917OQai1bNEpI0vGS/iBpST7+BNKvZD3PAFv1Ua6tSZdnRY+SfjkblX1Cfr0D8Laasr++Jr+Fhdc7AIfWHP8PwJZN5LUdqQZbzw7A1YU07wO6qf8ZNvN5LLS9qrCt9vNYyNrqbQNA0icLjRXf7yPvXuXMPxA9ej6LzUk1l3qfRa+/ZX4PC1spe62axsftcz6/yT8kz5Fq8TuRbk1ACuqvKCTxCuBFp+pY7b6e/S/UyfcNhXznNShevbyol17BE4XXy+qs93z/Pka6UvujpHmS3tNHmm2tHYNjIwuBh/Ovdc+yke2e2lev4YVy69uVpFrBFrY3Id3/EvX9CnhLH40hfyEFl6LtSfeCmin7j2rKPt72VwrHuOb4m2uOn2D7X5rMa+c+9h1fk+6Y/Otf61fAsX201P8F2K7m86r9POoN+dRwGCjbX/Kaxor3NTquSU+TLv/qNab1+lvm2tl2tFD2tQ4sND72Ubs0a76P81hzlUR+Pa+wb79i6zWp0Wat4Gf7t4V8926Qb728nrD9TN/vqn+2H7f9T7a3JjUofVfSrq2mOxQ6KTj+EXhB0scljZU0QtI+kg7O+58Adiz8Zx1Nur/0FLBS0vGk+yWNfIP0C3qxpB0AJG0j6RuS9iMF1t0lvUPSSElnAHuRLiX7cylwsqRjc7nHSHqTpG0bHP/znNc7JY3Ky8GSXt3g+Npzt5L0IUkbStpI0qF53/eB/114f6+UdGqDdH5ECqZXStpT0gaSNsu1uxNIjUBLgY/l8r2JdIthehNlrFyuDU4FviFp6/y5H5Z/NK8ATpR0lKRRpPuzy0mNIpWQtHd+HGeEpAmkhqLFpNo7pPuZH87fua1zmS7K+2aRavgfyH/Tc/P2G9exOJcAZ0vaS9ImwKcLebVE0tsK3+tnST8Aq/o4pW11THDM9yBPIt1be5hUM/ghqfEC4Mf532ck3WH7BeADpP8IzwLvAGb0kf4S0j2nLuA2SS+QHmV5Hpiff1VPIn1pnyFdPpxk++kmyr4QOJV0G+ApUtD5dxp8/rnsbybdT+25gf9VUrDvL68XSA0NJ+fzHgSOzLu/RfoMrs/v7w/AoQ3SWU5qYb6f1PL9N9IP1ObAbbZX5DyOJ/0tvgv8o+37+yvjIPoo6QmE20kt918lPY3wAOk+8P8llf1k0mNlK5pNWNL2xUvmJmwBXE76HBeQ7j2elO95Qmoc+1ku7z2kJxF+AJDLdRrwj6RGoPcApw2kvEW2rwX+A7iJ1BD4KPDZdUmrjoNJ/39eJH3XPpjbBzqO0i2NEEIIRR1TcwwhhMEUwTGE0DEkHSfpAUnzJZ1XZ/+HJd0r6W5Jv+65v573dUu6My8Nb7GtPj4uq0MInSB3Qf0z6Z76ItK95Cm27y0ccyTpnvhSSf9C6ol2Rt73olPX5aZEzTGE0CkOITWOLsiNUdNJDZ2r2b7J9tK82tMbbJ1EcAwhdIpt6P0g/iJ6P7hf62zgl4X1MZJm544hp/WXWUeMjFLPaG3oMVQ3kljXltWOUjbqiaX9H9Sqim+ZaPSoStN318r+D2pB16vGVZo+wMhlFf8N/lbd9+hlXmKFlzfqNNGUY48c72eWdDd17Jy7l88jPbjf40LbF65LvpLOJPVPP6KweQfbiyXtDNwoaa7tRr3JOjc4jmE8h+qoytJf/O7D+z+oBdt+fXal6QO4a50eg2vayK23qzT97sefrDT9v7zzoErTB9h8brV/g1HXV/c9us2/bjmNp5d0c9t1zV3ZjtrqoZdtT+rjkMWknkw9tqVODzVJR5MGvzgiP68LQE9PMNsLJM0i9QdvGBzjsjqEUCHT7VVNLU24HdhN0k5KY1lOpqZjh6QDSA/Pn2L7ycL2ibl3FEqDNb8OuJc+dGzNMYTQ/lLfwXJuLdhembtOXkcagm6q7XmSzgdm254BfI00CMaPc1f0x2yfQhrg4weSVpEqhV8ptnLXE8ExhFAZY7rc3D3HptJLw9LNrNn2mcLroxucdyuw70DyiuAYQqhUWTXHwRbBMYRQGQPdERxDCGFtnVpzLL21Og9VVG/7mbm/4zxJd0n6YR5Lrva4N0lqZozEEEKbM9BtN7W0m0GpOUo6jjRJ0PH5IcwRwFmkMe6eG4wyhBAGnzFdHVpzHKzL6k8BHy08hNlNGqUZWB08v0kaWbpjJ+QJIdQwdHdmbBy0h8D3Bu6ot0PSGNIsgicDB9F7EqkQQgfrmSOhmaXdDHoPGUn75vHUHsrzsOxJmjjrwTzT2qV9nHtO7jg+u4vljQ4LIbQN0d3k0m4GKzjOAw4EsD3X9v6k0TLGDiQR2xfanmR70qj+p1MJIQwxA6vc3NJuBis4fhm4oGa2vZ7AeD9p1sCeKTSnEEIYNjq15lhFg8w4SYsK69+w/Q1JrwR+mVuqnyPNsHad7ZclnQP8QtJS4LfARhWUK4QwyAx0uTPHtyk9ONr1PwnbFwMXN9h3LeneYwhhGEk9ZNqvVtiM6CETQqiMEd0dOjJiBMcQQqVWOWqOIYTQS1xWhxBCHUZ0uTPDTGeWOoTQMaLmGEIINWzRHY/yDK6uLcdXOkPgNl+5tbK0gUEZp2TJew6rNP0Nn6+2R+wz+1Y7u+H2n6v2bwyw7LRDKk2/2slxy7Eqao4hhNBbapCJmmMIIdSIy+oQQlhL6j44YqiLsU4iOIYQKhM9ZEIIoYFVcVkdQgi9RYNMCCHUYUR39K0OIYS1rerQmmMV81Z/Ks9NfXeeK+ZQSY9I2rxwzOq5qSW9S9JT+dh5kn4iaVzZ5QohDD5bdHlEU0u7KTU4SjoMOAk40PZ+wNHAwiZOvdz2/rb3BlYAZ5RZrhDC0DDQ7Q2aWtpN2ZfVWwFP214OYPtpAKm5ew6SRgLjgWdLLlcIYYh0aoNM2aW+HthO0p8lfVfSEU2ed4akO4HFwKbAz0ouVwhhCBixys0t7abU4Gj7ReAg4BzgKeBySe+i/jgLxW2X5+latwTmAv9eL/3ivNUrl75UZtFDCBXpZoOmlnZTeolsd9ueZfuzwLnA3wPPABMLh20KPF3nXJNqjW9skPbqeatHjhtfdtFDCCVLg91GgwyS9pC0W2HT/sCjwCzgnfmYEcCZwE0Nknk98FCZ5QohDA2Tesg0s7SbshtkJgD/V9ImwEpgPukSuwv4nqS7AAHXApcWzjtD0utJwXoR8K6SyxVCGCIxEjhgew7QaATadzQ45yLgojLLEUJoD7ZKrRVKOg74FjAC+KHtr9Ts/zDwXlLl7CngPbYfzfvOAj6dD/2i7Yv7yit6yIQQKlXWM4z5ltx3gGNIV5i3S5ph+97CYX8CJtleKulfgP8gXZluCnwWmES62p+Tz2342GD7XeiHEIYNk6ZJaGZpwiHAfNsLbK8ApgOn9srPvsn20rz6B2Db/PpY4AbbS3JAvAE4rq/MouYYQqiMEV2rmm6J3lzS7ML6hbYvLKxvQ+8ed4uAQ/tI72zgl32cu01fhYngGEKo1ACeYXza9qQy8pR0JukSutmOKGuJy+oQQmVK7iGzGChOSblt3taLpKOBTwGn9HRlbvbcogiOIYRKrWKDppYm3A7sJmknSaOBycCM4gGSDgB+QAqMTxZ2XQe8WdJESROBN+dtDXXsZfWoJ5ay7ddn93/gOhqMeaWrtunU3w91EVoy/sqhLkHrxv/6vkrTr3bm8NbZlDbYre2Vks4lBbURwFTb8ySdD8y2PQP4Gul56x/nAW8es32K7SWSvkAKsADn217SV34dGxxDCJ2hzEElbM8EZtZs+0zh9dF9nDsVmNpsXhEcQwiV6elb3YkiOIYQKpP6Vkf3wRBCqFFu98HBFMExhFCpJnu/tJ0IjiGEypTZWj3YIjiGECpjxMrmuw+2lQiOIYRKxWV1CCHU6OTW6lKakSTdJOnYmm0fkvRLSffk9XGSLpM0V9I9km6RNKGM/EMI7atTp0koq0TTSP0ciyYDXy6sfxB4wva+tvchDSfUVVL+IYR21OSgE+1YuywrOP4EODF3BkfSjsDW9B4/bSsKo2DYfsD2ckk7Sro/1yrvk/QTSeNKKlcIYQiVPNjtoColOOYO3H8Ejs+bJgNX0Hv8hqnAxyX9XtIXa2Yp3AP4ru1XA38D/rVePsV5q7v8chlFDyFUyMDKVRs0tbSbMktUvLSenNdXs30nsDNp1IxNSfM/vDrvXmj7d/n1paTpWddSnLd6lMaUWPQQQlXW98tqgJ8CR0k6EBiXZyLsxfaLtq+y/a+kIHhCz67aQ0ssVwhhiJQ82O2gKi042n4RuIl0+Tytdr+k1+VBJsn3JvcCHs27t5d0WH79DuCWssoVQhha6/U9x4JpwGuoExyBXYCbJc0lTZ84G+gZzvQB4P2S7gMmAt8ruVwhhKHgzr2sLvUhcNvXwJqfANuPAPvk15cAl9Sek0frXWn7zDLLEkIYep38EHj0kAkhVCb1rW6/luhmDHlwLNYuQwjDj6PmGEIIa2vHxpZmRHAMIVTGjnuOIYRQV1xWhxDCWkR3NMgMMht3rags+SXvOaz/g1qw6dTfV5r+YNCGG1ab/shqv57L3rhXpekDdI+tNjCMu+q2StNvVTzKE0II9Tjdd+xEERxDCJWK1uoQQqhhokEmhBDqaM9+082I4BhCqNSqVREcQwihFzsuq0MIoa64rA4hhDo69VGeSp5QlXSaJEvas7DtEEm/kfSApD9J+mFxlkFJ10j6QxXlCSEMHVtNLe2mqsf3p5CmOpgCIGkL4MfAx23vYfsA4Fpgo7x/E+AgYGNJO1dUphDCIIs5ZAokTSDNHng2a2YjfD9wse3VfeZs/8T2E3n1rcDPgOmFc0IInc7l1hwlHZevPudLOq/O/jdKukPSSkmn1+zrlnRnXmb0l1cVNcdTgWtt/xl4RtJBpMFs15qNsGAKad6Zafl1Xb3mrWZ5mWUOIVTFTS79kDQC+A5wPGmCvimSajvIPwa8C/ifOkkss71/Xk7pL78qguMUUg2Q/G/DYAerL7l3A27JAbVLUt2RwXvNW021gx6EEMpRYs3xEGC+7QW2V5Diy6m98/Ijtu8GVrVa7lJbqyVtCvwdsK8kAyNIvwkXk+4p/rTOaW8nzTj4cJ5s6xWkgPqpMssWQhgaJbZWbwMsLKwvAg4dwPljJM0GVgJfyRMCNlR2zfF04Ee2d7C9o+3tgIeBXwFnSVr9RiS9NdcapwDH5eN3JAXRuO8YwjDQ07e6yZrj5j23zfJyTsnF2cH2JOAdwDcl7dLXwWU/5zgF+GrNtitJwW4ycIGkV5GqvL8B7gd2AFY/wmP7YUnPSzrUdnsPVhdC6JvBzXcffDoHr0YWA9sV1rfN25orir04/7tA0izgAOChRseXPW/1kXW2fbuw+oY6p21T55wDyyxXCGEIlXdZfTuwm6SdSEFxMqkW2C9JE4GltpdL2hx4HfAffZ3TmeOXhxA6RHOX1M00yNheCZwLXAfcB1xhe56k8yWdAiDpYEmLgLcBP5A0L5/+amC2pLuAm0j3HO/tK7/oPhhCqFaJ3QdtzwRm1mz7TOH17aTL7drzbgX2HUheERxDCNWJUXlCCKGBDh14IoJjCKFaUXMMIYQ6ouY4uDR6FCO33q7/A9fRhs+33PtoyFU9r7SXV9u/XSNGVJr+io2rTR+oPDCM2L3P55hbokduaT0REzXHEEKop1MHu43gGEKoVgTHEEJYm2L2wRBCqNHkWI3tKIJjCKFCigaZEEKoK2qOIYRQRwTHEEKoo0ODY2lDluV5qi8trI+U9JSknxe2nSbpbkn3SZor6bSy8g8htCGn1upmlnZTZs3xJWAfSWNtLwOOoTBKr6TXABcAx+TRvncCbpC0IE+IE0IYjtb3mmM2Ezgxv+6ZbrXHR4Ev2X4Y0nQIwJeBfweQ9AFJ9+aa5XRCCGEIlR0cpwOTJY0B9gOKc8DszdpzV8/O2wHOAw6wvR/wvnqJF+etXtG9rNyShxAqITe3tJtSg2O+PN6RVGuc2ffRa7kbuEzSmaSpE+ulv3re6tEjxrZU1hDCILGaW9pMFXPIzCDdW5xWs/1e0rSrRQcBPXM8nAh8BzgQuF1StKSH0Ok8gKXNVBGApgLP2Z4r6U2F7RcAP5Z0o+1HJO0IfBI4XdIGwHa2b5J0C2lWsQnAcxWUL4QwiNSho/+VHhxtLwK+XWf7nZI+DvxM0iigC/hY3j4KuFTSxoCAb9uOwBjCcNCGtcJmlBYcbU+os20WMKuwfhVwVZ3juoDXl1WWEEIbWd+DYwgh1GrXluhmRHAMIVSrDVuimxHBMYRQqWiQCSGEeuKyOoQQasQ9xxBCaCCC4+By10q6H3+ysvSf2be6ObEBxl9ZafIAaGS1f96q55VetXRppem/sG0VHcR6m7C42htu3X9+qLK07ZLmJY/gGEIIa4vL6hBCqCeCYwgh1IgGmRBCaCCCYwgh1BHBMYQQehOde1ld/bMMIYT1W4mD3Uo6TtIDkuZLOq/O/jdKukPSSkmn1+w7S9KDeTmrv7yi5hhCqI7L61staQRptoBjgEWkGQNm2L63cNhjwLtIE/oVz90U+CwwKZWKOfncZxvlN+g1R0ndku6UdI+kH0sal7e/ONhlCSEMgvJqjocA820vsL2CNKHfqb2ysh/Jc1nVhuRjgRtsL8kB8QbguL4yG4rL6mW297e9D7CCBjMNhhCGhwHMPrh5z+yieTmnJqltgIWF9UV5WzMGfO5QX1b/ljSFawhhuGq+QeZp25MqLMmADFmDTJ5d8Hhg7gDOWT1vdZdfrq5wIYRylDv74GKgOOjBtnlbJecORXAcK+lOYDbp5ul/N3ticd7qURpTWQFDCOXRquaWJtwO7CZpJ0mjSbOUzmiyGNcBb5Y0UdJE4M15W0NDcVm9zPb+Q5BvCGEIlPWco+2Vks4lBbURwFTb8ySdD8y2PUPSwcDVwETgZEmft7237SWSvkAKsADn217SV35Dfc8xhDDclfgQuO2ZwMyabZ8pvL6ddMlc79ypwNRm82qnh8DHSVpUWD481AUKIbSo3HuOg2rQa4715rfO29spUIcQSqC8dKK4rA4hVKsNa4XNiOAYQqhUTM0aQgj1RM0xhBBqxEjgIYTQQATHEEJYW9QcB1nXq8bxl3ceVFn623/u1srSHizL3rhXpemv2Ljaeaurnld66wuq/xv/9SOHV5r+xC1eVVnaerqk8BDBMYQQapQ42O1gi+AYQqhW1BxDCKG3Tp5gK4JjCKFaERxDCGFtcmdGxwiOIYTqRINMCCE00JkVxwiOIYRqdWqDTOlP2UraUtJ0SQ9JmiNppqTd874PSXpZ0saF4zeTdJOkFyX9Z9nlCSEMsQ4d7LbU4ChJpPkbZtnexfZBwCeALfIhU0hzOLy1cNrLwP8HfLTMsoQQ2kCTc1a3Y+2y7JrjkUCX7e/3bLB9l+3fStoFmAB8mhQke/a/ZPsWUpAMIQw3HVpzLPue4z7AnAb7JgPTgd8Ce0jawvYTJecfQmgjArSqDSNfEwZz3pYpwHTbq4ArgbcNNAFJ50iaLWl299KXSi9gCKF8nXpZXXbNcR5weu1GSfsCuwE3pNuSjAYeBgbUAGP7QuBCgLFbbteGH2cIoZc2vWRuRtk1xxuBDSWd07NB0n7At4HP2d4xL1sDW0vaoeT8QwhtRquaW9pNqcHRtoG3AEfnR3nmAV8G3kRqxS66mnQfEkmPAN8A3pXnrK52IMIQwuCJBpnE9l+Atzdx3IcLr3csuxwhhDbgzm2QiR4yIYRKtWNjSzMiOIYQqhXBMYQQeovBbkMIoR47LR0ogmMIoVJRcwwhhDra8RnGZkRwDCFUx0A8yjO4Ri4zm89dUVn6y047pLK0Acb/+r5K0wfoHltx1/mKv/MTFldb5fjrRw6vNH2Arb5+a6Xpd1eYtr2ypITKSWawDebAEyGE9VCZA09IOk7SA5LmSzqvzv4NJV2e998mace8fUdJyyTdmZfv155bq2NrjiGEDlFSa7WkEcB3gGOARcDtkmbYvrdw2NnAs7Z3lTQZ+CpwRt73kO39m80vao4hhEqVWHM8BJhve4HtFaTxYU+tOeZU4OL8+ifAUXmGggGL4BhCqIxy3+pmliZsAywsrC/K2+oe43TT9Hlgs7xvJ0l/knSzpDf0l1lcVocQqtV8u9rmkmYX1i/MY7iW4a/A9rafkXQQcI2kvW3/rdEJERxDCJVS8/ccn7Y9qY/9i4HtCuvb5m31jlkkaSSwMfBMHk5xOYDtOZIeAnYHZtNAXFaHEKrT7FiOzcXP24HdJO0kaTRpPNgZNcfMAM7Kr08HbrRtSa/MDTpI2pk0M8GCvjKrPDhKsqRLC+sjJT0l6ed5/V15/U+SHpR0naTqH0ALIQwCr+lf3d/SX0rpHuK5wHXAfcAVtudJOl/SKfmw/wY2kzQf+DDQ87jPG4G7Jd1Jaqh5n+0lfeU3GJfVLwH7SBprexmpGb62Kny57XMBJB0JXCXpSNvVPykdQqhUmYPd2p4JzKzZ9pnC65epM3mf7StJE/s1bbAuq2cCJ+bXU4BpjQ60fRNpEq1zGh0TQugQjjlk+jMdmCxpDLAfcFs/x98B7Fl5qUII1SvpsnqwDUprte27czeeKdRUiRuo+9BmntXwHIANx2xSVvFCCFVqv7jXlMF8lGcGcAFpJsLN+j6UA0g3XHspzlu90cbbduhHHsL6ZQCP8rSVwQyOU4HnbM+V9KZGB0k6glQ7PHKwChZCqFAEx77ZXgR8u8HuMyS9HhgHPAz8fbRUh9D5ZKPuCI512Z5QZ9ssYFZ+fRFwUdXlCCEMkag5hhBCHREcQwihhhnIwBNtJYJjCKFS0VodQgj1RHAMIYQaNqzqzOvqCI4hhGp1ZmyM4BhCqFbccxxk+ttSRl3fcBDflo2qLOVkMH5Mx13V3/gerRmx+y6Vpt/954cqTX/iFq+qNH2odl7pjhHBMYQQahgocTzHwRTBMYRQoWiQCSGE+uKyOoQQasRldQgh1GNwXFaHEMLa4rI6hBBqxGV1CCE00KGt1X3OPijpJknH1mz7kKRfSronr4+TdJmkuZLukXSLpLUGuA0hrI+anHmwDS+9+5uadS2uKB0AAA0nSURBVBowuWbbZODLhfUPAk/Y3tf2PsDZQFd5RQwhdCyTao7NLG2mv+D4E+BESaMB8vSqWwMLC8dsBSzuWbH9gO3l+fhrJM2RNC9Pq0re/qKkr+Xtv5J0iKRZkhZIOqWctxZCaAvDseZoewnwR+D4vGkycAW9Z6KdCnxc0u8lfVHSboV977F9EDAJ+ICknilZxwM32t4beAH4InAM8Bbg/EblkXSOpNmSZnexvOk3GUIYQsMxOGbFS+vJeX0123cCOwNfAzYFbpf06rz7A5LuAv4AbAf0BM4VwLX59VzgZttd+fWOjQpi+0Lbk2xPGsWGTRQ9hDC0nFqrm1naTDOt1T8F/o+kA4Fxtufky+vVbL8IXAVcJWkVcIKkLYCjgcNsL5U0CxiTT+myV/9UrIJUDbS9SlK0oIcwXBjc3ZljE/Vbc8yB7ybS5fO02v2SXidpYn49GtgLeBTYGHg2B8Y9gdeWWfAQQofo0MvqZmtp04CrWbvlGmAX4HuSRAq2vwCuBEYD75N0H/AA6dI6hLA+Ge7TJNi+BlBh/RFgn/z6EuCSOqctZ01DTm16EwqvP9doXwhhGGjDWmEz4v5eCKFSHs41xxBCWCc2dEdwDCGEtXXokGXNPOcYQgjrxIBXuamlGZKOk/SApPmSzquzf0NJl+f9txUfO5T0ibz9gdoxI+qJ4BhCqI7zYLfNLP2QNAL4Dqmhdy9giqS9ag47m/QI4a7A/wG+ms/di/S0zd7AccB3c3oNRXAMIVSqxJrjIcB82wtsrwCmA6fWHHMqcHF+/RPgqPyY4anAdNvLbT8MzM/pNdSx9xxf4Nmnf+WfPDqAUzYHnq6qPMMg/YHn8UDF6Q/cwNJ/vOL0B67d0t+h1Qxf4NnrfrXqis2bPHyMpOJk9BfavrCwvg29B71ZBBxak8bqY2yvlPQ8sFne/oeac7fpqzAdGxxtv3Igx0uabXtSVeXp9PQHI49If3inX4/t4wYzvzLFZXUIoVMsJg1g02NbCsMl1h6Tx2nYGHimyXN7ieAYQugUtwO7Sdopj+MwGZhRc8wM4Kz8+nTS0IjO2yfn1uydSCOE/bGvzDr2snodXNj/Iet1+oORR6Q/vNOvVL6HeC5wHTACmGp7nqTzgdm2ZwD/DfxI0nxgCXk8iHzcFcC9wErg/bb7HC5I7tB+jyGEUKW4rA4hhDoiOIYQQh0RHEMIoY4IjqEXSQdXmPYGkl5RVfqFfD5UdR5h+BuWDTKSfkbvGRJ7sd3S9K+SxgBnAM8CPwM+BrwBeAj4gu2WeznkPP4VeD3pvdwCfM/2y62mXSevvYApeXmuzAeFJf0P8D6gm/QoxiuAb9n+Wll51MnzMdvbV5T2nrbvLymtPYBzgD3zpvuA/7I98L5H9dM/CfgCqafLSNKA1bZd+Q/UcDBcg+MR+eVbgS2BS/P6FOAJ2//WYvpXAF2kKWYnAveQguTrgf1tn9RK+oU8XmBN2d8BbGL7ba2mndPfkTUBsYv0H2hSHuW9NJLutL2/pH8ADgTOA+bY3q/MfGryXGh7u/6PXKe0Swm8kg4jTUr3A+BPpMB1APBPwFtttzytSH6c5a3AXA/H/+gVG5bPOdq+GUDS12tqQT+r6bu5rvayvU9+An+R7Z5gfG2eirYM+9gujjhyk6R7y0hY0u9JNbjpwN/bflDSw2UHxmyUpFHAacB/2u6SVPV/1JbSl/TtRruATVpJu+AzwBTbswrbrpF0I/BZGkwxMkALgXsiMK6bYRkcC8ZL2tn2AoD8ZPz4EtJdAasfSv1Lzb6y5qG8Q9Jre2oQkg4FygjsAE+QOt1vAbwSeJAWA0offgA8AtwF/EbSDsDfWk1U0gukMvfMbdRTfgFjW0z+3cBHyFMG15jSYto9dqkJjED6YZdU1sPaHwNmSrqZwnux/Y2S0h/WhuVldQ9Jx5F6BSwg/afZATjH9vUtpvskqdYl0r3H6T27gLfb3qKV9HMe9wF7AI/lTduTxsFZSbpv1NJlqaSNSZdcU0hdqTYBjrXdZ5eqVuXho0bYXlllPq3ItbdP2761zr6Hbe9UQh5zbB/UYN8dtg8sIY/rgReBuaT54QGw/flW014fDOvgCGlkYNbc8L7fdr3awEDTPKuv/bYv7mt/k3n0OVyU7YEM19ZfXq8C3k4KlNuXdb9O0j7Av5MGGAWYB1xge24JaY8hNfTsCtxN6kpWSsCVtCnwsu2lZaTXII+eH9i1dlHeD+w9tvdpNZ311bAMjvlxkS1sP5jX38aaS63rbD8xZIXrh6RxQJftrry+B3AC8Kjtq0rKYwywke2nara/CtjM9n0l5HEqcAHwZdbcDpgEfAL4qO2ftpj+5aSGpN+S7s89avuDraTZRJ7bAZPLaGkfpB/Y/wB+1eqV0vpquAbHC4FbbV+U1x8EfgmMA1bafl+L6dc+KmTSIKI32b60/llNp/0b4OzcSLIraeSQy0jDwv/R9idaST/ncSFwbW2wlfQW4M22/6WEPO4CTq1t5Mmt5D+1/ZoW059re9/8eiTps2n5UrROPq8E3kaqVW8NXG37o2XnU4V8X3Y86X5jF/Eoz4AM1waZg4F/Lqy/aPsDAJJuKSH9C+ps2xQ4U9I+ttea+GcAJvbUeElDL02z/b/yEE1zSDWvVh1k+5zajbavlvTFEtIHGFmv9dv2I7n1ulVdhTRXpluZ5ZC0Eel+7DuA3UmP3Oxke9sS86gdaquXVp/FzWls1Goa67PhGhxH1jy+8M7C65Yfxeh5VKhW/sLPIT3Lt87JF17/HfC1nOcKSWXNcTmuj31l9ZpaKWl7248VN+Z7qWXcG3yNpJ5WbwFj83oZtaMnSTX2TwO32HauVZfpMNKjNtOA21jT6l4qSRNJDW5jerbZ/k0VeQ03wzU4rpK0pe3HAWzfAyBpGwqtdmWz3V1CDeZuSReQRineFbgeQFJZz9cBPCnpkNqWaUmHAE81OGegPgv8StKXSD8YkO45ngd8vNXEbfc5c1yLPkEaB/C7wLR8f7NsWwLHkC7X3wH8gnSVMK+sDCS9F/ggadTrO4HXAr8n/eiG/tgedgtwJqmr2huBjfJyRN72jyWkv2mdZRfg88BlLaY9lhRAvgW8prD9cOCdJX0+h5CePfwccHJePg88DBxa4t/hNcAlpOA4B/hR8T21+wLsDHyS9CjMy6TnBnevIJ8NgXeRfpjOLTHduaQa4515fU/gqqH+XDtlGZYNMrD6GcdPkh4jMekxkq/Y/mUJaT9M7weQV5HmqZgFfNF2yw85V03SFqS+2/uw5vP5NXCG7fcPZdmGWm4I28L27wrb9iX9YB3hkmqt+TGzE0m1xx1JQ/lPtd3n3CYDSP922wdLupP0o7dc0jzbe/d7chi2l9XYvha4tna7pA/Z/maLyZ8BLLT915zmWcDfk36lW/pMJc2l70EzSumT7PQ402clHUj6z3kWqaZ9ZRnpD0aDQ4W+SU3Dl+25ebSfL5WRgaRLSD9MM4HPO9/6KdmifDvmGuAGSc8CpT0fO9wN25pjI2UMHCDpDuBo20skvZH0MO//AvYHXm379BbSnkn6D7iIOkHSJTz8LWl31gw68TRwOenZw5bnKS7k8RR9NDi4QaNWO+ipcTXYt/oRohbzWAW8lFeLf+dKHrfJg7FsTHqEa0WZaQ9Xw7bm2IcyWgVH2F6SX59Bmnz8SuDKfAnTiutILdRbAVeQbtL/qcU0a91Penj6JNvzASS1NFJRHZU3OFSor8avVvttA2C7srFUG3QkOJD0oHwExiatj4PdllFVHpEfPAY4CrixsK+lHxzb37J9GKkB6RlgqqT7JX021/jK8Fbgr6SRfv5L0lGU/CiJ7W7b19o+i9RKOh+YpTR7XLubLemfajfm1t85dY5vN9eS7mH23D/9Palx6f2SvjKE5eoow/KyujBiy1q7gLG2W70v+ClSl76nSQNCHGjb+Yt4se3XtZJ+nfwOAKYC+5XVGJDTHQ+cSqrd/R2pZflql9TdrOoGh6rkxqqrSaMvFR9DGg28xfkRsXZV03voC8Cmtt/f05GgjNsC64NhGRwHg6TXki59r7f9Ut62OzDB9h0lpD+S1Gd4Mql2Oot0WdpSn+Q+8ptI6iZ3hu2jSkiv2OAwvaIGh0pJOpL0HgDm2b6xr+PbhaS7exruJP0O+Jrta/L6XW6x6+b6IoJjm5HUc5/uBFIvjemkvsgv9XlimxnsBoewhqRLgcdJHQnOI3V9XJpbrm+O4NicCI5tRmkswf8BrrT97FCXJ3QeSWNJPWO2It3GuCtvP5w0yO6PhrJ8nSKCYwgh1LE+PsoTwrA2WB0JhrsIjiEMPwvpoyNBaE4ExxCGn8HoSDDsxT3HEIapPHbm5LyMJXXlnGb7z0NasA4RwTGE9UBVHQmGs/Wx+2AI6wVJIyWdLOky0hxKD5C6joYmRM0xhGFmuHQkGGoRHEMYZqIjQTkiOIYQQh1xzzGEEOqI4BhCCHVEcAwhhDoiOIYQQh0RHEMIoY7/H8MPoRNzH0eaAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 432x288 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "oK-RHJ11qDSH"
      },
      "source": [
        "#Granger Causality Analysis"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "UYTYQ0UrqGNT"
      },
      "source": [
        "#@title Prepping Averaged Data \n",
        "session11_lfp_corr = session11_lfp[:,correct_right_idx==1,:] # pulling correct trials only \n",
        "session11_lfp_incorr = session11_lfp[:,error_idx==1,:] #pulling incorrect trials only \n",
        "session11_lfp_corr_avg = np.mean(session11_lfp_corr,axis=1) #averaging \n",
        "session11_lfp_incorr_avg = np.mean(session11_lfp_incorr,axis=1)\n"
      ],
      "execution_count": 56,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "5s2A5xJ-tX2o"
      },
      "source": [
        "#@title Running Grangers \n",
        "TR = .01 #time res 100 hz = 10 ms  \n",
        "fb = [5,8] # theta band \n",
        "\n",
        "G_corr, g2_corr = granger_causality_analysis(session11_lfp_corr_avg,TR,freq_band=fb)\n",
        "G_incorr, g2_incorr = granger_causality_analysis(session11_lfp_incorr_avg,TR,freq_band=fb)"
      ],
      "execution_count": 57,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 320
        },
        "id": "qLsH9vJUu3yd",
        "outputId": "e0e90417-915b-4875-f692-7c2df0a06023"
      },
      "source": [
        "#@title Plotting Correct Right Granger\n",
        "fig01 = drawmatrix_channels(np.transpose(g2_corr), labels, size=[5., 5.], color_anchor=0)"
      ],
      "execution_count": 67,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAR8AAAEvCAYAAAB8LcXgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAZ70lEQVR4nO3de7hdVX3u8e+7dhIgBCKIcjOyDSLKTUw4oJYewUsRsN7wSOJDC0/bx8cKj4qXQzm0gjxU2ppDT6lWa5UColzkZlVA2yOhRAW5GAiRRi6JhehBUUhNSEiy1+/8McfamXtl7+y1k7nWWJf382Q9e84xxhxrrLWyf3vMMeYcSxGBmVmn1XI3wMwGk4OPmWXh4GNmWTj4mFkWDj5mloWDj5ll4eBjZlk4+JhZFn0bfCSNSFoqabmkByR9TFKtlH+0pMWSHpF0v6RvSzo8Z5vNBsm03A1oo/URcSSApBcDXwN2B86XtDdwHfC+iPhBKnMscCCwLFN7zQaK+vX2CklrI2JWaX8ucA+wF3AhUI+I83O1z2zQ9e1pV7OIeBwYAl4MHArcn7dFZoNtYILPtki6W9LDkv4ud1vMBsXABJ902jUC/BJYDsxr5EXEMcBfALPztM5s8AxE8JH0IuALwGejGOT6HHCGpNeXis3M0jizAdXPA84jFDNX04HNwFeASyKinvJfC/w1sD9Fb+hp4MKIuDdPi80GS98GHzPrbgNx2mVm3cfBx8yycPAxsywcfMwsCwcfM8vCwcfMsnDwMbMsHHzMLAsHHzPLoi8WE/uvdesrvUx79113UZX1mdnW3PMxsywcfFokae0E6adJerC0VvSXJL1gnHLHSfpW+1tq1hv64rQrF0lvBc4GToyI1ZKGgNOBvYFnszbOrMs5+OyY84CPR8RqgIgYAS5rZKbg9H+A54AlWVpo1qV82rVjJlwLWtLOwD8Bvw/MB/bpYLvMup6DT0UkHZ6+J+wxSacCrwRWRsQjafXEqzI30ayrOPjsmNG1oCNiWfqesFuBXbK2yqwHVB58BmxW6GJgkaSXlNIagec/gGFJB6b9hR1tmVmX68iAc5/MCs2U9GRp/5KIuCQtTn9rek3PAg8B34mIDZLeD3xb0nPAncBunW+2WXeqfA3n5m8KTWl3Ap+MiNsnOKZ5VmhuRLyt1ef0Fc5mvadTYz6eFTKzMTo+4OxZITODzgUfzwqZ2RidCj6eFTKzMdox2+VZITObVF98Y6lnu8x6j69wNrMsHHzMLAsHHzPLoi/GfF7+pzdW+iI+/YFjqqyO9756f48hmTVxz8fMsnDwyUjSeeku/wfTVd/HSFolaa9SmdG7/CWdIelXqexySddLmpnvFZhtPwefTCS9DngbMC8ijgDeDDzRwqHXRsSREXEosBE4tY3NNGsbr+Gcz77A0xHxPEBEPA0gtTY8JGkasCvwTLsaaNZO7vnk811gjqSfSvoHSW9o8bhTJS0FVgN7At9sWwvN2sjBJ5OIWEuxhMj7gV8B10o6Axhv5q6cdm26MXcfYBnwiTY31awtHHwyioiRiFgcEecDZwGnAL8G9igV2xN4epxjg6LX89870VazqrVjDWfP4LRA0sGSDiolHQn8DFgM/EEqMwScBoy7AiRwLPBYG5tp1jaVDjg3zeA8nwLOjBYOvTYizkp1fI1iBuefq2xbF5oF/H1aRH8z8CjFKdgm4POSHgAE3MbYBdZOlXQsxR+OJ4EzOtlos6pUPdvlGZwWRcR9wOsnyH7fBMdcDlzepiaZdVTVp12ewTGzllQafDyDY2atqnzA2TM4ZtaKSoOPZ3DMrFVVDzh7BsfMWlJp8PEMjpm1ylc4m1kWDj5mloWDj5ll0Rfr+bzrpFdUWt/RL9m90vqef+apSteY3mmPvb0mtPU893zMLAsHnz4g6XZJJzSlfUTSrZIeSvszJX1V0jJJD0laImlWnhabOfj0i6uBBU1pC4CLS/sfBp6KiMMj4jDgjymuvzLLwsGnP1wPnCxpBoCkYWA/xi5Ivy/FjbsARMSKtOzJsKT/SL2ihwdlPSXLz8GnD0TEb4AfASempAXAdYy9efcy4BxJP5R0UdNtMAcD/xARrwL+C/hgB5ptA87Bp3+UT70WpP1REbEUmAt8huLG3nskvSplPxER30/bV1HcX2fWVg4+/eMbwJskzQNmpltdxoiItRFxY0R8kCLInNTIai7a3qaaVRR8PNuSX1pL6XaK06urm/Ml/Y6kPdL2DOAQihUHAF6alsCF4h68Je1vsQ26qno+nm3pDlcDr2ac4AMcCNwhaRnwY+Be4IaUtwI4U9LDFOsufb4DbbUBV9UVztcDF0maEREbtzHb0vhLS0SsgNGZmduA+4B5wHLgDyPiuYraNjAi4maKJUsa+6uAw9L2lcCVzcek9bU3R8RpnWmlWaGSno9nW8xsqqoccPZsSw+KiFXpNNiso6oMPp5tMbOWVRZ8PNtiZlNR9XU+nm0xs5ZUvYazZ1vMrCW+wtnMssi+kmG5d2Rmg8M9HzPLQsU3FPe2NWufq/RF7PrMyiqro776p5XWF4ceV2l9O+22h9eEto5zz8fMsnDwsQlJeqekkPTKUtrRkv5d0gpJP5b0pfLKh5JulnRXnhZbL3HwsW1ZSHHB50IASXsDXwfOiYiDI+I1FDcF75byXwDMB2ZLmpunydYrHHxsXGmtpWMplj5p3LN3JnBFRPywUS4iro+Ip9Luu4FvAtew9RIrZmM4+NhE3gHcFhE/BX4taT7FJRFb3bNXspDi6var07bZhBx8bCILKXowpJ/bDCbplOwgYEkKWJsk+fotm1D2iwyt+0jaE3gjcLikAIYoVhq4gmJM5xvjHPZeivvyVqZbZnanCFjndaLN1nva0vPxLEnPew/wlYg4ICKGI2IOsBL4N+B0Scc0Ckp6d+r1LATemsoPUwQpj/vYhNp12uVZkt62ELipKe0GimCyAFiU/og8DJwAvBA4ABj94xERK4E15UBlVlb5Fc5plmQFcDzwzYg4WNKFABHxyQmO+SPgKOApYFNEfHoqz+krnHeMr3C2HNrR8/EsiZlNqh3Bx7MkZjapSme7PEtiZq2quufjWRIza0nVwcezJGbWEq/nMw7Pdpm1n2+vMLMsHHzMLAsHHzPLoi9uLN31N49VWt+KPz+30voO/sTZldZ386r1ldb39v2fr3TMbMae+3gMySblno+ZZeHgYx2RVjm4qrQ/TdKvJH2rlPZOSQ9KeljSMknvzNNa64S+OO2ynrAOOEzSLhGxHngLsLqRKenVwCLgLRGxUtLLgH+V9HhEPJinydZO7vlYJ90CnJy2GzcTN3wc+HS6yLRxsenFwCcAJH1I0k9Sz+garOc5+FgnXQMskLQzcARwdynvULZe+eDelA7wZ8BrIuII4APtbqi1n4OPdUw6fRqm6PXcMsXDHwS+Kuk0YHPFTbMMHHys0/6FYmzn6qb0n1DcVFw2H1ietk8GPgfMA+6R5PHKHldZ8PFshrXoMuBTEbGsKX0RcK6kYYD0838B/1tSDZgTEbcD5wCzgVkdaq+1SZV/PTybYZOKiCeBS8dJXyrpHOCbkqYDm4D/mdKnA1dJmg0IuDQinu1ow61yVXddG7MZ17NlNuN3U95WsxmSGrMZfyDpQxQDiZuBn0SE1/TpIxGxVU8lIhYDi0v7NwI3jlNuE8W3p1ofqXrMx7MZZtaSSoOPZzPMrFXtmO3ybIaZTaodv+CXAc9GxDJJx5XSFwFfl/S9iFhVms14T3k2Q9ISiiVXZwEeVDTrU5UHH89mmFkrKgs+ns0ws6nwFc5mloWDj5ll4eBjZln0xfd2Pb/m15W+iGnP/GeV1VGfvku19c16UaX1bZo+s9L6prVhBeedd9nF60L3Gfd8zCwLBx/rSZJGJC2V9JCkr0uamdLX5m6btcbBx3rV+og4MiIOAzbi+wF7joOP9YM7gZfnboRNjYOP9bR0D+CJQPPiZNblfPOm9apdJC1N23cCX87ZGJs6Bx/rVesj4sjcjbDt59MuM8ui48HHU6TWZjMlPVl6fDR3g2x8OU67RrvLkr5KMUV6SYZ2WA8bbxWFlO7efI/I/UF5itRsQGULPp4iNRtsOU67PEVqZnnHfMxscOUe8zGzAdVNwcdTpGYDpOOnXZ4iNTPorp6PmQ0QBx8zy8LBx8yy6Iu72mP6zpXWt3n2fpXWV9u4rtL6HlhT7Vrqu51/aqX1Db/79yqtD+DWD3+p0i8JOPHxpV6QPjP3fMwsCwcfs0TSPpKukfSYpPsk3SLpFSnvI5I2SJpdKv9CSbdLWivps/la3pscfMwASQJuAhZHxIERMR84F9g7FVkI3AO8u3TYBuAvgI93sq39wsHHrHA8sCkivtBIiIgHIuJOSQcCs4A/pwhCjfx1EbGEIgjZFDn4mBUOA+6bIG8BcA3FjdAHS9p7gnI2BQ4+ZpNbCFwTEXXgBuB/ZG5PX+iLqXazCiwH3tOcKOlw4CDgX4thIWYAKwEPMO+gyns+njGwHvU9YCdJ728kSDoCuBS4ICKG02M/YD9JB+RqaL+oNPh4xsB6VUQE8C7gzekP53LgYuA4iv/TZTdRjAMhaRXFGuRnpNUYDulYo3tc1add484YAJRmDD4InAf8c8pfByyR5LWcLauI+Dnw3hbKfbS0PdzONvWzqk+7PGNgZi3p5GyXZwzMbFTVp12eMTCzllTd8/GMgZm1pNLg4xkDM2tV5RcZesbAzFrh2yvMLAsHHzPLwsHHzLJQMUbc25b9fE2lL+KA2TOqrI4ZFYf42vo11db3xIOV1vfyS39TaX0A37roxErrO/jJOyqtb/p/e7vXhJ4i93zMLAsHH7M2kDQiaamk5ZIekPQxSbVS/tGSFkt6RNL9kr6dLsYdGF7Px6w91kfEkQCSXgx8DdgdOD/d13gd8L6I+EEqcyxwILAsU3s7zsHHrM0i4pfpqv97JF0AnAVc0Qg8qcySXO3LxaddZh0QEY8DQ8CLgUOB+/O2KD8HH7PMJN0t6WFJf5e7LZ3k4GPWAZLmAiPALylWf5jXyIuIYyhW85w9/tH9qe3Bx6P+NugkvQj4AvDZdPP15yhuon59qdjMLI3LqBMDzh71t0G0i6SlwHRgM/AVipUbiIj/J+lU4K8l7U/RG3oauDBXY3Po6GyXR/1tUETE0CT5dwFv6FBzulLHx3w86m9m0GUDzoM66m82iDoefDzqb2bQ4eDjUX8za+jEgLNH/c1sK20PPh71N7PxdNWAs5kNDgcfM8vCwcfMsuiLNZxveugXlb6I35v7giqrY6fHfjB5oSkYeemRlda3aXq1E4zTRzZUWh/A/129sdL65u0zq9L6Vr7rpErre93td/T9mtDu+ZhZFg4+Zj1A0tpx0i6Q9PGmtFWS9upcy7afg4+ZZeHgY2ZZeAF5s952tqTTSvv7ZWvJFDn4mPW2v42IRY0dSasytmVKfNplZll0+q72vhuxN7Pt49Mus94wU9KTpf1LsrWkIg4+Zj0gIlo6S4mI4TY3pTLdEnx6dsTezLZPtwSfnh2xN7Pt49kuM8vCwcfMsuj0aVffjdib2fbp9DeW9t2IvZltH592mVkWDj5mloWDj5ll0S3X+eyQk/cZqbS+er3a9YKfu3dxpfXNeNlRldZ333FvrLS++Yu/V2l9AG+cU+1/1T+58eFK6/vo7J0qre+3V5xf6brku53+qa5bE9o9HzPLwsHHbABJCklXlfanSfqVpG+l/TPS/o8lPSLpO5JeX2UbHHzMBtM64DBJu6T9twCrm8pcGxGviYiDgL8CbpT0qqoa4OBjNrhuAU5O2wuBqycqGBG3A18E3l/Vkzv4mA2ua4AFknYGjgDunqT8/cArq3pyBx+zARURDwLDFL2eW1o4pNIZs76Yajez7fYvwCLgOOCFk5R9DVDZNQpt7/l0w6i6mU3oMuBTEbFsW4UkvYFivOefqnriTvR8RkfVI2I9E4+qnwUg6XiKUfXjI6LaK8HMbIyIeBK4dILsUyUdC8wEVgKnVPk72akxn6yj6mY2VkTMGidtcUS8LW1fHhEviogjI+IVEXFCRHy/yjZ0KvhkHVU3s+7TkeCTe1TdzLpPJ2e7so2qm1n36WTwuQx4NiKWSTpuokKlUfXjO9UwM+u8jgWfnKPqZtZ92h58JhpVBxan7cuBy9vdDjPrLr69wsyycPAxsywcfMwsi764sXTot09VWt+q+u6V1rff9Grf5jUj1db3uttuqrQ+/fYXldYHcPcpf1hpfbO+/+TkhabgpWuWV1rfUxuqXZd8+nPrKl0TeueZu+7wtXju+ZhZFg4+ZpaFg4+ZZeHgY2ZZOPiYWRYOPmaWhYOPmWWxzeAj6XZJJzSlfUTSrZIeSvszJX1V0jJJD0laImmr+7nMzMom6/lcDSxoSlsAXFza/zDwVEQcHhGHAX8MbKquiWbWjyYLPtcDJ0uaASBpGNgPeKJUZl9KC8JHxIqIeD6Vv1nSfZKWSxpdk1nSWkmfSen/JuloSYslPS7p7dW8NDPrZtsMPhHxG+BHwIkpaQFwHVC+VPsy4BxJP5R0kaSDSnl/FBHzgaOAD0lqrGC4K/C9iDgU+C1wEcW3WrwLuHAHX5OZ9YBWBpzLp14LaPrmiYhYCswFPgPsCdxT+jL5D0l6ALgLmAM0AtNG4La0vQy4IyI2pe3h7XolZtZTWrlD8RvA30qaB8yMiPvS6deoiFgL3EjxfVt14CRJewNvBl4XEc9JWgzsnA7ZFBGN3lMdeD7VU5fUFze7mtm2TdrzSYHldorTq62+b0vS70jaI23PAA4BfgbMBp5JgeeVwGurbLiZ9bZWexlXAzex9cwXwIHA5yWJIph9G7gBmAF8QNLDwAqKUy8zM6DF4BMRN1P6Lq2IWAUclravBK4c57Dn2TJQ3VzfrNL2BRPlmVn/8hXOZpaFg4+ZZeHgY2ZZaMuMt5lZ57jnY2ZZOPiYWRYOPmaWhYOPmWXh4GNmWTj4mFkWfXEHeW32nGDzBkAgFfeBqLSd0knpxW1oSv+2pBflS2UhlR0tXt4oVVn65timckqbjUyVm1IqU6SPtqbYb9SkcfZHn0Jb7TeestT00ntROr65DU3tH7s/5mVsed/K6aXn3PK8Y8tsqbt8iUfaHr3sI1JSqUw5L+1rTNlx6okgynVFOX/89LHlm+qMSEWjlLzl+NHLVmJselF87PONqaeUXzStfBxbypeOGZM+2r7yU0T5HaFcbZ0YPWQ0r2k/tYimZhLA02z8TkS8lR3UF8GHzRuYccgpqDZUPIaGqNWGqE2bgWq1LemNvGkzqJXT0qM2fcbY/ZpQrQhWQ9NqqCZqSmk1UauJ2pC2lGvKHxqqoRpj6qnVxLRpNYZqGvOYVhMzptUYqjXlqfi50zjHDI0eM7Z8TWJIMH2oNrpdq02QnuqvCabXagzVUlqjntK+BLUUCGsp6DXqKNKKANooM9TYpnFcsT8kUSMg6ijqEPUimNQ3p+30qJfy6yOlsvUtx46MNNVTlI16HTZvJEZGoN5IK7Zj88bi58gI1EfGlk9lqY+MOba+eRMxUi8e9fro9simsen1tF3ftHk0vT5avnie+sbNRbnR9KA+UmdkU3k/0rFBjAQjG0dGy0V9S/7IxsbzFuVipNge2VxnJIKRYKufG+vbSg/qTJw3EvCP/GyvKn5tfdplZlk4+JhZFg4+ZpaFg4+ZZeHgY2ZZOPiYWRYOPmaWhYOPmWXh4GNmWTj4mFkWDj5mloWDj5ll4eBjZlk4+JhZFg4+ZpaFg4+ZZeHgY2ZZOPiYWRYOPmaWhYOPmWXh4GNmWTj4mFkWGv2uoR4m6SFgQ8XV7gU8XXGdNpbf4960b0Tsv6OV9Mf3dsGGiDiqygol3Vt1nTaW3+PeJGldFfX4tMvMsnDwMbMs+iX4fLFH6rSx/B73phurqKQvBpzNrPf0S8/HzHpM1wcfFS6V9KikByXNm6DcX0p6QtLapvSdJF2bjr9b0vAU6pwvaVkqd6kkpfQjJd0laamkeyUdXf0r722S3ippRXrv/myc/K0+l1LeuSl9haQTOtnuQSfpPEkb0+PWcfLPkvScpJC0qCnvi6VjJz+ljoiufgAnAbcCAl4L3D1BudcC+wJrm9I/CHwhbS8Arp1CnT9K+UrlT0zp3y1tnwQszv0+ddMDGAIeA+YCM4AHgEMm+1zS9iGp/E7Ay1I9Q7lf0yA8gOnAJuANwK7AeuD3m8ocC5ySPpdFpfS56di5wHDaHt7W83V9zwd4B3BlFO4CXiBp3+ZCEXFXRPxiguOvSNvXA29qpc60v3uqN4ArgXc2ng7YPW3PBn6+Yy+x7xwNPBoRj0fERuAaive8bKvPJfUs3wFcExHPR8RK4NFUn7XfGcCaiLgjItYBdwBnlgtExJKIuIHid6DsbOCR9JmvAh4BPratJ+uF4LM/8ERp/8mUNuXjI2IzsIYiMk9W5/4pfbwyHwE+I+kJYBFw7hTaMwha+czG+1xe2OKx1h4HM/aK81XAfi0eOwysLu3/PKVNqBeCTzf6U+DsiJhDEfG/nLk9Zj2nK4OPpDPTYO5S4BfAnFL2SxgbYSezGpjTqBM4APjPFupcndLHK3M6W651+Do+LWi2mtbe3zkAkqZRnL7+usVjrT1WUNxv1zBM60MKqxjbQ90vpU0s9yBXC4NgJzN2cPhHk5RvHnA+k7EDm9e1WidbDziflNIfBo5L228C7sv9PnXTg+KewccpBowbA86HTva5pO1DGTvg/DgecO7U57YTxUDx77JlwPntE5R9lPEHnIfZMuA8d5vPl/sFt/CGCPgcxej6MuCoUt7S0vbfUIwP1NPPC1L6zhS9k0dTMJk7hTqPAh5K5T7LlosyjwXuS78kdwPzc79P3fagmAX8aXrvzktpFzb+M4/3uZSOPS8dt4I0q+hHxz63TwIbU/D4bkq7A7gobZ8ObKYYcK5T3NTdOPaydOxG4MuTPZevcDazLLpyzMfM+p+Dj5ll4eBjZlk4+JhZFg4+ZpaFg4+ZZeHgY2ZZOPiYWRb/H35kQrPN/CN2AAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 360x360 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 320
        },
        "id": "u6lewrFqwEJ9",
        "outputId": "e1809e1c-cecf-4b79-8c73-59c58b7e8ce1"
      },
      "source": [
        "#@title Plotting Error Granger \n",
        "fig02 = drawmatrix_channels(np.transpose(g2_incorr), labels, size=[5., 5.], color_anchor=0)"
      ],
      "execution_count": 68,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAR8AAAEvCAYAAAB8LcXgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAaPklEQVR4nO3deZRc5X3m8e9TLbXQAmIzq8EyCvsSFh/AhBiI7cEsjtcEyYfYnHjGJ2N7sJ2QIR4mxuaQMDkmeMIYmzgZgjEYMJuNzeIkxxKDYoPZBEKWhVkPwgkgxGIJbV31mz/u263bpW6pJN2qt6vq+Rzq6C7vfestWnr63ve99y1FBGZmnVbL3QAz608OHzPLwuFjZlk4fMwsC4ePmWXh8DGzLBw+ZpaFw8fMsujZ8JFUl7RQ0mJJj0r6M0m10v5jJc2X9CtJD0u6Q9LhOdts1k8m5W5AG62OiCMBJO0GfBfYAbhQ0u7A94CPRcRPU5kTgdnAokztNesr6tXHKyStjIgZpfX9gAeAXYGLgEZEXJirfWb9rmcvu5pFxNPAALAbcCjwcN4WmfW3vgmfTZF0v6Qlkv4ud1vM+kXfhE+67KoDLwGLgaOH90XEccBfAjPztM6s//RF+Eh6C3Al8PUoOrmuAM6RdEKp2LQsjTPrU73c4VynGLmaDAwB3wEui4hG2n888DfA3hRnQ8uBiyLiwTwtNusvPRs+Zjax9cVll5lNPA4fM8vC4WNmWTh8zCwLh4+ZZeHwMbMsHD5mloXDx8yycPiYWRY9MZnY0LIlld6mPemtB6vK+sxsYz7zMbMsHD4tkrRynO1nS3qsNFf0P0racYxyJ0v6UftbatYdeuKyKxdJ7wO+AJwWES9IGgA+AewOvJa1cWYTnMNn21wAnBcRLwBERB24anhnCqf/DbwJLMjSQrMJypdd22bcuaAlbQf8A/B+4Bhgjw62y2zCc/hURNLh6XvCnpJ0FnAQ8ExE/CrNnnht5iaaTSgOn20zMhd0RCxK3xN2FzA1a6vMukDl4dNno0KXAJdKemtp23Dw/BKYJWl2Wp/b0ZaZTXAd6XDukVGhaZKWldYvi4jL0uT0d6XP9BrwOPDjiFgj6VPAHZLeBO4Ftu98s80mpsrncG7+ptC07V7gSxExb5xjmkeF9ouIM1t9T9/hbNZ9OtXn41EhMxul4x3OHhUyM+hc+HhUyMxG6VT4eFTIzEZpx2iXR4XMbLN64htLPdpl1n18h7OZZeHwMbMsHD5mlkVP9Pn8+tWVlX6I3VcsrrI6BmYf5z4ksyY+8zGzLBw+GUm6ID3l/1i66/s4Sc9K2rVUZuQpf0nnSHo5lV0s6WZJ0/J9ArOt5/DJRNI7gTOBoyPiCOA9wPMtHHpjRBwZEYcC64Cz2thMs7bxHM757Aksj4i1ABGxHEBqrXtI0iRgOvBquxpo1k4+88nnn4F9JD0h6RuSTmrxuLMkLQReAHYGfti2Fpq1kcMnk4hYSTGFyKeAl4EbJZ0DjDVyV952Y3owdw9gEfDnbW6qWVs4fDKKiHpEzI+IC4HPAh8BXgF2KhXbGVg+xrFBcdbzrk601axq7ZjD2SM4LZB0oKT9S5uOBJ4D5gN/lMoMAGcDY84ACZwIPNXGZpq1TaUdzk0jOGtT4Ay2cOiNEfHZVMd3KUZw/qnKtk1AM4D/kybRHwKepLgEWw98U9KjgIC7GT3B2lmSTqT4xbEMOKeTjTarStWjXR7BaVFEPAScMM7uj41zzNXA1W1qkllHVX3Z5REcM2tJpeHjERwza1XlHc4ewTGzVlQaPh7BMbNWVd3h7BEcM2uJ5/MZg+fzMWs/3+FsZlk4fMwsC4ePmWXRE/P5PLFiTaX11Xc6tNL6dl29utI+qalTp7oPybqez3zMLAuHTw+QNE/SqU3bPi/pLkmPp/Vpkq6TtEjS45IWSJqRp8VmDp9ecT0wp2nbHOCS0vrngBcj4vCIOAz4JMX9V2ZZOHx6w83AGZIGASTNAvZi9IT0e1I8uAtARCxN057MkvTLdFa0pF/mU7L8HD49ICJWAD8HTkub5gDfY/TDu1cB50v6maSLmx6DORD4RkQcDLwBfLoDzbY+5/DpHeVLrzlpfURELAT2A75K8WDvA5IOTrufj4h/S8vXUjxfZ9ZWDp/e8QPg3ZKOBqalycpGiYiVEXFrRHyaImROH97VXLS9TTWrKHw82pJfmktpHsXl1fXN+yX9jqSd0vIgcAjFjAMA+6YpcKGYRXFB+1ts/a6qMx+PtkwM1wO/zRjhA8wG7pG0CHgEeBC4Je1bCnxG0hKKeZe+2YG2Wp+r6g7nm4GLJQ1GxLpNjLYM/6YlIpbCyMjM3cBDwNHAYuDjEfFmRW3rGxHxfYopS4bXnwUOS8vXANc0H5Pm1x6KiLM700qzQiVnPh5tMbMtVWWHs0dbulBEPJsug806qsrw8WiLmbWssvDxaIuZbYmq7/PxaIuZtST7HM5ptOtH29LvMP+p5ZV+iNk7bVdldew6daDS+jyfj/UC3+FsZllkn8mwfC+KmfUPn/mYWRbZ+3yqsPb1Vyr9EDFQ7Qnhc6uq7aJ52w6TK61vO/chWQY+8zGzLBw+Ni5JH5QUkg4qbTtW0v+TtFTSI5L+sTzzoaTvS7ovT4utmzh8bFPmUtzwORdA0u7ATcD5EXFgRBxF8VDw9mn/jsAxwExJ++VpsnULh4+NKc21dCLF1CfDz+x9Bvh2RPxsuFxE3BwRL6bVDwM/BG5g4ylWzEZx+Nh4PgDcHRFPAK9IOobiloiNntkrmUtxd/v1adlsXA4fG89cijMY0p+bDJN0SbY/sCAF1npJvn/LxpX9JkObeCTtDPwecLikAAYoZhr4NkWfzg/GOOwPKZ7LeyZNULYDRWBd0Ik2W/dpy5mPR0m63keB70TE2yJiVkTsAzwD/CvwCUnHDReU9OF01jMXeF8qP4sipNzvY+Nq12WXR0m621zgtqZtt1CEyRzg0vRLZAlwKrAL8DZg5JdHRDwDvF4OKrOyyu9wTqMkS4FTgB9GxIGSLgKIiC+Nc8wfA+8AXgTWR8Rfb8l7+g7nbeM7nC2Hdpz5eJTEzDarHeHjURIz26xKry88SmJmrar6zMejJGbWkqrDx6MkZtYSz+czBo92mbWfH68wsywcPmaWhcPHzLLojQdLG0OVVqdatZm889TBSusjGpVWt2bVbyrtM9tu+vbuQ7LN8pmPmWXh8LGOSLMcXFtanyTpZUk/Km37oKTHJC2RtEjSB/O01jqhNy67rBusAg6TNDUiVgPvBV4Y3inpt4FLgfdGxDOS3g78i6SnI+KxPE22dvKZj3XSncAZaXn4YeJh5wF/nW4yHb7Z9BLgzwEknSvpF+nM6Aas6zl8rJNuAOZI2g44Ari/tO9QNp754MG0HeAvgKMi4gjgT9rdUGs/h491TLp8mkVx1nPnFh7+GHCdpLOBaoc3LQuHj3Xa7RR9O9c3bf8FxUPFZccAi9PyGcAVwNHAA5LcX9nlKgsfj2ZYi64CvhIRi5q2Xwp8UdIsgPTn/wD+VlIN2Cci5gHnAzOBGR1qr7VJlb89PJphmxURy4DLx9i+UNL5wA8lTQbWA/89bZ8MXCtpJiDg8oh4raMNt8pV9lS7pJUUf6kejoibJV1Dccr8uxFxpqTvAPMi4qrSMZ8ETo6IP5J0LkVH4hDwi4hoeU6fta++WO2j+ZOqvSP5jai2vu2rvuCo+I5p3+Fsrai6z8ejGWbWkkrDx6MZZtaqdox2eTTDzDarHf/ArwJei4hFkk4ubb8UuEnSTyLi2dJoxkfLoxmSFlBMuToDcKeiWY+qPHw8mmFmreiNOZw92rVtPNplGfgOZzPLwuFjZlk4fMwsCw9nj0GrX6+0vinTd6+0vnqltcFAxXNWVz0nNLgfqRf5zMfMsnD4WFeSVJe0UNLjkm6SNC1tX5m7bdYah491q9URcWREHAasw88Ddh2Hj/WCe4Hfyt0I2zIOH+tq6RnA04DmyclsgvNol3WrqZIWpuV7gf+bszG25Rw+1q1WR8SRuRthW8+XXWaWRcfDx0Ok1mbTJC0rvf40d4NsbDkuu0ZOlyVdRzFEelmGdlgXi4gxv70iInw23yVy/6A8RGrWp7KFj4dIzfpbjssuD5GaWd4+HzPrX7n7fMysT02k8PEQqVkf6fhll4dIzQwm1pmPmfURh4+ZZeHwMbMseuKpdjWqnVJdzy/efKEtUDu42gnkaxVPpf75aQdXWt8X/svRldaXVDop/YHfutUT0mfmMx8zy8LhY5ZI2kPSDZKekvSQpDslHZD2fV7SGkkzS+V3kTRP0kpJX8/X8u7k8DEDJAm4DZgfEbMj4hjgi8DwNfNc4AHgw6XD1gB/CZzXybb2CoePWeEUYH1EXDm8ISIejYh7Jc0GZgD/kyKEhveviogFFCFkW8jhY1Y4DHhonH1zgBsoHoQ+UFK1Iwh9yuFjtnlzgRsiogHcAvxB5vb0hJ4YajerwGLgo80bJR0O7A/8S9EtxCDwDOAO5m1U+ZmPRwysS/0EmCLpU8MbJB0BXA58OSJmpddewF6S3parob2i0vDxiIF1q4gI4EPAe9IvzsXAJcDJFH+ny26j6AdC0rMUc5Cfk2ZjOKRjje5yVV92jTliAFAaMfg0cAHwT2n/KmCBJM/lbFlFxK+BP2yh3J+Wlme1s029rOrLLo8YmFlLOjna5REDMxtR9WWXRwzMrCVVn/l4xMDMWlJp+HjEwMxaVflNhh4xMLNW+PEKM8vC4WNmWTh8zCwLFX3E3e35FSsr/RBVz5G8YnW1c0wfxIuV1jd/5Y6V1jdzu+qfV378xZWV1nf6/rtUWt8eO073nNBbyGc+ZpaFw8esDSTVJS2UtFjSo5L+TFKttP9YSfMl/UrSw5LuSDfj9g3P52PWHqsj4kgASbsB3wV2AC5MzzV+D/hYRPw0lTkRmA0sytTejnP4mLVZRLyU7vp/QNKXgc8C3x4OnlRmQa725eLLLrMOiIingQFgN+BQ4OG8LcrP4WOWmaT7JS2R9He529JJDh+zDpC0H1AHXqKY/WHkO6Uj4jiK2Txnjn10b2p7+LjX3/qdpLcAVwJfTw9fX0HxEPUJpWLTsjQuo050OLvX3/rRVEkLgcnAEPAdipkbiIj/kHQW8DeS9qY4G1oOXJSrsTl0dLTLvf7WLyJiYDP77wNO6lBzJqSO9/m419/MYIJ1OPdrr79ZP+p4+LjX38ygw+HjXn8zG9aJDmf3+pvZRtoePu71N7OxTKgOZzPrHw4fM8vC4WNmWfTEfD6rhxqV1jdzyia7qbbYfcteq7S+Aw6sdv7hQ6dOqbS+HSr+/wdwzzMrKq3vpTeHKq1v58Z/VDqP+ODOe/T8nNA+8zGzLBw+Zl1A0kZf3yHpy5LOa9r2rKRdO9eyrefwMbMsHD5mlkVPdDib9bEvSDq7tL5XtpZsIYePWXf7WkRcOrwi6dmMbdkivuwysyw6/VR7z/XYm9nW8WWXWXeYJmlZaf2ybC2piMPHrAtEREtXKRExq81NqcxECZ+u7bE3s60zUcKna3vszWzreLTLzLJw+JhZFp2+7Oq5Hnsz2zqd/sbSnuuxN7Ot48suM8vC4WNmWTh8zCyLiXKfzzbZd3q1092ujmrre3N9vdL6qFX7O+PV/zan0vrWvmXHSusD+M9f+Val9e1x0rmV1vfqlw6vtD5237faOaGP/9CEmxPaZz5mloXDx6wPSQpJ15bWJ0l6WdKP0vo5af0RSb+S9GNJJ1TZBoePWX9aBRwmaWpafy/wQlOZGyPiqIjYH/hfwK2SDq6qAQ4fs/51J3BGWp4LXD9ewYiYB3wL+FRVb+7wMetfNwBzJG0HHAHcv5nyDwMHVfXmDh+zPhURjwGzKM567mzhkEpHzHpiqN3MttrtwKXAycDmvof7KGBJVW/c9jOfidCrbmbjugr4SkQs2lQhSSdR9Pf8Q1Vv3Ikzn5Fe9YhYzfi96p8FkHQKRa/6KRFRWcqa2cYiYhlw+Ti7z5J0IjANeAb4SJX/JjvV55O1V93MRouIGWNsmx8RZ6blqyPiLRFxZEQcEBGnRsS/VdmGToVP1l51M5t4OhI+uXvVzWzi6eRoV7ZedTObeDoZPlcBr0XEIkknj1eo1Kt+SqcaZmad17HwydmrbmYTT9vDZ7xedWB+Wr4auLrd7TCzicWPV5hZFg4fM8vC4WNmWSii0qlis1j7m1cr/RCrmFJldUxvrK60vpfq1bavXvHfgZXrGpXWB7DP9pMrrW/G8qWV1je0076V1heD0yqt7+XV1c4jvs/OM7b5Xjyf+ZhZFg4fM8vC4WNmWTh8zCwLh4+ZZeHwMbMsHD5mlsUmw0fSPEmnNm37vKS7JD2e1qdJuk7SIkmPS1ogaaPnuczMyjZ35nM9MKdp2xzgktL654AXI+LwiDgM+CSwvrommlkv2lz43AycIWkQQNIsYC/g+VKZPSlNCB8RSyNibSr/fUkPSVosaWROZkkrJX01bf9XScdKmi/paUm/X81HM7OJbJPhExErgJ8Dp6VNc4DvAeX78a8Czpf0M0kXS9q/tO+PI+IY4B3AuZKGZzCcDvwkIg4FfgNcTPGtFh8CLtrGz2RmXaCVDufypdccmr55IiIWAvsBXwV2Bh4ofZn8uZIeBe4D9gGGg2kdcHdaXgTcExHr0/KsrfokZtZVWplM7AfA1yQdDUyLiIfS5deIiFgJ3ErxfVsN4HRJuwPvAd4ZEW9Kmg9slw5ZHxueaG0Aa1M9DUn+FlWzPrDZM58ULPMoLq82+r4tSb8jaae0PAgcAjwHzAReTcFzEHB8lQ03s+7W6lnG9cBtbDzyBTAb+KYkUYTZHcAtwCDwJ5KWAEspLr3MzADP5zMmz+ezbTyfz7bzfD5mZm3i8DGzLBw+ZpZFT/T5mFn38ZmPmWXh8DGzLBw+ZpaFw8fMsnD4mFkWDh8zy6InniCvzdwnGFoDCCQEoNJy2k7aXjyGpvTfhu1F+VJZSGVHipcXSlWW7jRvKqe0OLxT5aaUyhTbR1pTrA/XpDHWR95CG60Pv2Wp6aX/F6Xjm9vQ1P7R66M+xob/b+Xtpffc8L6jy2you3yLR1oeue0j0qZSmfK+tK5RZceoJ4Io1xXl/WNvH12+qc6IVDRKmzccP3LbSozeXhQf/X6j6intL5pWPo4N5UvHjNo+0r7yW0T5/wjlahvEyCEj+5rWU4toaiYBLGfdjyPifWyjnggfhtYweMhHUG2geA0MUKsNUJs0iGq1DduH900apFbell61yYOj12tCtSKsBibVUE3UlLbVRK0magPaUK5p/8BADdUYVU+tJiZNqjFQ06jXpJoYnFRjoNa0T8WfU8Y4ZmDkmNHlaxIDgskDtZHlWm2c7an+mmByrcZALW0brqe0LkEtBWEthd5wHcW2IkCHywwMLzN8XLE+IFEjIBooGhCNIkwaQ2k5vRql/Y16qWxjw7H1elM9RdloNGBoHVGvQ2N4W7EcQ+uKP+t1aNRHl09ladRHHdsYWk/UG8Wr0RhZrq8fvb2Rlhvrh0a2N0bKF+/TWDdUlBvZHjTqDerry+uRjg2iHtTX1UfKRWPD/vq64fctykW9WK4PNahHUA82+nNdY1Pbgwbj76sH/D3P7VrFP1tfdplZFg4fM8vC4WNmWTh8zCwLh4+ZZeHwMbMsHD5mloXDx8yycPiYWRYOHzPLwuFjZlk4fMwsC4ePmWXh8DGzLBw+ZpaFw8fMsnD4mFkWDh8zy8LhY2ZZOHzMLAuHj5ll4fAxsyw08l1DXUzS48AMYHnuttgW2RX/zLrRnhGx97ZW0hvf2wVrgDUR8Y7cDbHWSXrQP7PuI2lVFfX4ssvMsnD4mFkWvXLZ9a3cDbCt4p9bd7q1ikp6osPZzLqPL7vMLIsJHz4qXC7pNUlrJT0h6egxyt0t6Q1JayStkPT3kgbSvjslrZe0WtKLkk7v/CfpL5LeJ2mppCcl/cUY+6dIujHtv1/SrNK+L6btSyWd2sl29ztJF0hal153jbH/uvRvbLWk1yW9v7TvI5J+k/avkTRzk28WERP6BZwOPAjcBRwPLAbuH6PcDukl4BbgfmAOcArwNHB+Krdb7s/U6y9gAHgK2A8YBB4FDmkq82ngyrQ8B7gxLR+Syk8B3p7qGcj9mfrhBUwG1gMnAdOB1cD7m8rsXVq+GHg5LU9J5f8grf8WMHlT7zfhz3yADwCvAddExH0UneS7SNqzXCgi3oiIN9L+KenPAP4rsACop3IvdbDt/epY4MmIeDoi1gE3UPwcyz4AfDst3wy8W5LS9hsiYm1EPAM8meqz9jsHeD0i7omIVcA9wGfKBSLihdJq+czmfODfI+KmVO7JiFi/qTfrhvDZG5gKPJ/WlwEr0vZRJP0YeBP4T8ATFH+pDwD2Bf5K0kpJt0vaqRMN72N7s+HnBcXPrPnnNVImIoaA14FdWjzW2uNARt9x/iywV3OhdLm8nuIX+9y0+ai0b7mkNyXdsbk364bwaVlEnApsD9wBvBX4PYozoCeBacC7gXcBf5urjWbdLiLOiojJwJXpBcW/s32AE9Kf75J03qbqmZDhI+kzkpZJWg28E1hL8YGgCJWdgRfGOjYi1gC3UTxy8QGK35zXRUQ9Iu6n+A37zjZ/hH73Aht+XlD8zJp/XiNlJE2iOIV/pcVjrT2WUjxvN2wW8OtNlP8cMDstPw08FxFPRMQrFH2uJ2/qzSZk+ETEFRHx1oiYCnycoiP545KOB4aAFRHx78PlJc2QNFvSnukv8pkUAfVL4PvA+1O5A1Jdj3b2E/WdB4D9Jb1d0iBFh/LtTWVuBz6Rlj8K/CSKnsrbgTlpNOztwP7AzzvU7n53DTBT0u9Kmk7R8fyNcgFJ7ymtXkjRzQHwNWAvSbtImkJxGfbgJt8tdw97Cz3wAq6gOGNZS3EJ9Y60b2F67Q48kv5HrKH4Dfp1ilPBQYpUXkPRG/9Tiqdys3+2Xn5RjFI+QTFadUHadhHw+2l5O+Cm9PP8ObBf6dgL0nFLgdNyf5Z+egFfAtZRjHr9c9p2D3BxWn6k9G/p1eGfZ9r3jbRvDWOMSDe/fIezmWUxIS+7zKz3OXzMLAuHj5ll4fAxsywcPmaWhcPHzLJw+JhZFg4fM8vi/wPcELaqKiFdNQAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<Figure size 360x360 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 320
        },
        "id": "mZm_kEE_wLJK",
        "outputId": "19191b19-7a09-4cc9-bc49-ca4344ca3962"
      },
      "source": [
        "#@title Plotting Difference Granger \n",
        "g2_diff = np.subtract(g2_corr,g2_incorr)\n",
        "fig03 = drawmatrix_channels(np.transpose(g2_diff), labels, size=[5., 5.], color_anchor=0)"
      ],
      "execution_count": 70,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAR8AAAEvCAYAAAB8LcXgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAZ+klEQVR4nO3dfbRdVXnv8e9vn+SQhGDkRUAoNZCLIG9CcEBVuMBVi4hWq7YkDrwyyhhcL3h9qfai5SqWwSh1SGnL9YWrLUUBCYhvKC9qhwklLSIQAkkIUSA4CAwkIG8JCSfn7Of+seY+WefknGQnWXvP/fL7DPbIWnPNNfezOTlP5ppzrbkVEZiZtVstdwBm1p+cfMwsCycfM8vCycfMsnDyMbMsnHzMLAsnHzPLwsnHzLLo2eQjaUTSUkkrJN0v6dOSaqXjx0laJOk3kpZIulnSkTljNusnU3IH0EIbIuJoAEl7A98BXgVcKGkf4AbgQxHxn6nOCcAcYFmmeM36inr18QpJ6yJiZmn/IOBuYC/gIqAeERfmis+s3/XsZdd4EfEoMADsDRwOLMkbkVl/65vkszWS7pK0UtI/5Y7FrF/0TfJJl10jwNPACmBu41hEHA98HpiVJzqz/tMXyUfSa4ArgK9EMcj1VeAsSW8pVZuRJTizPtXLA84jFDNXU4Fh4Grgsoiop+N/BHwJ2J+iN/QMcFFE3JMnYrP+0rPJx8w6W19cdplZ53HyMbMsnHzMLAsnHzPLwsnHzLJw8jGzLJx8zCwLJx8zy8LJx8yy6InFxB586oVKb9M+bN9ZqrI9M9uSez5mloWTT5MkrZuk/ExJD5TWiv5nSa+eoN7Jkn7S+kjNukNPXHblIumdwKeA0yLiCUkDwEeAfYDnswZn1uGcfHbOBcBnIuIJgIgYAa5sHEzJ6R+Bl4HFWSI061C+7No5k64FLWka8E3gPcCxwL5tjMus4zn5VETSkel7wh6RdAZwKLA6In6TVk+8JnOIZh3FyWfnjK4FHRHL0veE3QpMzxqVWReoPPn02azQJcClkv6gVNZIPA8BsyXNSfvz2xqZWYdry4Bzj8wKzZC0prR/WURclhanvzV9pueB5cBPI2KjpHOAmyW9DNwB7Nb+sM06U+VrOI//ptBUdgfwhYhYOMk542eFDoqIdzf7nr7D2az7tGvMx7NCZjZG2wecPStkZtC+5ONZITMbo13Jx7NCZjZGK2a7PCtkZtvUE99Y6tkus+7jO5zNLAsnHzPLwsnHzLLoifV8RLVDNFcvebzSMaQPzz3AY0hm47jnY2ZZOPlkJOmC9JT/A+mu7+MlPSZpr1Kd0af8JZ0laW2qu0LSjZJm5PsEZjvOyScTSW8G3g3MjYijgLcDjzdx6vURcXREHA4MAWe0MEyzlumJMZ8u9VrgmYh4BSAingGQmhsekjQF2BV4rlUBmrWSez75/Aw4QNKvJX1N0klNnneGpKXAE8AewI9bFqFZCzn5ZBIR6yiWEDkHWAtcL+ksYKKZtnLZ9enB3H2BZcBftThUs5Zw8skoIkYiYlFEXAh8DPgA8Cywe6naHsAzE5wbFL2e/9qOWM2q1oo1nD2D0wRJh0g6uFR0NPBbYBHw4VRnADgTmHAFSOAE4JEWhmnWMpUOOI+bwXklJZzBJk69PiI+ltr4DsUMzr9WGVsHmgn837SI/jDwMMUl2Cbg65LuBwTcxtgF1s6QdALFPxxrgLPaGbRZVaqe7fIMTpMi4l7gLZMc/tAk51wFXNWikMzaqurLLs/gmFlTKk0+nsExs2ZVPuDsGRwza0alycczOGbWrKoHnD2DY2ZN6Yk1nFc+9WKlH+KeJ1+osjmv52M2Ad/hbGZZOPmYWRZOPmaWRU+M+QyvWVnph6hPq/Y7C6tub9pMf6+YdT/3fMwsCyefHiBpoaRTx5V9UtKtkpan/RmSrpW0TNJySYslzcwTsZmTT6+4Dpg3rmwecElp/xPA7yLiyIg4Ajib4v4rsyycfHrDjcDpkgYBJM0G9mPsgvSvpXhwF4CIWJWWPZkt6aHUK1rZL+spWX5OPj0gIn4P/Ao4LRXNA25g7MO7VwLnS7pT0sXjHoM5BPhaRLwBeBE4tw1hW59z8ukd5UuveWl/VEQsBQ4CvkzxYO/dkt6QDj8eEf+Rtq+heL7OrKWcfHrHj4C3SZoLzEiLlY0REesi4vsRcS5FknlX49D4qq0N1ayi5OPZlvzSWkoLKS6vrht/XNJbJe2etgeBwyhWHAD4w7QELhSrKC5ufcTW76rq+Xi2pTNcB7yRCZIPMAe4XdIy4D7gHuB76dgq4DxJKynWXfp6G2K1PlfVkho3AhdLGoyIoa3MtjT+pSUiVsHozMxtwL3AXGAF8N8j4uWKYusbEfFDiiVLGvuPAUek7W8D3x5/TlpfezgizmxPlGaFSno+nm0xs+1V5YCzZ1u6UEQ8li6DzdqqyuTj2RYza1plycezLWa2Paq+z8ezLWbWlOzr+aTZrp/szLiD1/Mx6z6+w9nMsqj6q3O2W/leFDPrH+75mFkW2cd8qrD2xfWVfojdKu4P1pbcXGl79bmnV9retBm7egzJ2s49HzPLwsnHJiXpfZJC0qGlsuMk/bukVZLuk/TP5ZUPJf1Q0i/zRGzdxMnHtmY+xQ2f8wEk7QN8Fzg/Ig6JiGMoHgreLR1/NXAsMEvSQXlCtm7h5GMTSmstnUCx9Enjmb3zgG9FxJ2NehFxY0T8Lu2+H/gxsIAtl1gxG8PJxybzXuC2iPg18KykYyluidjimb2S+RR3t1+Xts0m5eRjk5lP0YMh/bnVZJIuyQ4GFqeEtUmS79+ySWW/ydA6j6Q9gP8GHCkpgAGKlQa+RTGm86MJTvtziufyVqcFyl5FkbAuaEfM1n1a0vPxLEnX+yBwdUS8LiJmR8QBwGrg34CPSDq+UVHS+1OvZz7wzlR/NkWS8riPTapVl12eJelu84EfjCv7HkUymQdcmv4RWQmcCuwJvA4Y/ccjIlYDL5QTlVlZ5Xc4p1mSVcApwI8j4hBJFwFExBcmOecvgDcBvwM2RcTfbs97+g7nneM7nC2HVvR8PEtiZtvUiuTjWRIz26ZKLzA8S2Jmzaq65+NZEjNrStXJx7MkZtaUSi+7IuKUCcouL+2eOMFp+09wztwq4zKzzuPHK8wsCycfM8vCycfMsuiJNZzXvbyh0g8x/emHqmyOULU5fs3Map9A2XdGtbd0T5s+3XdM2za552NmWTj5WFukVQ6uKe1PkbRW0k9KZe+T9ICklZKWSXpfnmitHbyej7XLeuAISdMjYgPwDuCJxkFJbwQuBd4REaslHQj8XNKjEfFAnpCtldzzsXa6BWg8kt94mLjhM8DfpptMGzebXgL8FYCkj0t6MPWMFmBdz8nH2mkBME/SNOAo4K7SscPZcuWDe1I5wGeBYyLiKOCjrQ7UWs/Jx9omXT7Npuj13LKdpz8AXCvpTGC44tAsAycfa7ebKMZ2rhtX/iDFQ8VlxwIr0vbpwFeBucDdkjxe2eUqSz6ezbAmXQn8TUQsG1d+KfA5SbMB0p9/Dfy9pBpwQEQsBM4HZgEz2xSvtUiV/3p4NsO2KSLWAJdPUL5U0vnAjyVNBTYB/zuVTwWukTQLEHB5RDzf1sCtcpXd4SxpHcVfqiURcaOkb1N0mU+MiHdLuhpYGBFXls45Gzg5Ij4s6eMUA4nDwIMR0fSaPr7Deef4DmfLoeoxH89mmFlTKk0+ns0ws2a1YrbLsxlmtk2t+AW/Eng+IpZJOrlUfinwXUm/iIjHSrMZHyzPZkhaTLHk6kzAg4pmPary5OPZDDNrhtfzmYBnu3aOZ7usGb7D2cyycPIxsyycfMwsi54Y83nmxfUd/SFeNfRcpe3FLrtW2t7wwLRK22vFD2PmDI8j9Rr3fMwsCycf60qSRiQtlbRc0nclzUjl63LHZs1x8rFutSEijo6II4Ah/Dxg13HysV5wB/Bfcgdh28fJx7paegbwNGD84mTW4fzwpnWr6ZKWpu07gH/JGYxtPycf61YbIuLo3EHYjvNll5ll0fbk4ylSa7EZktaUXn+ZOyCbWI7LrtHusqRrKaZIL8sQh3WxiJjw2ysiwr35LpH7B+UpUrM+lS35eIrUrL/luOzyFKmZ5R3zMbP+lXvMx8z6VCclH0+RmvWRtl92eYrUzKCzej5m1kecfMwsCycfM8uiJ55qH65X296Sp6p9zOxt+0yttL3n64OVtveHJ/+vStv743POrrQ9gLe+fq9K16X/7CkHe0H6zNzzMbMsnHzMEkn7Slog6RFJ90q6RdLr07FPStooaVap/p6SFkpaJ+kr+SLvTk4+ZoAkAT8AFkXEnIg4FvgcsE+qMh+4G3h/6bSNwOeBz7Qz1l7h5GNWOAXYFBFXNAoi4v6IuEPSHGAm8H8oklDj+PqIWEyRhGw7OfmYFY4A7p3k2DxgAcWD0IdI2meSerYdnHzMtm0+sCAi6sD3gD/LHE9P6ImpdrMKrAA+OL5Q0pHAwcDPi2EhBoHVgAeYd1LlPR/PGFiX+gWwi6RzGgWSjgIuB74YEbPTaz9gP0mvyxVor6g0+XjGwLpVRATwp8Db0z+cK4BLgJMp/k6X/YBiHAhJj1GsQX5WWo3hsLYF3eWqvuyacMYAoDRjcC5wAfCv6fh6YLEkr+VsWUXEk8CfN1HvL0vbs1sZUy+r+rLLMwZm1pR2znZ5xsDMRlV92eUZAzNrStU9H88YmFlTKk0+njEws2ZVfpOhZwzMrBl+vMLMsnDyMbMsnHzMLIueeLB096h2zeV37L6+0va0ZnWl7W3Y77hK23v2UwdU2h7vr/7bsHXz5ZW2N/LIyZWuCT0w53ivCb2d3PMxsyycfMxaQNKIpKWSVki6X9KnJdVKx4+TtEjSbyQtkXRzuhm3b/TEZZdZB9oQEUcDSNob+A7wKuDC9FzjDcCHIuI/U50TgDnAskzxtp2Tj1mLRcTT6a7/uyV9EfgY8K1G4kl1FueKLxdfdpm1QUQ8CgwAewOHA0vyRpSfk49ZZpLukrRS0j/ljqWdnHzM2kDSQcAI8DTF6g9zG8ci4niK1TxnTXx2b2p58vGov/U7Sa8BrgC+kh6+/irFQ9RvKVWbkSW4jNox4OxRf+tH0yUtBaYCw8DVFCs3EBFPSToD+JKk/Sl6Q88AF+UKNoe2znZ51N/6RUQMbOP4L4GT2hROR2r7mI9H/c0MOmzAuV9H/c36UduTj0f9zQzanHw86m9mDe0YcPaov5ltoeXJx6P+ZjaRjhpwNrP+4eRjZlk4+ZhZFj2xno/qI5W2FwODlbb35HeurrS93S94c6XtDf7xRyptL0Y2VtoewKZXKm5zeLja5n57f6VrQk953Rt7fk1o93zMLAsnH7MuIGmLr2iR9EVJnxlX9pikvdoX2Y5z8jGzLJx8zCyLnhhwNutjn5J0Zml/v2yRbCcnH7Pu9g8RcWljR9JjGWPZLr7sMrMs2v1Ue8+N2JvZjvFll1l3mCFpTWn/smyRVMTJx6wLRERTVykRMbvFoVSmU5JP147Ym9mO6ZTk07Uj9ma2YzzbZWZZOPmYWRbtvuzquRF7M9sx7f7G0p4bsTezHePLLjPLwsnHzLJw8jGzLDrlPp+dUp+xe6XtaeOLlbb3/MNPVdreqyttDT5750uVtnfn8ocrbQ9g0bkfrbS98/Y4vtL2/sfjSyttb5enXqh0TejD9p3VcWtCu+djZlk4+Zj1IUkh6ZrS/hRJayX9JO2flfbvk/QbST+V9JYqY3DyMetP64EjJE1P++8AnhhX5/qIOCYiDgb+Dvi+pDdUFYCTj1n/ugU4PW3PB66brGJELAS+AZxT1Zs7+Zj1rwXAPEnTgKOAu7ZRfwlwaFVv7uRj1qci4gFgNkWv55YmTql0xqwnptrNbIfdBFwKnAzsuY26xwArq3rjlvd8OmFU3cwmdSXwNxGxbGuVJJ1EMd7zzareuB09n9FR9YjYwOSj6h8DkHQKxaj6KRFRWZY1sy1FxBrg8kkOnyHpBGAGsBr4QJW/k+0a88k6qm5mY0XEzAnKFkXEu9P2VRHxmog4OiJeHxGnRsR/VBlDu5JP1lF1M+s8bUk+uUfVzazztHO2K9uoupl1nnYmnyuB5yNimaSTJ6tUGlU/pV2BmVn7tS355BxVN7PO0/LkM9moOrAobV8FXNXqOMyss/jxCjPLwsnHzLJw8jGzLBRR6VKxWWxc/1K1H6JW7VBYbd3aStsj6pU2N/Dyc5W2V1/7eKXtAdQ3rK+0vYEDj6q0vedmza60vapvdJv1UrU/k6n7H7rTIbrnY2ZZOPmYWRZOPmaWhZOPmWXh5GNmWTj5mFkWTj5mlsVWk4+khZJOHVf2SUm3Slqe9mdIulbSMknLJS2WtMXzXGZmZdvq+VwHzBtXNg+4pLT/CeB3EXFkRBwBnA1sqi5EM+tF20o+NwKnSxoEkDQb2A8o3y75WkoLwkfEqoh4JdX/oaR7Ja2QNLoms6R1kr6cyv9N0nGSFkl6VNKfVPPRzKyTbTX5RMTvgV8Bp6WiecANQPlxhiuB8yXdKeliSQeXjv1FRBwLvAn4uKTGCoa7Ar+IiMOBl4CLKb7V4k+Bi3byM5lZF2hmwLl86TWPcd88ERFLgYOALwN7AHeXvkz+45LuB34JHAA0EtMQcFvaXgbcHhGb0vbsHfokZtZVmnmC8kfAP0iaC8yIiHvT5deoiFgHfJ/i+7bqwLsk7QO8HXhzRLwsaREwLZ2yKTY/0VoHXknt1CX5W1TN+sA2ez4psSykuLza4vu2JL1V0u5pexA4DPgtMAt4LiWeQ4E/qjJwM+tuzfYyrgN+wJYzXwBzgK9LEkUyuxn4HjAIfFTSSmAVxaWXmRng9Xwm5vV8dorX89l5Xs/HzKxFnHzMLAsnHzPLoifGfMys+7jnY2ZZOPmYWRZOPmaWhZOPmWXh5GNmWTj5mFkWPfEEeW3WAcHwRkAgFbemq7SdyknlxWNoSv9tLi/ql+pCqjtavbxRarJ0p/m4ekqbjYMqh1KqU5SPRlPsN1rSBPujb6Et9htvWQq99P+idP74GMbFP3Z/zMfY/P+tXF56z83vO7bO5rbLt3ik7dHbPiIVleqUj6V9jak7QTsRRLmtKB+fuHxs/XFtRqSqUSrefP7obSsxtryoPvb9xrRTOl6EVj6PzfVL54wpH42v/BZR/j9Cudk6MXrK6LFx+ykixoVJAM8w9NOIeCc7qSeSD8MbGTzsA6g2ULwGBqjVBqhNGUS12ubyxrEpg9TKZelVmzo4dr8mVCuS1cCUGqqJmlJZTdRqojagzfXGHR8YqKEaY9qp1cSUKTUGahrzmlITg1NqDNTGHVPx5y4TnDMwes7Y+jWJAcHUgdrodq02SXlqvyaYWqsxUEtljXZK+xLUUiKspaTXaKMoKxJoo85AY5vGecX+gESNgKijqBfPqkWg+nDaTq966Xh9pFS3vvnckZFx7RR1o16H4SFiZATqjbJiO4aHij9HRqA+MrZ+qkt9ZMy59eFNxEi9eNXro9sjm8aW19N2fdPwaHl9tH7xPvWh4aLeaHlQH6kzsqm8H+ncIEaCkaGR0XpR33x8ZKjxvkW9GCm2R4brjEQwEmzx51B9a+VBncmPjQT8P367VxW/tr7sMrMsnHzMLAsnHzPLwsnHzLJw8jGzLJx8zCwLJx8zy8LJx8yycPIxsyycfMwsCycfM8vCycfMsnDyMbMsnHzMLAsnHzPLwsnHzLJw8jGzLJx8zCwLJx8zy8LJx8yycPIxsyycfMwsC41+11AXk7QcmAk8kzsW2y574Z9ZN3ptROy/s430xvd2wUZgY0S8KXcg1jxJ9/hn1n0kra+iHV92mVkWTj5mlkWvXHZ9I3cAtkP8c+tO36+ikZ4YcDaz7uPLLjPLouOTjwpLJQ1J2iDpQxPUmSHpaUkb0+s+SQPp2FGSnpX0iqTfSzq8/Z+i/0h6p6RVkh6W9NkJju8i6fp0/C5Js0vHPpfKV0k6tZ1x9ztJF6TftSFJt05w/Nr0O7ZB0guS3lM6dlvp3L/e5ptFREe/gC8AawEBZwPrJqgzAzg9bQ9S3DtySdp/CLghbX8TWJ77M/X6CxgAHgEOSj+P+4HDxtU5F7gibc8Drk/bh6X6uwAHpnYGcn+mfngBU4FNwEnArsAG4D3j6uxf2r4YWJu235Pq7wacmNqZurX36/ieD8VfzAVR+BdgUNIbyxUi4uWIuLmxC7wI7JH2DwS+lLa/ALyhDTH3u+OAhyPi0YgYAhYA7x1X573At9L2jcDbJCmVL4iIVyJiNfBwas9a7yzghYi4PSLWA7cD55UrRMQTpd1Zpe3zgNsj4qWIuAN4IbU3qW5IPnsCD5b21wFHTVRR0k8pekl7A19uFAMnpO03AzVJe7YmVEv2Bx4v7a9JZRPWiYhhir+sezZ5rrXGIYy94/wxYL/xldLl8ibgfwLzU/F+wKOlas+m9ibVDclne5wO3EXRVT8wlW0ATpJ0H0V3MoCRPOGZdb+IOCMipgJXpNcO6cjkkzLrBkkbgOcoxgEaZgIPTHLqN4BVwN+zuZv/FHBeRBwDXA6MRMTzrYnckieAA0r7f5DKJqwjaQpFF/7ZJs+11lhF8bxdw2zgya3U/wQwJ20/STHG17Bnam9SHZl8UmadHhHTKcYL5qVZr7OBoYi4v1xf0kxJ/0jxF/jTFD2gh9Lhn7P52vObwK/a8Rn63N3AwZIOlDRIMW5307g6NwEfSdsfBH4RxcjlTRQ/710kHQgcjH9m7fJtYJakEyXtSnGl8LVyBUlvL+1eCLyctr9GcYWxm6QTKX4Xr9rqu+UeYW9iBF7AMorR843Ah0vHNqQ/j6a4nNpIcZm1FjgnHTsr/Q8aovgXdN/cn6kfXsC7gF9TXAJfkMouAv4kbU8DvksxoPwr4KDSuRek81YBp+X+LP30opiUGUq/bz9LZbcDF6ft+0q/Z881fp7p2M/SeUPA57f1Xr7D2cyy6MjLLjPrfU4+ZpaFk4+ZZeHkY2ZZOPmYWRZOPmaWhZOPmWXh5GNmWfx/Rq+3uFS5EEgAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 360x360 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    }
  ]
}