{
"cells": [
{
"cell_type": "markdown",
"id": "32647951-7dcf-4553-b42a-5d8508d95cf9",
"metadata": {},
"source": [
"# Data Explorer\n",
"\n",
"
"
]
},
{
"cell_type": "markdown",
"id": "f07451ff-cf56-4763-b4ac-bd099ba1f145",
"metadata": {},
"source": [
"\n",
"# Table of Contents\n",
"\n",
"- [Introduction](#intro)\n",
"- [Setup](#setup)\n",
"- [Get your bouts and trials here](#one)\n",
"- [CAP analysis](#two)\n"
]
},
{
"cell_type": "markdown",
"id": "38b445ff-7cf7-442e-bb25-83afe6d5a3fb",
"metadata": {},
"source": [
"\n",
"# Giant Fiber System Activity\n",
"\n",
"In this lab you measured the activity of the giant fiber system using extracellular differential electrodes. The activity that you will observe in your recording is referred to as a **complex action potential** (CAP). You will apply your knowledge from the previous two labs to interpret this signal. One of the predominant analysis frameworks you will use is comparing *trials* across *bouts*."
]
},
{
"cell_type": "markdown",
"id": "488d476b-9387-43d7-b4bf-0528d9f4950a",
"metadata": {},
"source": [
"\n",
"# Setup\n",
"\n",
"[toc](#toc)"
]
},
{
"cell_type": "markdown",
"id": "ac5e6fe8-5ee8-4a44-a353-785ed3d50412",
"metadata": {},
"source": [
"## Import and define functions"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "b907ad44-4d78-45f0-8775-8cdaf6584ceb",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\" }\n",
"\n",
"#@markdown Run this code cell to import packages and define functions \n",
"import numpy as np\n",
"import pandas as pd\n",
"import plotly.express as px\n",
"import plotly.graph_objects as go\n",
"from plotly.subplots import make_subplots\n",
"from scipy import ndimage\n",
"from scipy.signal import hilbert,medfilt,resample, find_peaks, unit_impulse\n",
"import seaborn as sns\n",
"from datetime import datetime,timezone,timedelta\n",
"pal = sns.color_palette(n_colors=15)\n",
"pal = pal.as_hex()\n",
"import matplotlib.pyplot as plt\n",
"import random\n",
"\n",
"from pathlib import Path\n",
"\n",
"from matplotlib.ticker import (AutoMinorLocator, MultipleLocator)\n",
"from ipywidgets import widgets, interact, interactive\n",
"%config InlineBackend.figure_format = 'retina'\n",
"plt.style.use(\"https://raw.githubusercontent.com/NeuromatchAcademy/course-content/master/nma.mplstyle\")\n",
"\n",
"print('Task completed at ' + str(datetime.now(timezone(-timedelta(hours=5)))))"
]
},
{
"cell_type": "markdown",
"id": "85f5b456-2a84-416a-8f1d-f480c7a10ad5",
"metadata": {},
"source": [
"## Mount Google Drive"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d84a46ac-5cf3-4765-94fd-07e19f0881d2",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\" }\n",
"\n",
"#@markdown Run this cell to mount your Google Drive.\n",
"\n",
"from google.colab import drive\n",
"drive.mount('/content/drive')\n",
"\n",
"print('Task completed at ' + str(datetime.now(timezone(-timedelta(hours=5)))))"
]
},