Alexis - Data Explorer
Contents
Motor Unit Population Stimulation¶
Comparing electrically stimulated motor unit activity to voluntary motor unit activity. And comparing electrically stimulated motor unit activity across conditions (for example temperature).
Setup¶
Import and define functions¶
#@title {display-mode: "form" }
#@markdown Run this code cell to import packages and define functions
import numpy as np
import pandas as pd
import plotly.express as px
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from scipy import ndimage
from scipy.signal import hilbert,medfilt,resample, find_peaks, unit_impulse
from scipy.optimize import curve_fit
import seaborn as sns
from datetime import datetime,timezone,timedelta
pal = sns.color_palette(n_colors=15)
pal = pal.as_hex()
import matplotlib.pyplot as plt
import random
from numpy import NaN
from pathlib import Path
from matplotlib.ticker import (AutoMinorLocator, MultipleLocator)
from ipywidgets import widgets, interact, interactive,interactive_output
%config InlineBackend.figure_format = 'retina'
plt.style.use("https://raw.githubusercontent.com/NeuromatchAcademy/course-content/master/nma.mplstyle")
def hyper_fit(t,r,c):
return r + ((r*c)/t)
print('Task completed at ' + str(datetime.now(timezone(-timedelta(hours=5)))))
Mount Google Drive¶
#@title {display-mode: "form" }
#@markdown Run this cell to mount your Google Drive.
from google.colab import drive
drive.mount('/content/drive')
print('Task completed at ' + str(datetime.now(timezone(-timedelta(hours=5)))))
Data Processing¶
In this section, you will load the raw data aquired and saved using Bonsai-rx as a .bin file. A single trial consists of a stimulus pulse (from the TENS). Each bout consists of a set of trials at a given frequency.
Import data¶
Import data digitized with Nidaq USB6211 and recorded using Bonsai-rx as a .bin file
#@title {display-mode: "form" }
#@markdown Specify the file path
#@markdown to your recorded data in the colab runtime (find the filepath in the colab file manager):
filepath = "full filepath goes here" #@param
# filepath = '/Users/kperks/Downloads/Alexis Capstone/hot_temp_0.bin' #@param
#@markdown Specify the sampling rate and number of channels recorded.
sampling_rate = None #@param
number_channels = None #@param
# sampling_rate = 40000 #@param
# number_channels = 1 #@param
# downsample = False #@param
# newfs = 10000 #@param
#@markdown After you have filled out all form fields,
#@markdown run this code cell to load the data.
filepath = Path(filepath)
# No need to edit below this line
#################################
data = np.fromfile(Path(filepath), dtype = np.float64)
data = data.reshape(-1,number_channels)
data_dur = np.shape(data)[0]/sampling_rate
print('duration of recording was %0.2f seconds' %data_dur)
fs = sampling_rate
# if downsample:
# # newfs = 10000 #downsample emg data
# chunksize = int(sampling_rate/newfs)
# data = data[0::chunksize,:]
# fs = int(np.shape(data)[0]/data_dur)
time = np.linspace(0,data_dur,np.shape(data)[0])
print('Data upload completed at ' + str(datetime.now(timezone(-timedelta(hours=5)))))
For a more extensive RAW Data Explorer than the one provided in the above figure, use the DataExplorer.py application found in the howto section of the course website.
Define bout and stimulus times¶
The time between stimulus onset and action potential, and the time between two stimulus pulses are critical parameters of the data on each trial.
Our first task in processing and analyzing data from this experiment is to figure out the stimulus onset times. You can then segment the data in to separate bouts if the raw recording was not one continuous successful protocol.
Define stimulus times¶
#@title {display-mode: "form"}
#@markdown Run this cell to create an interactive plot with a slider to scroll
#@markdown through the signal
#@markdown and set an appropriate event detection threshold
#@markdown (you can do so based on level crossing or peaks).
slider_xrange = widgets.FloatRangeSlider(
min=0,
max=data_dur,
value=(0,1),
step= 0.05,
readout=True,
continuous_update=False,
description='Time Range (s)',
style = {'description_width': '200px'})
slider_xrange.layout.width = '600px'
# slider_yrange = widgets.FloatRangeSlider(
# min=np.min(stim)-0.5,
# max=np.max(stim)+0.5,
# value=[np.min(stim),np.max(stim)],
# step=0.05,
# continuous_update=False,
# readout=True,
# description='yrange',
# style = {'description_width': '200px'})
# slider_yrange.layout.width = '600px'
select_channel = widgets.Select(
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)),
value=0,
#rows=10,
description='Channel used to detect events',
style = {'description_width': '200px'},
disabled=False
)
slider_threshold = widgets.FloatSlider(
min=-2,
max=2,
value=0.2,
step=0.001,
readout_format='.3f',
continuous_update=False,
readout=True,
description='event detection threshold',
style = {'description_width': '200px'})
slider_threshold.layout.width = '600px'
detect_type_radio = widgets.RadioButtons(
options=['peak', 'level crossing'],
value='peak', # Defaults to 'level crossing'
layout={'width': 'max-content'}, # If the items' names are long
description='Type of event detection',
style = {'description_width': '200px'},
disabled=False
)
radio_polarity = widgets.RadioButtons(
options=[1, -1],
value=-1,
description='peaks polarity',
disabled=False,
style = {'description_width': '200px'}
)
iei_text = widgets.Text(
value='0.005',
placeholder='0.005',
description='min IEI (seconds)',
style = {'description_width': '200px'},
disabled=False
)
def update_plot(chan_,xrange,thresh_,detect_type,polarity,iei):
fig, ax = plt.subplots(figsize=(10,5),num=1); #specify figure number so that it does not keep creating new ones
signal = data[:,chan_]
signal = signal-np.median(signal)
iei = float(iei)
if iei>0.001:
d = iei*fs #minimum time allowed between distinct events
if detect_type == 'peak':
r = find_peaks(signal*polarity,height=thresh_,distance=d)
trial_times = r[0]/fs
# print(r[1])
ax.scatter(trial_times,[r[1]['peak_heights']*polarity],marker='^',s=300,color='purple',zorder=3)
if detect_type == 'level crossing':
# get the changes in bool value for a bool of signal greater than threshold
# if polarity == 1:
threshold_crossings = np.diff(signal*polarity > thresh_, prepend=False)
# get indices where threshold crossings are true
tcross = np.argwhere(threshold_crossings)[:,0]
# get a mask for only positive level crossings
mask_ = [signal[t]-signal[t-1] > 0 for t in tcross]
# if polarity == -1:
# threshold_crossings = np.diff(signal*polarity < thresh_*polarity, prepend=False)
# # get indices where threshold crossings are true
# tcross = np.argwhere(threshold_crossings)[:,0]
# # get a mask for only positive level crossings
# mask_ = [signal[t]-signal[t-1] > 0 for t in tcross]
# trial times are positive level crossings
trial_times = tcross[mask_]/fs
ax.scatter(trial_times,[thresh_*polarity]*len(trial_times),marker='^',s=300,color='purple',zorder=3)
starti = int(xrange[0]*fs)+1
stopi = int(xrange[1]*fs)-1
ax.plot(time[starti:stopi], signal[starti:stopi], color='black')
# ax.plot(tmp,color='black')
ax.hlines(thresh_*polarity, time[starti],time[stopi],linestyle='--',color='green')
# ax.set_ylim(yrange[0],yrange[1])
ax.set_xlim(xrange[0],xrange[1])
ax.xaxis.set_minor_locator(AutoMinorLocator(5))
return trial_times
w_trials_ = interactive(update_plot, chan_=select_channel,
xrange=slider_xrange,
thresh_=slider_threshold, detect_type = detect_type_radio,
polarity=radio_polarity, iei = iei_text);
display(w_trials_)
#@title {display-mode: "form"}
#@markdown Run this cell to finalize the list of event times
#@markdown after settling on a channel and threshold in the interactive plot. <br>
#@markdown This stores the event times in an array called 'event_times'. <br>
trial_times = w_trials_.result
Define Bouts¶
#@title {display-mode: "form"}
#@markdown For this experiment, the entire file should be one long bout,
#@markdown but if there were regions that something got messed up or that you want to exclude, you can specify bouts with good data.
#@markdown Specify the list of bout ranges as follows: [[start of bout 0, end of bout 0],[start 1, end 1],...]] <br>
bouts_list = [[None,None]] #@param
# bouts_list = [[5,10],[10,15],[15,20],[20,27],[27,35]] #@param #control condition
# bouts_list = [[12,17],[17,23],[23,30],[30,35],[35,42],[42,47]] # cold condition
# bouts_list = [[20,25],[25,30],[30,35],[35,40]] # hot condition
#@markdown Then run this code cell to programatically define the list of bouts as 'bouts_list'.
Analyze Data¶
Measure the raw data¶
Obtain necessary information from the raw signal time-locked to each trial start time (stimulus time).
#@title {display-mode:"form"}
#@markdown Run this code cell to create an interactive plot to
#@markdown examine the raw signal time-locked to individual events (event_times).
#@markdown You can overlay multple channels by selecting more than one.
#@markdown You can overlay multiple event times by selecting more than one.
#@markdown (To select more than one item from an option menu, press the control/command key
#@markdown while mouse clicking or shift while using up/down arrows)
slider_xrange = widgets.FloatRangeSlider(
min=-0.01,
max=0.05,
value=(-0.001,0.03),
step=0.0005,
continuous_update=False,
readout=True,
readout_format='.3f',
description='xrange (s)'
)
slider_xrange.layout.width = '600px'
slider_yrange = widgets.FloatRangeSlider(
min=-5,
max=5, # normal range for earthworm experiments
value=(-0.5,1),
step=0.05,
continuous_update=False,
readout=True,
description='yrange'
)
slider_yrange.layout.width = '600px'
ui_range = widgets.VBox([slider_xrange, slider_yrange])
# trials in bout 0 to start...
trials_t = trial_times[(trial_times>bouts_list[0][0]) & (trial_times<bouts_list[0][1])]
select_channels = widgets.SelectMultiple(
options=np.arange(np.shape(data)[1]), # start with a single trial on a single bout... it will update when runs ,
value=[0],
#rows=10,
description='Channels',
disabled=False
)
select_bouts = widgets.Select(
options=np.arange(len(bouts_list)), # start with a single trial on a single bout... it will update when runs ; old: np.arange(len(trial_times)),
value=0,
#rows=10,
description='Bouts',
disabled=False
)
select_trials = widgets.SelectMultiple(
options=np.arange(len(trials_t)), # start with a single trial on a single bout... it will update when runs ,
value=[0],
#rows=10,
description='Events',
disabled=False
)
ui_trials = widgets.HBox([select_channels, select_trials, select_bouts])
slider_threshold = widgets.FloatSlider(
min=0,
max=1,
value=0.25,
step=0.001,
readout_format='.3f',
continuous_update=False,
readout=True,
description='peak detection threshold',
style = {'description_width': '250px'})
slider_threshold.layout.width = '600px'
detect_chan_radio = widgets.RadioButtons(
options=['none'] + [v.astype(str) for v in np.arange(np.shape(data)[1])],
value='none', # Defaults to 'none'
layout={'width': 'max-content'}, # If the items' names are long
description='detect delay to peaks on channel: ',
style = {'description_width': '250px'},
disabled=False
)
slider_inwindow = widgets.FloatRangeSlider(
min= -0.01,
max= 0.1,
value=(0.004,0.015),
step=0.0001,
continuous_update=False,
readout=True,
description='window to detect peaks within',
style = {'description_width': '250px'},
readout_format='.4f'
)
slider_inwindow.layout.width = '600px'
radio_polarity = widgets.RadioButtons(
options=[1, -1],
value=1,
description='peaks polarity',
disabled=False,
style = {'description_width': '200px'}
)
iei_text = widgets.Text(
value='0.001',
placeholder='0.001',
description='minimum interval (s)',
style = {'description_width': '200px'},
disabled=False
)
ui_peaks = widgets.VBox([detect_chan_radio,slider_inwindow, radio_polarity, slider_threshold, iei_text])
trial_abs_readout = widgets.Label(
value=f'time of this event is (sec): {NaN}'
)
trial_abs_readout.layout.width = '600px'
trial_readout = widgets.Label(
value=f'time since last event is: {NaN}'
)
trial_readout.layout.width = '600px'
lagging_time_readout = widgets.Label(
value=f'lagging peak times are: {NaN}'
)
lagging_time_readout.layout.width = '600px'
lagging_amp_readout = widgets.Label(
value=f'lagging peak amplitudes are: {NaN}'
)
lagging_amp_readout.layout.width = '600px'
def update_plot(chan_list,trial_list,bout_,yrange,xrange,lagging_chan_,detect_inwin,thresh_,iei,polarity):
fig, ax = plt.subplots(figsize=(8,4))# ,ncols=1, nrows=1); #specify figure number so that it does not keep creating new ones
iei = float(iei)
win_0 = int(xrange[0]*fs)
win_1 = int(xrange[1]*fs)
xtime = np.linspace(xrange[0],xrange[1],(win_1 - win_0))
trials_t = trial_times[(trial_times>bouts_list[bout_][0]) & (trial_times<bouts_list[bout_][1])]
trials_init_ = np.arange(len(trials_t))
select_trials.options = trials_init_
trial_list = [t_try for t_try in trial_list if t_try in trials_init_]
select_trials.value = trial_list
lagging_time_readout.value=f'lagging peak times are: {NaN}'
lagging_amp_readout.value=f'lagging peak amplitudes are: {NaN}'
trial_abs_readout.value=f'time of this event is: {NaN}'
trial_readout.value=f'time since last event is: {NaN}'
channel_colors = ['purple','green','blue','orange']
for chan_ in chan_list:
this_chan = data[:,chan_]
for trial_ in trial_list:
if trial_ in trials_init_:
t_ = trials_t[trial_]
if ((int(fs*t_)+win_0)>0) & ((int(fs*t_)+win_1))<len(this_chan):
data_sweep = this_chan[(int(fs*t_)+win_0):(int(fs*t_)+win_1)]
ax.plot(xtime,data_sweep,color=channel_colors[chan_],linewidth=2,alpha=0.5)
if iei > 1/fs:
d = iei*fs
if (lagging_chan_ != 'none') & (len(trial_list)==1):
ax.hlines(thresh_*polarity, xrange[0],xrange[1],linestyle='--',color='green')
ax.vlines(detect_inwin[0],yrange[0],yrange[1],linestyle='--',color='blue')
ax.vlines(detect_inwin[1],yrange[0],yrange[1],linestyle='--',color='blue')
lagging_chan_ = int(lagging_chan_)
lagging_signal = data[(int(fs*t_)+win_0):(int(fs*t_)+win_1),lagging_chan_]
r = find_peaks(lagging_signal*polarity,height=thresh_,distance=d)
mask = (((r[0]/fs)>detect_inwin[0]) & ((r[0]/fs)<detect_inwin[1]))
# lagging_peak_times = [np.round(lt/fs,5) for lt in r[0][mask]]#r[0]/fs
lagging_peak_times = (r[0][mask]/fs) + xrange[0]
lagging_peak_amp = r[1]['peak_heights'][mask]*polarity
lagging_peak_amp = [np.round(a,2) for a in lagging_peak_amp]
lagging_time_readout.value=f'lagging peak times are (ms): {[t*1000 for t in lagging_peak_times]}'
lagging_amp_readout.value=f'lagging peak amplitudes are (V): {lagging_peak_amp}'
if trial_list[0] == 0:
trial_readout.value=f'time since last event is: {NaN}'
trial_abs_readout.value=f'time of this event is: {NaN}'
if trial_list[0] > 0:
iti = 1000*(trials_t[trial_list[0]] - trials_t[trial_list[0]-1])
trial_readout.value=f'time since last event is (msec): {iti:.2f}'
trial_abs_readout.value=f'time of this event is (sec): {trials_t[trial_list[0]]}'
for lt_ in lagging_peak_times:
ax.scatter(lagging_peak_times,lagging_peak_amp,color='black',s=50,zorder=3)
ax.set_ylim(yrange[0],yrange[1]);
ax.set_xlabel('seconds')
ax.xaxis.set_minor_locator(AutoMinorLocator(10))
# # Turn grid on for both major and minor ticks and style minor slightly
# # # differently.
ax.grid(which='major', color='gray', linestyle='-')
ax.grid(which='minor', color='gray', linestyle=':')
# def on_button_export_data(b):
# if len(chan_list)==1:
# df = pd.DataFrame()
# for b_, bout_ in enumerate(bouts_list):
# trials_t = trial_times[(trial_times>bouts_list[bout_][0]) & (trial_times<bouts_list[bout_][1])]
# with output:
# print("Button clicked.")
# lagging_chan_ = int(lagging_chan_)
# lagging_signal = data[(int(fs*t_)+win_0):(int(fs*t_)+win_1),lagging_chan_]
# r = find_peaks(lagging_signal*polarity,height=thresh_,distance=d)
# mask = (((r[0]/fs)>detect_inwin[0]) & ((r[0]/fs)<detect_inwin[1]))
# # lagging_peak_times = [np.round(lt/fs,5) for lt in r[0][mask]]#r[0]/fs
# lagging_peak_times = (r[0][mask]/fs) + xrange[0]
# lagging_peak_amp = r[1]['peak_heights'][mask]*polarity
# lagging_peak_amp = [np.round(a,2) for a in lagging_peak_amp]
# df_ = pd.DataFrame({'amp':amp,'time':time,'bout',b_})
# df.append(df_)
# display(w)
w = interactive_output(update_plot, {'chan_list':select_channels,
'trial_list':select_trials,
'bout_':select_bouts,
'yrange':slider_yrange,
'xrange':slider_xrange,
'lagging_chan_':detect_chan_radio,
'detect_inwin':slider_inwindow,
'thresh_':slider_threshold,
'iei':iei_text,
'polarity':radio_polarity});
display(trial_abs_readout,trial_readout,lagging_time_readout,lagging_amp_readout,ui_trials,ui_peaks,w,ui_range)
#@title {display-mode:"form"}
#@markdown Use the parameters you optimized using the interactive plot to automatically measure peaks
#@markdown across all trials and bouts. The results will be saved as a .csv file.
#@markdown All trials in all bouts are included in the analysis - remove bouts from bouts_list if you do not want them.
condition = 'hot' #@param
channel_use = 0 #@param
polarity_use = -1 #@param
threshold_use = 0.04 #@param
inwindow_use = [0.0050,0.0130] #@param
iei_use = 0.005 #@param
onset_samps = int(inwindow_use[0]*fs)
df = pd.DataFrame()
synwav = []
for b_, bout_ in enumerate(bouts_list):
trials_t = trial_times[(trial_times>bout_[0]) & (trial_times<bout_[1])]
peak_amp = []
peak_time = []
for t_ in trials_t:
lagging_signal = data[(int(fs*(t_+inwindow_use[0]))):(int(fs*(t_+inwindow_use[1]))),channel_use]
r = find_peaks(lagging_signal*polarity_use,height=threshold_use,distance=int(iei_use*fs))
if len(r[0])==0:
peak_time.append(np.NaN)
peak_amp.append(np.NaN)
if len(r[0])>0:
this_i = int(fs*(t_+inwindow_use[0])) + r[0][0]
peak_time.append((r[0][0]/fs)+inwindow_use[0])
peak_amp.append(r[1]['peak_heights'][0]*polarity_use)
synwav.append(data[this_i-int(0.005*fs):this_i+int(0.005*fs),channel_use])
df_ = pd.DataFrame({'amp':peak_amp,'time':peak_time,'bout':b_,'n':np.arange(0,len(peak_amp)),'condition':condition})
df = pd.concat([df,df_])
df.to_csv('evoked-emg_' + condition + '.csv')
#@title {display-mode:"form"}
#@markdown Run this code cell to get a plot of the average EMG waveform (centered at the peak of each event).
#@markdown The variable values from the last code cell will automatically be used and the plot will include 5msec before and after the peak.
#@markdown <br> You can change the y-axis limits so that the y-scale is the same across conditions.
ylim = [-0.25,0.15] #@param
#@markdown Finally, specify whether you want to save the resulting figure or not to the workspace.
save_fig = False #@param
pre_post_dur = 0.005
synwav = []
for b_, bout_ in enumerate(bouts_list):
trials_t = trial_times[(trial_times>bout_[0]) & (trial_times<bout_[1])]
peak_amp = []
peak_time = []
for t_ in trials_t:
lagging_signal = data[(int(fs*(t_+inwindow_use[0]))):(int(fs*(t_+inwindow_use[1]))),channel_use]
r = find_peaks(lagging_signal*polarity_use,height=threshold_use,distance=int(iei_use*fs))
if len(r[0])==0:
peak_time.append(np.NaN)
peak_amp.append(np.NaN)
if len(r[0])>0:
this_i = int(fs*(t_+inwindow_use[0])) + r[0][0]
peak_time.append((r[0][0]/fs)+inwindow_use[0])
peak_amp.append(r[1]['peak_heights'][0]*polarity_use)
synwav.append(data[this_i-int(0.005*fs):this_i+int(0.005*fs),channel_use])
synwav = np.asarray(synwav).T
wav_u = np.mean(synwav,1)
wav_std = np.std(synwav,1)
hfig,ax = plt.subplots(figsize=(5,4))
x = np.linspace(-pre_post_dur,pre_post_dur,np.shape(wav_u)[0])
ax.plot(x,wav_u,color = 'black');
ax.fill_between(x, wav_u-wav_std, wav_u+wav_std, alpha = 0.25, color = 'purple');
ax.set_ylim(ylim[0],ylim[1]);
ax.set_ylabel('Volts (recorded)')
ax.set_xlabel('time (seconds)');
if save_fig == True:
plt.savefig(condition + '_avg.png',format='png',facecolor='white')
Visualize the processed data¶
Plot the processed data from one condition¶
Draw a scatter plot to show the relationship between x and y for different subsets of the data using the hue parameter to identify different bouts (or other subset groupings of the data if you want to change it). It is possible to show up to three dimensions independently by using all three semantic types, but this style of plot can be hard to interpret and is often ineffective. Using redundant semantics (i.e. both hue and style for the same variable) can be helpful for making graphics more accessible.
See the tutorial for more information.
The default treatment of the hue (and to a lesser extent, size) semantic, if present, depends on whether the variable is inferred to represent “numeric” or “categorical” data. In particular, numeric variables are represented with a sequential colormap by default, and the legend entries show regular “ticks” with values that may or may not exist in the data. This behavior can be controlled through various parameters, as described and illustrated below.
#@title {display-mode:"form"}
#@markdown Specify the filepath to a csv file
filepath = 'evoked-emg_hot.csv' #@param
#@markdown Specify the header name of the column you want for your x points.
#@markdown If more than one header is specified (separated by commas), each will be plotted overlaid in a different color for a violin plot.
x_column = 'column header for x axis' #@param
#@markdown Specify the header name of the column you want for your y points.
#@markdown If more than one header is specified (separated by commas), each will be plotted overlaid in a different color for a scatter plot
y_column = 'column header for y axis' #@param
hue_column = 'bout' #@param
#@markdown Specify the plot type ('scatter', 'box' or 'point').
plot_type = 'scatter' #@param
#@markdown Finally, specify whether you want to save the resulting figure or not to the workspace.
save_fig = False #@param
df = pd.read_csv(filepath)
hfig,ax = plt.subplots(figsize=(5,4))
if plot_type == 'scatter':
if hue_column == 'none':
g = sns.scatterplot(data=df,x=x_column,y=y_column,alpha=0.75);
else:
g = sns.scatterplot(data=df,x=x_column,y=y_column,hue=hue_column,style=hue_column,alpha=0.75);
if plot_type == 'point':
g = sns.pointplot(data=df,x=x_column,y=y_column,alpha=0.75);
if plot_type == 'box':
if hue_column == 'none':
g = sns.boxplot(data=df,x=x_column,y=y_column);
else:
g = sns.boxplot(data=df,x=x_column,y=y_column,hue=hue_column);
g.legend(loc='center left', bbox_to_anchor=(1, 0.5), ncol=1);
# ax.set_xlim(x_lim[0],x_lim[1])
if save_fig == True:
plt.savefig(condition + '_' + y_column + '.png',format='png',facecolor='white')
Compare data across conditions¶
Show point estimates and errors using dot marks.
A point plot represents an estimate of central tendency for a numeric variable by the position of the dot and provides some indication of the uncertainty around that estimate using error bars.
Point plots can be more useful than bar plots for focusing comparisons between different levels of one or more categorical variables. They are particularly adept at showing interactions: how the relationship between levels of one categorical variable changes across levels of a second categorical variable. The lines that join each point from the same hue level allow interactions to be judged by differences in slope, which is easier for the eyes than comparing the heights of several groups of points or bars.
It is important to keep in mind that a point plot shows only the mean (or other estimator) value, but in many cases it may be more informative to show the distribution of values at each level of the categorical variables. In that case, other approaches such as a box or violin plot may be more appropriate.
In the following code cell, you can plot any column from the csv files against the stimulus number within each bout. The “hue value” will be condition (so that each point in the plot shows the average across bouts for each stimulus number within a condition).
#@title {display-mode: "form"}
#@markdown List the full path to each csv file created using the "adaptation analysis tool."
file_list = ['evoked-emg_control.csv', 'evoked-emg_cold.csv', 'evoked-emg_hot.csv'] #@param
# file_list = ['adaptation_2_b_.csv','adaptation_5_b_.csv','adaptation_10_b_.csv']
#@markdown Then run this code cell to create a dataframe combining the data from all listed files. <br>
#@markdown And make a "point plot" across conditions for the specified csv column.
y_column = 'time' #@param
#@markdown Finally, specify whether you want to save the resulting figure or not to the workspace.
save_fig = False #@param
df = pd.DataFrame({})
for f in file_list:
df_ = pd.read_csv(f)
df = pd.concat([df,df_])
hfig,ax = plt.subplots(figsize=(8,4))
g = sns.pointplot(data=df,x='n',y=y_column,hue='condition',alpha=0.75);
g.legend(loc='center left', bbox_to_anchor=(1, 0.5), ncol=1);
if save_fig == True:
plt.savefig('all-conditions_' + y_column + '.png',format='png',facecolor='white')
Written by Dr. Krista Perks for courses taught at Wesleyan University.