From b618c2d9cab5ca62e8955a353ff966aed30560ee Mon Sep 17 00:00:00 2001 From: efish Date: Thu, 29 Nov 2018 16:47:30 +0100 Subject: [PATCH 1/2] zwischendrin --- code/base_chirps.py | 6 ++- code/base_eod.py | 18 +++++--- code/base_spikes.py | 48 ++++++++++++-------- code/func_chirp.py | 53 ++++++++++++++--------- code/func_spike.py | 13 ++++-- code/order_eff.py | 14 +++--- code/plot_spikesduringbaselineactivity.py | 1 + 7 files changed, 96 insertions(+), 57 deletions(-) diff --git a/code/base_chirps.py b/code/base_chirps.py index 1b26e6c..9ead914 100644 --- a/code/base_chirps.py +++ b/code/base_chirps.py @@ -8,7 +8,7 @@ from IPython import embed data_dir = "../data" -dataset = "2018-11-13-ah-invivo-1" +dataset = "2018-11-13-ad-invivo-1" data = ["2018-11-09-ad-invivo-1", "2018-11-09-ae-invivo-1", "2018-11-09-ag-invivo-1", "2018-11-13-aa-invivo-1", "2018-11-13-ac-invivo-1", "2018-11-13-ad-invivo-1", "2018-11-13-ah-invivo-1", "2018-11-13-ai-invivo-1", "2018-11-13-aj-invivo-1", "2018-11-13-ak-invivo-1", "2018-11-13-al-invivo-1", "2018-11-14-aa-invivo-1", "2018-11-14-ac-invivo-1", "2018-11-14-ad-invivo-1", "2018-11-14-af-invivo-1", "2018-11-14-ag-invivo-1", "2018-11-14-ah-invivo-1", "2018-11-14-ai-invivo-1", "2018-11-14-ak-invivo-1", "2018-11-14-al-invivo-1", "2018-11-14-am-invivo-1", "2018-11-14-an-invivo-1", "2018-11-20-aa-invivo-1", "2018-11-20-ab-invivo-1", "2018-11-20-ac-invivo-1", "2018-11-20-ad-invivo-1", "2018-11-20-ae-invivo-1", "2018-11-20-af-invivo-1", "2018-11-20-ag-invivo-1", "2018-11-20-ah-invivo-1", "2018-11-20-ai-invivo-1"] @@ -41,7 +41,9 @@ for i in sort_df: chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) df_map = map_keys(chirp_spikes) sort_df = sorted(df_map.keys()) -dct_phase = plot_std_chirp(sort_df, df_map, chirp_spikes, chirp_mods) +example = [-50, 200, 400] + +dct_phase = plot_std_chirp(example, df_map, chirp_spikes, chirp_mods) plt.show() plt.close('all') diff --git a/code/base_eod.py b/code/base_eod.py index 86a5ec2..8789b9e 100644 --- a/code/base_eod.py +++ b/code/base_eod.py @@ -14,10 +14,16 @@ time,eod = read_baseline_eod(os.path.join(data_dir, dataset)) zeit = np.asarray(time) -plt.plot(zeit[0:1000], eod[0:1000]) -plt.title('A.lepto EOD', fontsize = 18)#Plottitelk -plt.xlabel('time [ms]', fontsize = 16)#Achsentitel -plt.ylabel('amplitude[mv]', fontsize = 16)#Achsentitel -plt.xticks(fontsize = 14) -plt.yticks(fontsize = 14) + +inch_factor = 2.54 +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) + +plt.plot(zeit[0:1000], eod[0:1000], color = 'darkblue') +plt.title('A.lepto EOD', fontsize = 24)#Plottitel +plt.xlabel('time [ms]', fontsize = 22)#Achsentitel +plt.ylabel('amplitude[mv]', fontsize = 22)#Achsentitel +plt.tick_params(axis='both', which='major', labelsize = 22) +ax.spines['right'].set_visible(False) +ax.spines['top'].set_visible(False) +fig.tight_layout() plt.show() diff --git a/code/base_spikes.py b/code/base_spikes.py index 32ffd68..014b7b8 100644 --- a/code/base_spikes.py +++ b/code/base_spikes.py @@ -12,25 +12,28 @@ data_base = ("2018-11-09-ab-invivo-1", "2018-11-09-ad-invivo-1", "2018-11-13-aa- data_chirps = ("2018-11-09-ad-invivo-1", "2018-11-09-ae-invivo-1", "2018-11-09-ag-invivo-1", "2018-11-13-aa-invivo-1", "2018-11-13-ac-invivo-1", "2018-11-13-ad-invivo-1", "2018-11-13-ah-invivo-1", "2018-11-13-ai-invivo-1", "2018-11-13-aj-invivo-1", "2018-11-13-ak-invivo-1", "2018-11-13-al-invivo-1", "2018-11-14-aa-invivo-1", "2018-11-14-ac-invivo-1", "2018-11-14-ad-invivo-1", "2018-11-14-af-invivo-1", "2018-11-14-ag-invivo-1", "2018-11-14-ah-invivo-1", "2018-11-14-ai-invivo-1", "2018-11-14-ak-invivo-1", "2018-11-14-al-invivo-1", "2018-11-14-am-invivo-1", "2018-11-14-an-invivo-1", "2018-11-20-aa-invivo-1", "2018-11-20-ab-invivo-1", "2018-11-20-ac-invivo-1", "2018-11-20-ad-invivo-1", "2018-11-20-ae-invivo-1", "2018-11-20-af-invivo-1", "2018-11-20-ag-invivo-1", "2018-11-20-ah-invivo-1", "2018-11-20-ai-invivo-1") dataset = "2018-11-13-ad-invivo-1" - +inch_factor = 2.54 #for dataset in data_base: spike_times = read_baseline_spikes(os.path.join(data_dir, dataset)) spike_iv = np.diff(spike_times) x = np.arange(0.001, 0.01, 0.0001) -plt.hist(spike_iv,x) + +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) +plt.hist(spike_iv,x, color = 'darkblue') mu = np.mean(spike_iv) sigma = np.std(spike_iv) cv = sigma/mu -plt.title('A.lepto ISI Histogramm', fontsize = 18) -plt.xlabel('duration ISI[ms]', fontsize = 16) -plt.ylabel('number of ISI', fontsize = 16) - -plt.xticks(fontsize = 14) -plt.yticks(fontsize = 14) +plt.title('A.lepto ISI Histogramm', fontsize = 24) +plt.xlabel('duration ISI[ms]', fontsize = 22) +plt.ylabel('number of ISI', fontsize = 22) +plt.tick_params(axis='both', which='major', labelsize = 22) +ax.spines['right'].set_visible(False) +ax.spines['top'].set_visible(False) +plt.tight_layout() plt.show() @@ -47,6 +50,7 @@ sort_df = sorted(df_map.keys()) dct_rate, over_r = spike_rates(sort_df, df_map, chirp_spikes) plt.figure() +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) ls_mean = plot_df_spikes(sort_df, dct_rate) plt.show() @@ -54,14 +58,17 @@ plt.show() #mittlere Feuerrate einer Frequenz auf Frequenz: -plt.figure() -plt.plot(np.arange(0,len(ls_mean),1),ls_mean) +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) +plt.plot(np.arange(0,len(ls_mean),1),ls_mean, color = 'darkblue') plt.scatter(np.arange(0,len(ls_mean),1), np.ones(len(ls_mean))*over_r, color = 'green') -plt.title('Mean firing rate of a cell for a range of frequency differences', fontsize = 18) +plt.title('Mean firing rate of a cell for a range of frequency differences', fontsize = 24) plt.xticks(np.arange(1,len(sort_df),1), (sort_df)) -plt.xlabel('Range of frequency differences [Hz]', fontsize = 16) -plt.ylabel('Mean firing rate of the cell', fontsize = 16) -plt.tick_params(axis='both', which='major', labelsize = 14) +plt.xlabel('Range of frequency differences [Hz]', fontsize = 22) +plt.ylabel('Mean firing rate of the cell', fontsize = 22) +plt.tick_params(axis='both', which='major', labelsize = 18) +ax.spines['right'].set_visible(False) +ax.spines['top'].set_visible(False) +plt.tight_layout() plt.show() @@ -70,12 +77,15 @@ plt.show() #wie viel Prozent der Anfangsrate macht die Adaption von Zellen aus? adapt = adaptation_df(sort_df, dct_rate) -plt.figure() +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) plt.boxplot(adapt) -plt.title('Adaptation of cell firing rate during a trial', fontsize = 18) -plt.xlabel('Cell', fontsize = 16) -plt.ylabel('Adaptation size [Hz]', fontsize = 16) -plt.tick_params(axis='both', which='major', labelsize = 14) +plt.title('Adaptation of cell firing rate during a trial', fontsize = 24) +plt.xlabel('Cell', fontsize = 22) +plt.ylabel('Adaptation size [Hz]', fontsize = 22) +plt.tick_params(axis='both', which='major', labelsize = 18) +ax.spines['right'].set_visible(False) +ax.spines['top'].set_visible(False) +plt.tight_layout() plt.show() diff --git a/code/func_chirp.py b/code/func_chirp.py index b358e3a..5e2bddb 100644 --- a/code/func_chirp.py +++ b/code/func_chirp.py @@ -2,8 +2,10 @@ from read_baseline_data import * from read_chirp_data import * from utility import * import matplotlib.pyplot as plt +import math import numpy as np +inch_factor = 2.54 def chirp_eod_plot(df_map, eod, times): #die äußere Schleife geht für alle Keys durch und somit durch alle dfs @@ -11,7 +13,7 @@ def chirp_eod_plot(df_map, eod, times): for i in df_map.keys(): freq = list(df_map[i]) - fig,axs = plt.subplots(2, 2, sharex = True, sharey = True) + fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) for idx, k in enumerate(freq): ct = times[k] @@ -19,22 +21,28 @@ def chirp_eod_plot(df_map, eod, times): zeit = e1[0] eods = e1[1] - if idx <= 3: - axs[0, 0].plot(zeit, eods, color= 'blue', linewidth = 0.25) - axs[0, 0].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) - elif 4<= idx <= 7: - axs[0, 1].plot(zeit, eods, color= 'blue', linewidth = 0.25) - axs[0, 1].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) - elif 8<= idx <= 11: - axs[1, 0].plot(zeit, eods, color= 'blue', linewidth = 0.25) - axs[1, 0].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) - else: - axs[1, 1].plot(zeit, eods, color= 'blue', linewidth = 0.25) - axs[1, 1].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) + if idx <= 1: + ax.plot(zeit, eods, color= 'darkblue') + ax.scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) + #elif 4<= idx <= 7: + # axs[0, 1].plot(zeit, eods, color= 'blue', linewidth = 0.25) + # axs[0, 1].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) + #elif 8<= idx <= 11: + # axs[1, 0].plot(zeit, eods, color= 'blue', linewidth = 0.25) + # axs[1, 0].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) + else: + continue + #axs[1, 1].plot(zeit, eods, color= 'blue', linewidth = 0.25) + #axs[1, 1].scatter(np.asarray(ct), np.ones(len(ct))*np.mean(eods), color = 'green', s= 22) + + ax.set_ylabel('Amplitude [mV]', fontsize = 22) + ax.set_xlabel('Time [ms]', fontsize = 22) + ax.tick_params(axis='both', which='major', labelsize = 18) + ax.spines['right'].set_visible(False) + ax.spines['top'].set_visible(False) + fig.suptitle('EOD for chirps', fontsize = 24) + fig.tight_layout() - axs[0,1].set_ylabel('Amplitude [mV]') - axs[1,0].set_xlabel('Time [ms]') - fig.suptitle('EOD for chirps', fontsize = 16) @@ -69,7 +77,7 @@ def cut_chirps(freq, eod, times): def plot_std_chirp(sort_df, df_map, chirp_spikes, chirp_mods): - plt.figure() + fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) dct_phase = {} num_bin = 12 phase_vec = np.arange(0, 1+1/num_bin, 1/num_bin) @@ -82,10 +90,13 @@ def plot_std_chirp(sort_df, df_map, chirp_spikes, chirp_mods): dct_phase[i].append(phase[1]) for i in sort_df: - plt.scatter(dct_phase[i], chirp_mods[i], label = i) - plt.title('Change of std depending on the phase where the chirp occured') - plt.xlabel('Phase') - plt.ylabel('Standard deviation of the amplitude modulation') + norm = np.asarray(dct_phase[i]) *2*math.pi + plt.scatter(norm, chirp_mods[i], label = i, s = 22) + plt.title('Change of std depending on the phase where the chirp occured', fontsize = 24) + plt.xlabel('Phase', fontsize = 22) + plt.ylabel('Standard deviation of the amplitude modulation', fontsize = 22) + plt.xticks([0, math.pi/2, math.pi, math.pi*1.5, math.pi*2], ('0', '$\pi$ /2', '$\pi$', '1.5 $\pi$', '2$\pi$')) + plt.tick_params(axis='both', which='major', labelsize = 18) plt.legend() return(dct_phase) diff --git a/code/func_spike.py b/code/func_spike.py index ac47c8f..458f9d5 100644 --- a/code/func_spike.py +++ b/code/func_spike.py @@ -43,6 +43,8 @@ def spike_rates(sort_df, df_map, chirp_spikes): def plot_df_spikes(sort_df, dct_rate): #gibt die Feuerrate gegen die Frequenz aufgetragen + inch_factor = 2.54 + fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) ls_mean = [] for h in sort_df: mean = np.mean(dct_rate[h]) @@ -50,10 +52,13 @@ def plot_df_spikes(sort_df, dct_rate): plt.plot(np.arange(0,len(dct_rate[h]),1),dct_rate[h], label = h) plt.legend() - plt.title('Firing rate of the cell for all trials, sorted by df', fontsize = 18) - plt.xlabel('# of trials', fontsize = 16) - plt.ylabel('Instant firing rate of the cell', fontsize = 16) - plt.tick_params(axis='both', which='major', labelsize = 14) + plt.title('Firing rate of the cell for all trials, sorted by df', fontsize = 24) + plt.xlabel('# of trials', fontsize = 22) + plt.ylabel('Instant firing rate of the cell', fontsize = 22) + plt.tick_params(axis='both', which='major', labelsize = 18) + ax.spines['right'].set_visible(False) + ax.spines['top'].set_visible(False) + plt.tight_layout() return(ls_mean) diff --git a/code/order_eff.py b/code/order_eff.py index 35f14e1..1657ef3 100644 --- a/code/order_eff.py +++ b/code/order_eff.py @@ -12,6 +12,7 @@ data_chirps = ("2018-11-09-ad-invivo-1", "2018-11-09-ae-invivo-1", "2018-11-09-a +inch_factor = 2.54 data_rate_dict = {} for dataset in data_chirps: @@ -31,14 +32,17 @@ for dataset in data_chirps: rate = len(spikes)/ 1.2 data_rate_dict[dataset].append(rate) - +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) for dataset in data_rate_dict: plt.plot(data_rate_dict[dataset]) -plt.title('Test for sequence effects', fontsize = 20) -plt.xlabel('Number of stimulus presentations', fontsize = 18) -plt.ylabel('Firing rates of cells', fontsize = 18) -plt.tick_params(axis='both', which='major', labelsize = 16) +plt.title('Test for sequence effects', fontsize = 24) +plt.xlabel('Number of stimulus presentations', fontsize = 22) +plt.ylabel('Firing rates of cells', fontsize = 22) +plt.tick_params(axis='both', which='major', labelsize = 22) +ax.spines['right'].set_visible(False) +ax.spines['top'].set_visible(False) +fig.tight_layout() plt.show() diff --git a/code/plot_spikesduringbaselineactivity.py b/code/plot_spikesduringbaselineactivity.py index a13c16a..d9c2e49 100644 --- a/code/plot_spikesduringbaselineactivity.py +++ b/code/plot_spikesduringbaselineactivity.py @@ -17,6 +17,7 @@ interspikeintervals = np.diff(spikes)*1000 fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) plt.hist(interspikeintervals, bins=np.arange(0, np.max(interspikeintervals), 0.0001), color='darkblue') +#Titel fehlt!! plt.xlabel("time [ms]", fontsize = 22) plt.xticks(fontsize = 18) plt.ylabel("Number of \n Interspikeinterval", fontsize = 22) From a7217fb4c651b4796a018505c616594b9ce05a8e Mon Sep 17 00:00:00 2001 From: Ramona Date: Thu, 29 Nov 2018 17:06:17 +0100 Subject: [PATCH 2/2] wuhu --- code/plot_eodform_spikehist.py | 4 +- code/response_beat.py | 61 +++++++++++++++----- code/spikes_analysis.py | 56 +++++++++--------- code/spikes_beat.py | 62 ++++++++++++++++++++ code/spikes_chirp.py | 102 +++++++++++++++++++++++++++++++++ 5 files changed, 242 insertions(+), 43 deletions(-) create mode 100644 code/spikes_beat.py create mode 100644 code/spikes_chirp.py diff --git a/code/plot_eodform_spikehist.py b/code/plot_eodform_spikehist.py index 8452399..9512ca6 100644 --- a/code/plot_eodform_spikehist.py +++ b/code/plot_eodform_spikehist.py @@ -9,8 +9,8 @@ from IPython import embed inch_factor = 2.54 sampling_rate = 40000 data_dir = '../data' -dataset = '2018-11-09-ad-invivo-1' -#dataset = '2018-11-13-aa-invivo-1' +#dataset = '2018-11-09-ad-invivo-1' +dataset = '2018-11-14-ad-invivo-1' # read eod and time of baseline time, eod = read_baseline_eod(os.path.join(data_dir, dataset)) diff --git a/code/response_beat.py b/code/response_beat.py index 2d28c12..e899e6b 100644 --- a/code/response_beat.py +++ b/code/response_beat.py @@ -8,29 +8,33 @@ from IPython import embed # define data path and important parameters data_dir = "../data" sampling_rate = 40 #kHz -cut_window = 40 +cut_window = 100 cut_range = np.arange(-cut_window * sampling_rate, 0, 1) window = 1 +''' # norm: -150, 150, 300 aa, #ac, aj?? -data = ["2018-11-13-aa-invivo-1"]#, "2018-11-13-ad-invivo-1", "2018-11-13-ah-invivo-1", "2018-11-13-ai-invivo-1", - #"2018-11-13-ak-invivo-1", "2018-11-13-al-invivo-1"] +data = ["2018-11-13-aa-invivo-1", "2018-11-13-ad-invivo-1", "2018-11-13-ah-invivo-1", "2018-11-13-ai-invivo-1", + "2018-11-13-ak-invivo-1", "2018-11-13-al-invivo-1"] -''' # norm: -50 data = ["2018-11-20-aa-invivo-1", "2018-11-20-ab-invivo-1", "2018-11-20-ac-invivo-1","2018-11-20-ad-invivo-1", "2018-11-20-ae-invivo-1", "2018-11-20-af-invivo-1", "2018-11-20-ag-invivo-1", "2018-11-20-ah-invivo-1", "2018-11-20-ai-invivo-1"] - -data = ["2018-11-14-aa-invivo-1", "2018-11-14-ac-invivo-1", "2018-11-14-ad-invivo-1", "2018-11-14-af-invivo-1", - "2018-11-14-ag-invivo-1", "2018-11-14-ah-invivo-1", "2018-11-14-ai-invivo-1", "2018-11-14-ak-invivo-1", - "2018-11-14-al-invivo-1", "2018-11-14-am-invivo-1", "2018-11-14-an-invivo-1"] ''' + +data = ["2018-11-14-ad-invivo-1", "2018-11-14-ak-invivo-1", "2018-11-14-am-invivo-1"] + +#data = ["2018-11-14-ac-invivo-1", "2018-11-14-ad-invivo-1", "2018-11-14-ak-invivo-1"] #data = ["2018-11-09-ad-invivo-1", "2018-11-14-af-invivo-1"] +#data = ["2018-11-20-ad-invivo-1", "2018-11-13-ad-invivo-1"] +#data = ["2018-11-09-ad-invivo-1"] + rates = {} for dataset in data: + rates[dataset] = {} print(dataset) # read baseline spikes base_spikes = read_baseline_spikes(os.path.join(data_dir, dataset)) @@ -66,21 +70,48 @@ for dataset in data: # also save as binary, 0 no spike, 1 spike binary_spikes = np.isin(cut_range, spikes_idx) * 1 smoothed_data = smooth(binary_spikes, window, 1 / sampling_rate) - train = smoothed_data[window:beat_window+window] - norm_train = train*1000/spikerate - rep_rates.append(np.std(norm_train))#/spikerate) + #train = smoothed_data[window:beat_window+window] + #norm_train = train*1000/spikerate + #df_rate = np.std(norm_train) + #rates[dataset][df] = [df_rate] + rep_rates.append(np.std(smoothed_data))#/spikerate) + ''' + if df in rates[dataset].keys(): + rates[dataset][df].append(np.std(norm_train)) + else: + rates[dataset][df] = [np.std(norm_train)] + ''' break + #break df_rate = np.mean(rep_rates) + #df_rate = rep_rates + rates[dataset][df] = df_rate + #embed() #exit() + ''' if df in rates.keys(): - rates[df].append(df_rate) + rates[dataset][df].append(df_rate) else: - rates[df] = [df_rate] + rates[dataset][df] = [df_rate] + ''' + +colors = ['royalblue', 'red', 'green', 'violet', 'orange', 'black', 'gray'] + +fig, ax = plt.subplots() +for i, cell in enumerate(rates.keys()): + for j, df in enumerate(sorted(rates[cell].keys())): + ax.plot(df, rates[cell][df], 'o', color=colors[i]) +#ax.legend(sorted(rates.keys()), loc='upper left', bbox_to_anchor=(1.04, 1)) +fig.tight_layout() +plt.show() +''' fig, ax = plt.subplots() -for i, k in enumerate(sorted(rates.keys())): - ax.plot(np.ones(len(rates[k]))*k, rates[k], 'o') +for i, cell in enumerate(rates.keys()): + for j, df in enumerate(sorted(rates[cell].keys())): + ax.plot(np.ones(len(rates[cell][df]))*df, rates[cell][df], 'o', color=colors[i]) #ax.legend(sorted(rates.keys()), loc='upper left', bbox_to_anchor=(1.04, 1)) fig.tight_layout() plt.show() +''' \ No newline at end of file diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index ea8648e..53750d8 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -8,7 +8,8 @@ from IPython import embed # define sampling rate and data path sampling_rate = 40 #kHz data_dir = "../data" -dataset = "2018-11-13-ah-invivo-1" +dataset = "2018-11-13-al-invivo-1" +#dataset = "2018-11-09-ad-invivo-1" ''' data = ["2018-11-09-ad-invivo-1", "2018-11-09-ae-invivo-1", "2018-11-09-ag-invivo-1", "2018-11-13-aa-invivo-1", @@ -90,14 +91,14 @@ for deltaf in df_map.keys(): # make dictionaries for csi and beat -csi_trains = {} +#csi_trains = {} csi_rates = {} -beat = {} +#beat = {} # for plotting and calculating iterate over delta f and phases for df in df_phase_time.keys(): - csi_trains[df] = [] - csi_rates[df] = [] - beat[df] = [] + #csi_trains[df] = {} + csi_rates[df] = {} + #beat[df] = [] beat_duration = int(abs(1/df*1000)*sampling_rate) #steps beat_window = 0 # beat window is at most 20 ms long, multiples of beat_duration @@ -123,9 +124,12 @@ for df in df_phase_time.keys(): std_chirp = np.std(np.mean(train_chirp, axis=0)) std_beat = np.std(np.mean(train_beat, axis=0)) - beat[df].append(std_beat) + #beat[df].append(std_beat) csi_spikerate = (std_chirp - std_beat) / (std_chirp + std_beat) + csi_rates[df][phase] = np.mean(csi_spikerate) + + ''' rcs = [] rbs = [] for i, train in enumerate(train_chirp): @@ -145,9 +149,7 @@ for df in df_phase_time.keys(): # add the csi to the dictionaries with the correct df and phase csi_trains[df].append(csi_train) - csi_rates[df].append(np.mean(csi_spikerate)) - ''' # plot plot_trials = df_phase_time[df][phase] plot_trials_binary = np.mean(df_phase_binary[df][phase], axis=0) @@ -170,37 +172,39 @@ for df in df_phase_time.keys(): plt.show() ''' + +colors = ['k', 'k', 'k', + 'k', 'k', 'k', + 'k', 'k', 'k', + 'k', 'k', 'firebrick'] + +sizes = [12, 12, 12, + 12, 12, 12, + 12, 12, 12, + 12, 12, 18] + upper_limit = np.max(sorted(csi_rates.keys()))+30 lower_limit = np.min(sorted(csi_rates.keys()))-30 fig, ax = plt.subplots() -for i, k in enumerate(sorted(csi_rates.keys())): - ax.scatter(np.ones(len(csi_rates[k]))*k, csi_rates[k], s=20) - #ax.plot(i, np.mean(csi_rates[k]), 'o', markersize=15) -#ax.legend(sorted(csi_rates.keys()), loc='upper left', bbox_to_anchor=(1.04, 1)) ax.plot([lower_limit, upper_limit], np.zeros(2), 'silver', linewidth=2, linestyle='--') -#ax.set_xticklabels(sorted(csi_rates.keys())) +for i, df in enumerate(sorted(csi_rates.keys())): + for j, phase in enumerate(sorted(csi_rates[df].keys())): + ax.plot(df, csi_rates[df][phase], 'o', color=colors[j], ms=sizes[j]) fig.tight_layout() plt.show() ''' fig, ax = plt.subplots() -for i, k in enumerate(sorted(csi_trains.keys())): - ax.plot(np.ones(len(csi_trains[k]))*i, csi_trains[k], 'o') - #ax.plot(i, np.mean(csi_trains[k]), 'o', markersize=15) -ax.legend(sorted(csi_trains.keys()), loc='upper left', bbox_to_anchor=(1.04, 1)) -ax.plot(np.arange(-1, len(csi_trains.keys())+1), np.zeros(len(csi_trains.keys())+2), 'silver', linewidth=2, linestyle='--') -#ax.set_xticklabels(sorted(csi_trains.keys())) +for i, k in enumerate(sorted(beat.keys())): + ax.plot(np.ones(len(beat[k]))*k, beat[k], 'o') fig.tight_layout() plt.show() -''' -''' fig, ax = plt.subplots() -for i, k in enumerate(sorted(beat.keys())): - ax.plot(np.ones(len(beat[k]))*i, beat[k], 'o') -ax.legend(sorted(beat.keys()), loc='upper left', bbox_to_anchor=(1.04, 1)) -#ax.set_xticklabels(sorted(csi_trains.keys())) +for i, k in enumerate(sorted(csi_trains.keys())): + ax.plot(np.ones(len(csi_trains[k]))*i, csi_trains[k], 'o') +ax.plot(np.arange(-1, len(csi_trains.keys())+1), np.zeros(len(csi_trains.keys())+2), 'silver', linewidth=2, linestyle='--') fig.tight_layout() plt.show() ''' \ No newline at end of file diff --git a/code/spikes_beat.py b/code/spikes_beat.py new file mode 100644 index 0000000..c9f124c --- /dev/null +++ b/code/spikes_beat.py @@ -0,0 +1,62 @@ +import matplotlib.pyplot as plt +import numpy as np +from read_chirp_data import * +from read_baseline_data import * +from utility import * +from IPython import embed + +# define data path and important parameters +data_dir = "../data" +sampling_rate = 40 #kHz +cut_window = 100 +cut_range = np.arange(-cut_window * sampling_rate, 0, 1) +window = 1 +#dataset = "2018-11-13-ad-invivo-1" +#dataset = "2018-11-13-aj-invivo-1" +#dataset = "2018-11-13-ak-invivo-1" #al +#dataset = "2018-11-14-ad-invivo-1" +dataset = "2018-11-20-af-invivo-1" + +base_spikes = read_baseline_spikes(os.path.join(data_dir, dataset)) +base_spikes = base_spikes[1000:2000] +spikerate = len(base_spikes) / base_spikes[-1] +print(spikerate) + +# read spikes during chirp stimulation +spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) +df_map = map_keys(spikes) + +rates = {} +# iterate over df +for deltaf in df_map.keys(): + rates[deltaf] = {} + beat_duration = int(abs(1 / deltaf) * 1000) + beat_window = 0 + while beat_window + beat_duration <= cut_window/2: + beat_window = beat_window + beat_duration + for x, repetition in enumerate(df_map[deltaf]): + for phase in spikes[repetition]: + # get spikes some ms before the chirp first chirp + spikes_to_cut = np.asarray(spikes[repetition][phase]) + spikes_cut = spikes_to_cut[(spikes_to_cut > -cut_window) & (spikes_to_cut < 0)] + spikes_idx = np.round(spikes_cut * sampling_rate) + # also save as binary, 0 no spike, 1 spike + binary_spikes = np.isin(cut_range, spikes_idx) * 1 + smoothed_data = smooth(binary_spikes, window, 1 / sampling_rate) + #train = smoothed_data[window*sampling_rate:beat_window*sampling_rate+window*sampling_rate] + modulation = np.std(smoothed_data) + rates[deltaf][x] = modulation + break + +fig, ax = plt.subplots() +for i, df in enumerate(sorted(rates.keys())): + for j, rep in enumerate(rates[df].keys()): + if j == 15: + farbe = 'royalblue' + gro = 18 + else: + farbe = 'k' + gro = 12 + ax.plot(df, rates[df][rep], marker='o', color=farbe, ms=gro) +fig.tight_layout() +plt.show() diff --git a/code/spikes_chirp.py b/code/spikes_chirp.py new file mode 100644 index 0000000..fcdda70 --- /dev/null +++ b/code/spikes_chirp.py @@ -0,0 +1,102 @@ +import matplotlib.pyplot as plt +import numpy as np +from read_chirp_data import * +from utility import * +from IPython import embed + +# define sampling rate and data path +sampling_rate = 40 #kHz +data_dir = "../data" +dataset = "2018-11-13-al-invivo-1" + +# parameters for binning, smoothing and plotting +cut_window = 20 +chirp_duration = 14 #ms +neuronal_delay = 5 #ms +chirp_start = int((-chirp_duration / 2 + neuronal_delay + cut_window * 2) * sampling_rate) #index +chirp_end = int((chirp_duration / 2 + neuronal_delay + cut_window * 2) * sampling_rate) #index +number_bins = 12 +window = 1 #ms +time_axis = np.arange(-cut_window*2, cut_window*2, 1/sampling_rate) #steps +spike_bins = np.arange(-cut_window*2, cut_window*2) #ms + +colors = ['k', 'k', 'k', + 'k', 'k', 'k', + 'k', 'k', 'k', + 'k', 'k', 'firebrick'] + +sizes = [12, 12, 12, + 12, 12, 12, + 12, 12, 12, + 12, 12, 18] + +# differentiate between phases +phase_vec = np.arange(0, 1 + 1 / number_bins, 1 / number_bins) +cut_range = np.arange(-cut_window*2*sampling_rate, cut_window*2*sampling_rate, 1) + +df_phase_binary = {} + +spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) +df_map = map_keys(spikes) + +for deltaf in df_map.keys(): + df_phase_binary[deltaf] = {} + for rep in df_map[deltaf]: + chirp_size = int(rep[-1].strip('Hz')) + if chirp_size == 150: + continue + for phase in spikes[rep]: + for idx in np.arange(number_bins): + # check the phase + if phase[1] > phase_vec[idx] and phase[1] < phase_vec[idx+1]: + + # get spikes between 40 ms before and after the chirp + spikes_to_cut = np.asarray(spikes[rep][phase]) + spikes_cut = spikes_to_cut[(spikes_to_cut > -cut_window*2) & (spikes_to_cut < cut_window*2)] + spikes_idx = np.round(spikes_cut*sampling_rate) + # save as binary, 0 no spike, 1 spike + binary_spikes = np.isin(cut_range, spikes_idx)*1 + + # add the spikes to the dictionary with the correct df and phase + if idx in df_phase_binary[deltaf].keys(): + df_phase_binary[deltaf][idx] = np.vstack((df_phase_binary[deltaf][idx], binary_spikes)) + else: + df_phase_binary[deltaf][idx] = binary_spikes + +csi_rates = {} + +for df in df_phase_binary.keys(): + csi_rates[df] = {} + beat_duration = int(abs(1/df*1000)*sampling_rate) #steps + beat_window = 0 + # beat window is at most 20 ms long, multiples of beat_duration + while beat_window+beat_duration <= cut_window*sampling_rate: + beat_window = beat_window+beat_duration + for phase in df_phase_binary[df].keys(): + # csi calculation + trials_binary = df_phase_binary[df][phase] + + train_chirp = [] + train_beat = [] + for i, trial in enumerate(trials_binary): + smoothed_trial = smooth(trial, window, 1/sampling_rate) + train_chirp.append(smoothed_trial[chirp_start:chirp_end]) + train_beat.append(smoothed_trial[chirp_start-beat_window:chirp_start]) + + std_chirp = np.std(np.mean(train_chirp, axis=0)) + std_beat = np.std(np.mean(train_beat, axis=0)) + csi_spikerate = (std_chirp - std_beat) / (std_chirp + std_beat) + + csi_rates[df][phase] = np.mean(csi_spikerate) + + +upper_limit = np.max(sorted(csi_rates.keys()))+30 +lower_limit = np.min(sorted(csi_rates.keys()))-30 + +fig, ax = plt.subplots() +ax.plot([lower_limit, upper_limit], np.zeros(2), 'silver', linewidth=2, linestyle='--') +for i, df in enumerate(sorted(csi_rates.keys())): + for j, phase in enumerate(sorted(csi_rates[df].keys())): + ax.plot(df, csi_rates[df][phase], 'o', color=colors[j], ms=sizes[j]) +fig.tight_layout() +plt.show()