
{
"cell_type": "markdown",
"id": "dfd155f9-cf42-4d7d-98a7-3f2432cf2e5e",
"metadata": {},
"source": [
"## Import data \n",
"\n",
"Import data digitized with *Nidaq USB6211* and recorded using *Bonsai-rx* as a *.bin* file\n",
"\n",
"> If you would like to explore the analysis for this lab, but do not have data, you can download examples for the following experiments using the linked shared file:\n",
" - [CAP measured across different inter-electrode distances (Part II. Experiment 1)](https://drive.google.com/file/d/1SyTJtg0GgHtXdJYQV51QIuu2zCRcmbI3/view?usp=sharing). If you are using this example file, the sample rate was 30000 with three channels (channel 0 was the nerve signal and channel 2 was the stimulus monitor; channel 1 was unused). The inter-electrode distances were: 1, 2, 3 and 4. Bouts approximately: [[0,15],[15,35],[35,60],[70,100]]. The stimulus anode was 3cm from the closest measurement electrode."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "6fc27eac-7157-4510-9239-466ca324c58d",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\" }\n",
"\n",
"#@markdown Specify the file path \n",
"#@markdown to your recorded data in the colab runtime (find the filepath in the colab file manager):\n",
"\n",
"filepath = \"full filepath goes here\" #@param \n",
"# filepath = '/Volumes/Untitled/BIOL247/data/earthworm-giant-fiber-cap/diff_cv_0.bin' #@param \n",
"# filepath = '/Users/kperks/OneDrive - wesleyan.edu/Teaching/Neurophysiology_FA22/data/earthworm-cap/stim_2chan.bin'\n",
"\n",
"#@markdown Specify the sampling rate and number of channels recorded.\n",
"\n",
"sampling_rate = 30000 #@param\n",
"number_channels = 3 #@param\n",
"\n",
"# downsample = False #@param\n",
"# newfs = 10000 #@param\n",
"\n",
"#@markdown After you have filled out all form fields, \n",
"#@markdown run this code cell to load the data. \n",
"\n",
"filepath = Path(filepath)\n",
"\n",
"# No need to edit below this line\n",
"#################################\n",
"data = np.fromfile(Path(filepath), dtype = np.float64)\n",
"data = data.reshape(-1,number_channels)\n",
"data_dur = np.shape(data)[0]/sampling_rate\n",
"print('duration of recording was %0.2f seconds' %data_dur)\n",
"\n",
"fs = sampling_rate\n",
"# if downsample:\n",
"# # newfs = 10000 #downsample emg data\n",
"# chunksize = int(sampling_rate/newfs)\n",
"# data = data[0::chunksize,:]\n",
"# fs = int(np.shape(data)[0]/data_dur)\n",
"\n",
"time = np.linspace(0,data_dur,np.shape(data)[0])\n",
"\n",
"print('Data upload completed at ' + str(datetime.now(timezone(-timedelta(hours=5)))))"
]
},
{
"cell_type": "markdown",
"id": "caaf4c27-c930-4542-9f51-aba353b37a12",
"metadata": {},
"source": [
"## Plot raw data"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "cb4451d0-cff7-484e-bf5d-3e65c425b147",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\"}\n",
"\n",
"#@markdown Run this code cell to plot the imported data.
\n",
"#@markdown Use the range slider to scroll through the data in time.\n",
"#@markdown Use the channel slider to choose which channel to plot\n",
"#@markdown Be patient with the range refresh... the more data you are plotting the slower it will be. \n",
"\n",
"slider_xrange = widgets.FloatRangeSlider(\n",
" min=0,\n",
" max=data_dur,\n",
" value=(0,1),\n",
" step= 1,\n",
" readout=True,\n",
" continuous_update=False,\n",
" description='Time Range (s)')\n",
"slider_xrange.layout.width = '600px'\n",
"\n",
"slider_chan = widgets.IntSlider(\n",
" min=0,\n",
" max=number_channels-1,\n",
" value=0,\n",
" step= 1,\n",
" continuous_update=False,\n",
" description='channel')\n",
"slider_chan.layout.width = '300px'\n",
"\n",
"# a function that will modify the xaxis range\n",
"def update_plot(x,chan):\n",
" fig, ax = plt.subplots(figsize=(10,5),num=1); #specify figure number so that it does not keep creating new ones\n",
" starti = int(x[0]*fs)\n",
" stopi = int(x[1]*fs)\n",
" ax.plot(time[starti:stopi], data[starti:stopi,chan])\n",
"\n",
"w = interact(update_plot, x=slider_xrange, chan=slider_chan);"
]
},
{
"cell_type": "markdown",
"id": "69eef312-a1ec-4dc6-ba47-a45be1f91e58",
"metadata": {},
"source": [
"For a more extensive ***RAW*** Data Explorer than the one provided in the above figure, use the [DataExplorer.py](https://raw.githubusercontent.com/neurologic/Neurophysiology-Lab/main/howto/Data-Explorer.py) application found in the [howto section](https://neurologic.github.io/Neurophysiology-Lab/howto/Dash-Data-Explorer.html) of the course website."
]
},
{
"cell_type": "markdown",
"id": "be92808a-e067-4662-a7ae-57707d90765a",
"metadata": {},
"source": [
"\n",
"# Part I. Get your bouts and trials here\n",
"\n",
"The presentation of each stimulus marks the start of each ***trial*** in your experiment. Changes in stimulus parameters mark the start/end of different ***bouts*** in your experiment. \n",
"\n",
"Therefore, our first task in processing and analyzing data from the experiment is to figure out the trial times. You will do this using the same algorithm that you used to detect EOD events in the last lab. In this analysis, you will be processing the stimulus monitor signal directly."
]
},