diff --git a/actions/stimulus-response-fwhm/attributes.yaml b/actions/stimulus-response-fwhm/attributes.yaml new file mode 100644 index 000000000..5adfa43fc --- /dev/null +++ b/actions/stimulus-response-fwhm/attributes.yaml @@ -0,0 +1,5 @@ +registered: '2021-02-01T11:52:58' +data: + results: results.csv + notebook: 10-calculate-stimulus-response-fwhm.ipynb + html: 10-calculate-stimulus-response-fwhm.html diff --git a/actions/stimulus-response-fwhm/data/10-calculate-stimulus-response-fwhm.html b/actions/stimulus-response-fwhm/data/10-calculate-stimulus-response-fwhm.html new file mode 100644 index 000000000..562ecf06a --- /dev/null +++ b/actions/stimulus-response-fwhm/data/10-calculate-stimulus-response-fwhm.html @@ -0,0 +1,14226 @@ + + +
+ +Only PSTHs with a significance level below 0.05 are used
+ +%load_ext autoreload
+%autoreload 2
+
import matplotlib.pyplot as plt
+%matplotlib inline
+import spatial_maps as sp
+import septum_mec.analysis.data_processing as dp
+import septum_mec.analysis.registration
+import expipe
+import os
+import pathlib
+import numpy as np
+import exdir
+import pandas as pd
+import optogenetics as og
+import quantities as pq
+import shutil
+from distutils.dir_util import copy_tree
+import seaborn as sns
+from functools import reduce
+from septum_mec.analysis.stimulus_response import stimulus_response_latency, compute_response
+import scipy.stats
+from tqdm import tqdm_notebook as tqdm
+from tqdm._tqdm_notebook import tqdm_notebook
+from septum_mec.analysis.statistics import load_data_frames, make_paired_tables, make_statistics_table
+tqdm_notebook.pandas()
+
+%matplotlib widget
+
std_gaussian_kde = 0.04
+window_size = 0.03
+tmin = 0. # delay to compute FWHM after stimulus onset
+
data_loader = dp.Data()
+actions = data_loader.actions
+project = data_loader.project
+
output = pathlib.Path('output/stimulus-response-fwhm')
+(output / 'data').mkdir(parents=True, exist_ok=True)
+(output / 'figures').mkdir(parents=True, exist_ok=True)
+(output / 'statistics').mkdir(parents=True, exist_ok=True)
+
colors = ['#d95f02','#e7298a']
+labels = [
+ '11 Hz',
+ '30 Hz'
+]
+queries = [
+ 'frequency==11 and stim_location=="ms"',
+ 'frequency==30 and stim_location=="ms"']
+
data, labels, colors, queries = load_data_frames(queries, labels, colors)
+
# use baseline==False and p_e_peak < 0.05
+session_units_sig = data.query('p_e_peak == p_e_peak')
+
session_units_sig
+
def compute_fwhm(psth, times, t_min=None):
+ """
+ Compute PSTH. The function handles a single peak or multiple peaks. In the latter case, the peak containing the max peak is used
+ """
+ times_accepted = times > tmin
+ max_idx = np.argmax(psth[times_accepted])
+ idxs_greater = np.where(psth[times_accepted] > 0.5 * np.ptp(psth[times_accepted]))[0]
+
+ if np.all(np.diff(idxs_greater) == 1):
+ fwhm = times[times_accepted][idxs_greater[-1]] - times[times_accepted][idxs_greater[0]]
+ else:
+ # deal with multiple peaks. when multiple peaks are found, the one containing max_idx is used
+ diff_splits = np.where((np.diff(idxs_greater) == 1) == False)[0]
+ idxs_list = []
+ for i, ds in enumerate(diff_splits):
+ if i == 0:
+ idxs_list.append(idxs_greater[:ds+1])
+ elif i < len(diff_splits):
+ idxs_list.append(idxs_greater[diff_splits[i-1]+1:ds+1])
+ idxs_list.append(idxs_greater[ds+1:])
+
+ for idxs in idxs_list:
+ if max_idx in idxs:
+ fwhm = times[times_accepted][idxs[-1]] - times[times_accepted][idxs[0]]
+
+ return fwhm, psth[times_accepted][max_idx]
+
times = np.arange(-0.005, window_size, 1e-4)
+tmin = 0.005
+psths = []
+def process(row):
+ action_id = row['action']
+ channel_group = row['channel_group']
+ unit_name = row['unit_name']
+ name = f'{action_id}_{channel_group}_{unit_name}'
+ spike_times = data_loader.spike_train(action_id, channel_group, unit_name)
+
+ spike_times = np.array(spike_times)
+
+ stim_times = data_loader.stim_times(action_id)
+
+ stim_times = np.array(stim_times)
+
+ _, spikes, kernel, p_e, p_i = stimulus_response_latency(
+ spike_times, stim_times, window_size, std_gaussian_kde)
+
+ psths.append(kernel(times))
+ fwhm, peak = compute_fwhm(kernel(times), times, tmin)
+ return pd.Series({
+ 'fwhm': fwhm,
+ 'peak': peak
+ })
+
session_units_sig = session_units_sig.merge(
+ session_units_sig.progress_apply(process, axis=1),
+ left_index=True, right_index=True)
+
# Quality control
+fig = plt.figure()
+for psth in psths:
+ plt.plot(times, psth, color='C0', alpha=0.3)
+
keys = [
+ 'fwhm',
+ 'peak'
+]
+
+results, labels = make_paired_tables(session_units_sig, keys, queries=queries, labels=labels, cell_types=['gridcell', 'ns_inhibited', 'ns_not_inhibited'])
+
plt.rc('axes', titlesize=12)
+plt.rcParams.update({
+ 'font.size': 12,
+ 'figure.figsize': (1.7, 3),
+ 'figure.dpi': 150
+})
+
def violinplot2(data, xticks, colors):
+ pos = [i * 0.6 for i in range(len(data))]
+ print(pos)
+ violins = plt.violinplot(data, pos, showmeans=True, showextrema=False)
+
+ for i, b in enumerate(violins['bodies']):
+ b.set_color(colors[i])
+ b.set_alpha (0.8)
+
+
+
+ # for i, body in enumerate(violins['cbars']):
+ # body.set_color('C{}'.format(i))
+
+
+
+ for category in ['cbars', 'cmins', 'cmaxes', 'cmedians', 'cmeans']:
+ if category in violins:
+ violins[category].set_color(['k', 'k'])
+ violins[category].set_linewidth(2.0)
+ plt.xticks(pos, xticks, rotation=45)
+ plt.gca().spines['top'].set_visible(False)
+ plt.gca().spines['right'].set_visible(False)
+
for cell_type, result, in results.items():
+ fig = plt.figure()
+ violinplot2(
+ [result['fwhm']['11 Hz'].dropna().values, result['fwhm']['30 Hz'].dropna().values],
+ colors=colors,
+ xticks=["11 Hz ", " 30 Hz"],
+ )
+ plt.title(cell_type)
+ figname = f'{cell_type}-stim-response-fwhm-11-30Hz'
+ fig.savefig(
+ output / 'figures' / f'{figname}.png',
+ bbox_inches='tight', transparent=True)
+ fig.savefig(
+ output / 'figures' / f'{figname}.svg',
+ bbox_inches='tight', transparent=True)
+
for cell_type, result, in results.items():
+ fig = plt.figure()
+ violinplot2(
+ [result['peak']['11 Hz'].dropna().values, result['peak']['30 Hz'].dropna().values],
+ colors=colors,
+ xticks=["11 Hz ", " 30 Hz"],
+ )
+ plt.title(cell_type)
+ figname = f'{cell_type}-stim-response-peak-11-30Hz'
+ fig.savefig(
+ output / 'figures' / f'{figname}.png',
+ bbox_inches='tight', transparent=True)
+ fig.savefig(
+ output / 'figures' / f'{figname}.svg',
+ bbox_inches='tight', transparent=True)
+
stats = {}
+for cell_type, result in results.items():
+ stats[cell_type], _ = make_statistics_table(result, labels)
+
stats['gridcell']
+
for cell_type, stat in stats.items():
+ stat.to_latex(output / "statistics" / f"statistics_{cell_type}.tex")
+ stat.to_csv(output / "statistics" / f"statistics_{cell_type}.csv")
+
action = project.require_action("stimulus-response-fwhm")
+
action.modules['parameters'] = {
+ 'window_size': window_size,
+ 'std_gaussian_kde': std_gaussian_kde,
+ 'tmin': tmin
+}
+
action.data['results'] = 'results.csv'
+session_units_sig.to_csv(action.data_path('results'), index=False)
+
copy_tree(output, str(action.data_path()))
+
septum_mec.analysis.registration.store_notebook(action, "10-calculate-stimulus-response-fwhm.ipynb")
+
+
\n", + " | action | \n", + "baseline | \n", + "entity | \n", + "frequency | \n", + "i | \n", + "ii | \n", + "session | \n", + "stim_location | \n", + "stimulated | \n", + "tag | \n", + "... | \n", + "bs_ctrl | \n", + "ns_inhibited | \n", + "ns_not_inhibited | \n", + "gridcell | \n", + "bs_not_gridcell | \n", + "label | \n", + "label_num | \n", + "query | \n", + "color | \n", + "cell_type | \n", + "
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
13 | \n", + "1839-120619-4 | \n", + "False | \n", + "1839 | \n", + "30.0 | \n", + "False | \n", + "True | \n", + "4 | \n", + "ms | \n", + "True | \n", + "stim ii | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "30 Hz | \n", + "1.0 | \n", + "frequency==30 and stim_location==\"ms\" | \n", + "#e7298a | \n", + "bs_not_gridcell | \n", + "
14 | \n", + "1839-120619-4 | \n", + "False | \n", + "1839 | \n", + "30.0 | \n", + "False | \n", + "True | \n", + "4 | \n", + "ms | \n", + "True | \n", + "stim ii | \n", + "... | \n", + "NaN | \n", + "True | \n", + "False | \n", + "True | \n", + "False | \n", + "30 Hz | \n", + "1.0 | \n", + "frequency==30 and stim_location==\"ms\" | \n", + "#e7298a | \n", + "gridcell | \n", + "
15 | \n", + "1839-120619-4 | \n", + "False | \n", + "1839 | \n", + "30.0 | \n", + "False | \n", + "True | \n", + "4 | \n", + "ms | \n", + "True | \n", + "stim ii | \n", + "... | \n", + "NaN | \n", + "True | \n", + "False | \n", + "False | \n", + "False | \n", + "30 Hz | \n", + "1.0 | \n", + "frequency==30 and stim_location==\"ms\" | \n", + "#e7298a | \n", + "ns_inhibited | \n", + "
16 | \n", + "1839-120619-4 | \n", + "False | \n", + "1839 | \n", + "30.0 | \n", + "False | \n", + "True | \n", + "4 | \n", + "ms | \n", + "True | \n", + "stim ii | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "30 Hz | \n", + "1.0 | \n", + "frequency==30 and stim_location==\"ms\" | \n", + "#e7298a | \n", + "bs_not_gridcell | \n", + "
17 | \n", + "1839-120619-4 | \n", + "False | \n", + "1839 | \n", + "30.0 | \n", + "False | \n", + "True | \n", + "4 | \n", + "ms | \n", + "True | \n", + "stim ii | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "30 Hz | \n", + "1.0 | \n", + "frequency==30 and stim_location==\"ms\" | \n", + "#e7298a | \n", + "bs_not_gridcell | \n", + "
... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "... | \n", + "
1279 | \n", + "1833-010719-2 | \n", + "False | \n", + "1833 | \n", + "11.0 | \n", + "True | \n", + "False | \n", + "2 | \n", + "ms | \n", + "True | \n", + "stim i | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "11 Hz | \n", + "0.0 | \n", + "frequency==11 and stim_location==\"ms\" | \n", + "#d95f02 | \n", + "bs_not_gridcell | \n", + "
1280 | \n", + "1833-010719-2 | \n", + "False | \n", + "1833 | \n", + "11.0 | \n", + "True | \n", + "False | \n", + "2 | \n", + "ms | \n", + "True | \n", + "stim i | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "11 Hz | \n", + "0.0 | \n", + "frequency==11 and stim_location==\"ms\" | \n", + "#d95f02 | \n", + "bs_not_gridcell | \n", + "
1281 | \n", + "1833-010719-2 | \n", + "False | \n", + "1833 | \n", + "11.0 | \n", + "True | \n", + "False | \n", + "2 | \n", + "ms | \n", + "True | \n", + "stim i | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "11 Hz | \n", + "0.0 | \n", + "frequency==11 and stim_location==\"ms\" | \n", + "#d95f02 | \n", + "bs_not_gridcell | \n", + "
1282 | \n", + "1833-010719-2 | \n", + "False | \n", + "1833 | \n", + "11.0 | \n", + "True | \n", + "False | \n", + "2 | \n", + "ms | \n", + "True | \n", + "stim i | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "11 Hz | \n", + "0.0 | \n", + "frequency==11 and stim_location==\"ms\" | \n", + "#d95f02 | \n", + "bs_not_gridcell | \n", + "
1283 | \n", + "1833-010719-2 | \n", + "False | \n", + "1833 | \n", + "11.0 | \n", + "True | \n", + "False | \n", + "2 | \n", + "ms | \n", + "True | \n", + "stim i | \n", + "... | \n", + "NaN | \n", + "False | \n", + "False | \n", + "False | \n", + "True | \n", + "11 Hz | \n", + "0.0 | \n", + "frequency==11 and stim_location==\"ms\" | \n", + "#d95f02 | \n", + "bs_not_gridcell | \n", + "
647 rows × 68 columns
\n", + "