
{
"cell_type": "markdown",
"id": "5e30e1b5-5812-484f-96fa-c9a722877b89",
"metadata": {},
"source": [
"## Define trial times"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d0256f0d-3f55-4fcb-9325-e71513f61157",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\"}\n",
"\n",
"#@markdown Run this cell to create an interactive plot with a slider to scroll \n",
"#@markdown through the signal\n",
"#@markdown and set an appropriate event detection threshold \n",
"#@markdown (you can do so based on level crossing or peaks). \n",
"\n",
"slider_xrange = widgets.FloatRangeSlider(\n",
" min=0,\n",
" max=data_dur,\n",
" value=(0,1),\n",
" step= 0.05,\n",
" readout=True,\n",
" continuous_update=False,\n",
" description='Time Range (s)',\n",
" style = {'description_width': '200px'})\n",
"slider_xrange.layout.width = '600px'\n",
"\n",
"# slider_yrange = widgets.FloatRangeSlider(\n",
"# min=np.min(stim)-0.5,\n",
"# max=np.max(stim)+0.5,\n",
"# value=[np.min(stim),np.max(stim)],\n",
"# step=0.05,\n",
"# continuous_update=False,\n",
"# readout=True,\n",
"# description='yrange',\n",
"# style = {'description_width': '200px'})\n",
"# slider_yrange.layout.width = '600px'\n",
"\n",
"select_channel = widgets.Select(\n",
" options=np.arange(np.shape(data)[1]), # start with a single trial on a single bout... it will update when runs ; old: np.arange(len(trial_times)),\n",
" value=0,\n",
" #rows=10,\n",
" description='Channel used to detect events',\n",
" style = {'description_width': '200px'},\n",
" disabled=False\n",
")\n",
"\n",
"slider_threshold = widgets.FloatSlider(\n",
" min=-2,\n",
" max=2,\n",
" value=0.2,\n",
" step=0.001,\n",
" readout_format='.3f',\n",
" continuous_update=False,\n",
" readout=True,\n",
" description='event detection threshold',\n",
" style = {'description_width': '200px'})\n",
"slider_threshold.layout.width = '600px'\n",
"\n",
"detect_type_radio = widgets.RadioButtons(\n",
" options=['peak', 'level crossing'],\n",
" value='peak', # Defaults to 'level crossing'\n",
" layout={'width': 'max-content'}, # If the items' names are long\n",
" description='Type of event detection',\n",
" style = {'description_width': '200px'},\n",
" disabled=False\n",
")\n",
"\n",
"radio_polarity = widgets.RadioButtons(\n",
" options=[1, -1],\n",
" value=-1,\n",
" description='peaks polarity',\n",
" disabled=False,\n",
" style = {'description_width': '200px'}\n",
")\n",
"\n",
"iei_text = widgets.Text(\n",
" value='0.005',\n",
" placeholder='0.005',\n",
" description='min IEI (seconds)',\n",
" style = {'description_width': '200px'},\n",
" disabled=False\n",
")\n",
"\n",
"def update_plot(chan_,xrange,thresh_,detect_type,polarity,iei):\n",
" fig, ax = plt.subplots(figsize=(10,5),num=1); #specify figure number so that it does not keep creating new ones\n",
" \n",
" signal = data[:,chan_]\n",
" signal = signal-np.median(signal)\n",
" \n",
" iei = float(iei)\n",
" \n",
" if iei>0.001:\n",
" d = iei*fs #minimum time allowed between distinct events\n",
" \n",
" if detect_type == 'peak':\n",
" r = find_peaks(signal*polarity,height=thresh_,distance=d)\n",
" trial_times = r[0]/fs\n",
" # print(r[1])\n",
" ax.scatter(trial_times,[r[1]['peak_heights']*polarity],marker='^',s=300,color='purple',zorder=3)\n",
" \n",
" if detect_type == 'level crossing':\n",
" # get the changes in bool value for a bool of signal greater than threshold\n",
" # if polarity == 1:\n",
" threshold_crossings = np.diff(signal*polarity > thresh_, prepend=False)\n",
" # get indices where threshold crossings are true\n",
" tcross = np.argwhere(threshold_crossings)[:,0]\n",
" # get a mask for only positive level crossings\n",
" mask_ = [signal[t]-signal[t-1] > 0 for t in tcross]\n",
" # if polarity == -1:\n",
" # threshold_crossings = np.diff(signal*polarity < thresh_*polarity, prepend=False)\n",
" # # get indices where threshold crossings are true\n",
" # tcross = np.argwhere(threshold_crossings)[:,0]\n",
" # # get a mask for only positive level crossings\n",
" # mask_ = [signal[t]-signal[t-1] > 0 for t in tcross]\n",
" \n",
" # trial times are positive level crossings\n",
" trial_times = tcross[mask_]/fs\n",
" ax.scatter(trial_times,[thresh_*polarity]*len(trial_times),marker='^',s=300,color='purple',zorder=3)\n",
"\n",
" starti = int(xrange[0]*fs)+1\n",
" stopi = int(xrange[1]*fs)-1\n",
" ax.plot(time[starti:stopi], signal[starti:stopi], color='black')\n",
" \n",
" # ax.plot(tmp,color='black')\n",
" ax.hlines(thresh_*polarity, time[starti],time[stopi],linestyle='--',color='green')\n",
" \n",
" # ax.set_ylim(yrange[0],yrange[1])\n",
" ax.set_xlim(xrange[0],xrange[1])\n",
" \n",
"\n",
" ax.xaxis.set_minor_locator(AutoMinorLocator(5))\n",
"\n",
" \n",
" return trial_times\n",
"\n",
"w_trials_ = interactive(update_plot, chan_=select_channel, \n",
" xrange=slider_xrange, \n",
" thresh_=slider_threshold, detect_type = detect_type_radio, \n",
" polarity=radio_polarity, iei = iei_text);\n",
"display(w_trials_)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "0ebea80c-6fc4-4ad0-900a-8c4bf314f6eb",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\"}\n",
"\n",
"#@markdown Run this cell to finalize the list of trial times \n",
"#@markdown after settling on a channel and threshold in the interactive plot.
\n",
"#@markdown This stores the trial times in an array called 'trial_times'.\n",
"trial_times = w_trials_.result\n",
"\n"
]
},
{
"cell_type": "markdown",
"id": "2c61993b-bad0-4f65-853d-efa2eee9f37f",
"metadata": {},
"source": [
"## Define Bouts"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "5c081e0b-7544-43e8-9b7b-c8ed8ff3eece",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode: \"form\"}\n",
"\n",
"#@markdown Use the interactive plot of the stimulus monitor channel in the last section to determine \n",
"#@markdown the start and stop time for each of the bouts in your experiment. \n",
"#@markdown Specify the list of bout ranges as follows: [[start of bout 0, end of bout 0],[start 1, end 1],...]]
\n",
"\n",
"bouts_list = [[0,30],[35,45]] #@param\n",
"# bouts_list = [[2,10],[10,20],[20,30],[30,45],[45,55],[55,70],[70,85],[85,100],[100,120]]\n",
"\n",
"#@markdown Then run this code cell to programatically define the list of bouts as 'bouts_list'."
]
},
{
"cell_type": "markdown",
"id": "d9823a59-8623-454e-b18f-5380d150d174",
"metadata": {},
"source": [
"\n",
"# Part II. CAP analysis"
]
},
{
"cell_type": "markdown",
"id": "cdda2714-862d-478b-a6df-be178d761892",
"metadata": {},
"source": [
"## Visualize Trials"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "60f56f25-65ac-43bf-a986-56d48dab703a",
"metadata": {
"tags": [
"hide-input"
]
},
"outputs": [],
"source": [
"#@title {display-mode:\"form\"}\n",
"\n",
"#@markdown Run this code cell to create an interactive plot to \n",
"#@markdown examine the CAP on individual trials for each bout.\n",
"#@markdown You can overlay multple channels by selecting more than one.\n",
"#@markdown You can overlay multiple trials by selecting more than one. \n",
"#@markdown (To select more than one item from an option menu, press the control/command key \n",
"#@markdown while mouse clicking or shift while using up/down arrows)\n",
"\n",
"slider_xrange = widgets.FloatRangeSlider(\n",
" min=-0.01,\n",
" max=0.05,\n",
" value=(-0.001,0.03),\n",
" step=0.0005,\n",
" continuous_update=False,\n",
" readout=True,\n",
" readout_format='.4f',\n",
" description='xrange (s)'\n",
")\n",
"slider_xrange.layout.width = '600px'\n",
"\n",
"slider_yrange = widgets.FloatRangeSlider(\n",
" min=-1,\n",
" max=1, # normal range for earthworm experiments\n",
" value=(-0.5,0.5),\n",
" step=0.01,\n",
" continuous_update=False,\n",
" readout=True,\n",
" description='yrange'\n",
")\n",
"slider_yrange.layout.width = '600px'\n",
"\n",
"# trials in bout 0 to start...\n",
"trials_t = trial_times[(trial_times>bouts_list[0][0]) & (trial_timesbouts_list[bout_][0]) & (trial_times0) & ((int(fs*t_)+win_1)) Note that there must be a list of trials for every bout.\n",
"trials_all = [[0,1,2],[0,1,2],[0,1,2],[0,1,2]] #@param\n",
"#@markdown Specify the time range around each trial time that you want to include in the plot. \n",
"xrange = [0,0.006] #@param\n",
"\n",
"#@markdown Now run this code cell to create a plot that shows the trial-averaged response on each bout\n",
"\n",
"win_0 = int(xrange[0]*fs)\n",
"win_1 = int(xrange[1]*fs)\n",
"xtime = np.linspace(xrange[0],xrange[1],(win_1 - win_0))\n",
"\n",
"fig, ax = plt.subplots(figsize=(15,8),nrows=len(channel_list),ncols=1)\n",
"\n",
"for j,chan_ in enumerate(channel_list):\n",
" \n",
" this_chan = data[:,chan_]\n",
" \n",
" for i,bout_ in enumerate(bouts_all):\n",
" data_avg = []\n",
" trial_list = trials_all[i]\n",
" trials_t = trial_times[(trial_times>bouts_list[bout_][0]) & (trial_times0) & ((int(fs*t_)+win_1))1:\n",
" ax[j].plot(xtime,np.mean(data_avg,1),linewidth=2,label = f'bout {bout_}')\n",
" if len(channel_list)==1:\n",
" ax.plot(xtime,np.mean(data_avg,1),linewidth=2,label = f'bout {bout_}')\n",
"\n",
" \n",
" if len(channel_list)>1:\n",
"\n",
" # ax.set_ylim(yrange[0],yrange[1]);\n",
" ax[j].set_xlabel('seconds')\n",
" ax[j].legend()\n",
" # ax.vlines(0,yrange[0],yrange[1],color='green')\n",
" ax[j].xaxis.set_major_locator(MultipleLocator(0.002))\n",
" ax[j].xaxis.set_minor_locator(AutoMinorLocator(10))\n",
" ax[j].grid(which='major', color='gray', linestyle='-')\n",
" ax[j].grid(which='minor', color='gray', linestyle=':')\n",
"\n",
"if len(channel_list)==1:\n",
" # ax.set_ylim(yrange[0],yrange[1]);\n",
" ax.set_xlabel('seconds')\n",
" ax.legend()\n",
" # ax.vlines(0,yrange[0],yrange[1],color='green')\n",
" ax.xaxis.set_major_locator(MultipleLocator(0.002))\n",
" ax.xaxis.set_minor_locator(AutoMinorLocator(10))\n",
" ax.grid(which='major', color='gray', linestyle='-')\n",
" ax.grid(which='minor', color='gray', linestyle=':')"
]
},
{
"cell_type": "markdown",
"id": "d451b557-6777-4492-955e-a41130ad39b9",
"metadata": {},
"source": [
"
\n",
"Written by Dr. Krista Perks for courses taught at Wesleyan University."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "5c707588-d5f6-4142-8289-da7ba6a97dbc",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"id": "57a8a68c-055e-4af0-ba4d-67d0d67148f1",
"metadata": {},
"source": [
""
]
},
{
"cell_type": "markdown",
"id": "a350712e-e145-4475-9588-e3f19469b5ce",
"metadata": {},
"source": [
""
]
},
{
"cell_type": "markdown",
"id": "86fd2b4a-b890-465f-bf33-c2057738789e",
"metadata": {},
"source": [
""
]
},
{
"cell_type": "markdown",
"id": "b05de192-2baf-4d46-8848-6db1e3b50bff",
"metadata": {},
"source": [
""
]
},
{
"cell_type": "markdown",
"id": "e6519277-fc86-44fd-a75b-f74cad0efcf1",
"metadata": {},
"source": [
""
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.13"
},
"widgets": {
"application/vnd.jupyter.widget-state+json": {
"state": {},
"version_major": 2,
"version_minor": 0
}
}
},
"nbformat": 4,
"nbformat_minor": 5
}