diff --git a/code/base_chirps.py b/code/base_chirps.py index f19b18e..9ead914 100644 --- a/code/base_chirps.py +++ b/code/base_chirps.py @@ -8,45 +8,50 @@ from IPython import embed data_dir = "../data" -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") +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"] +#for dataset in data: -for dataset in data: - print(dataset) - eod = read_chirp_eod(os.path.join(data_dir, dataset)) - times = read_chirp_times(os.path.join(data_dir, dataset)) - df_map = map_keys(eod) - sort_df = sorted(df_map.keys()) +eod = read_chirp_eod(os.path.join(data_dir, dataset)) +times = read_chirp_times(os.path.join(data_dir, dataset)) +df_map = map_keys(eod) +sort_df = sorted(df_map.keys()) - chirp_eod_plot(df_map, eod, times) +eods = chirp_eod_plot(df_map, eod, times) +plt.show() +plt.close('all') - chirp_mods = [] - beat_mods = [] - for i in sort_df: - freq = list(df_map[i]) - ls_mod, beat_mod = cut_chirps(freq, eod, times) - chirp_mods.append(ls_mod) - beat_mods.append(beat_mod) - +chirp_mods = {} +beat_mods = [] +for i in sort_df: + chirp_mods[i] = [] + freq = list(df_map[i]) + ls_mod, beat_mod = cut_chirps(freq, eod, times) + chirp_mods[i].append(ls_mod) + beat_mods.append(beat_mod) #Chirps einer Phase zuordnen - zusammen plotten - chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) - df_map = map_keys(chirp_spikes) - sort_df = sorted(df_map.keys()) - - #plot_std_chirp(sort_df, df_map, chirp_spikes, chirp_mods) - +chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) +df_map = map_keys(chirp_spikes) +sort_df = sorted(df_map.keys()) +example = [-50, 200, 400] +dct_phase = plot_std_chirp(example, df_map, chirp_spikes, chirp_mods) +plt.show() +plt.close('all') +''' #Vatriablen speichern, die man für die Übersicht aller Zellen braucht - name = str(dataset.strip('invivo-1')) + name = str(dataset.replace('-invivo-1', '')) + print('saving ../results/Chirpcut/Cc_' + name + '.dat') f = open('../results/Chirpcut/Cc_' + name + '.dat' , 'w') f.write(str(sort_df)) f.write(str(df_map)) @@ -56,3 +61,4 @@ for dataset in data: #f.write(str(chirp_mods)) #f.write(str(beat_mods)) f.close() +''' diff --git a/code/base_eod.py b/code/base_eod.py index bf1e757..8789b9e 100644 --- a/code/base_eod.py +++ b/code/base_eod.py @@ -8,14 +8,22 @@ from IPython import embed #Funktionen importieren data_dir = "../data" dataset = "2018-11-09-aa-invivo-1" #data = ("2018-11-09-aa-invivo-1", "2018-11-09-ab-invivo-1", "2018-11-09-ac-invivo-1", "2018-11-09-ad-invivo-1", "2018-11-13-aa-invivo-1", "2018-11-13-ab-invivo-1", "2018-11-13-ad-invivo-1", "2018-11-09-af-invivo-1", "2018-11-09-ag-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-ab-invivo-1", "2018-11-14-ac-invivo-1", "2018-11-14-ad-invivo-1", "2018-11-14-ae-invivo-1", "2018-11-14-af-invivo-1", "2018-11-14-ag-invivo-1", "2018-11-14-aa-invivo-1", "2018-11-14-aj-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") + + 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')#Plottitelk -plt.xlabel('time [ms]', fontsize = 12)#Achsentitel -plt.ylabel('amplitude[mv]', fontsize = 12)#Achsentitel -plt.xticks(fontsize = 12) -plt.yticks(fontsize = 12) + +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 9f8bbc2..014b7b8 100644 --- a/code/base_spikes.py +++ b/code/base_spikes.py @@ -3,7 +3,7 @@ from read_chirp_data import * from func_spike import * import matplotlib.pyplot as plt import numpy as np -from IPython import embed #Funktionen imposrtieren +from IPython import embed #Funktionen importieren @@ -11,77 +11,89 @@ data_dir = "../data" data_base = ("2018-11-09-ab-invivo-1", "2018-11-09-ad-invivo-1", "2018-11-13-aa-invivo-1", "2018-11-13-ab-invivo-1", "2018-11-13-ad-invivo-1", "2018-11-13-af-invivo-1", "2018-11-13-ag-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-ab-invivo-1", "2018-11-14-ac-invivo-1", "2018-11-14-ad-invivo-1", "2018-11-14-ae-invivo-1", "2018-11-14-af-invivo-1", "2018-11-14-ag-invivo-1", "2018-11-14-aj-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-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_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: -''' -for dataset in data_base: - - print(dataset) - 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) +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) - mu = np.mean(spike_iv) - sigma = np.std(spike_iv) - cv = sigma/mu +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) +plt.hist(spike_iv,x, color = 'darkblue') - plt.title('A.lepto ISI Histogramm', fontsize = 14) - plt.xlabel('duration ISI[ms]', fontsize = 12) - plt.ylabel('number of ISI', fontsize = 12) +mu = np.mean(spike_iv) +sigma = np.std(spike_iv) +cv = sigma/mu - plt.xticks(fontsize = 12) - plt.yticks(fontsize = 12) -''' +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() -for dataset in data_chirps: +#for dataset in data_chirps: #Nyquist-Theorem Plot: - print(dataset) - chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) - times = read_chirp_times(os.path.join(data_dir, dataset)) - eod = read_chirp_eod(os.path.join(data_dir, dataset)) - df_map = map_keys(chirp_spikes) - sort_df = sorted(df_map.keys()) - - dct_rate, over_r = spike_rates(sort_df, df_map, chirp_spikes) +chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) +times = read_chirp_times(os.path.join(data_dir, dataset)) +eod = read_chirp_eod(os.path.join(data_dir, dataset)) +df_map = map_keys(chirp_spikes) +sort_df = sorted(df_map.keys()) - plt.figure() - ls_mean = plot_df_spikes(sort_df, dct_rate) +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() #mittlere Feuerrate einer Frequenz auf Frequenz: - plt.figure() - plt.plot(np.arange(0,len(ls_mean),1),ls_mean) - plt.scatter(np.arange(0,len(ls_mean),1), np.ones(len(ls_mean))*over_r) - plt.title('Mean firing rate of a cell for a range of frequency differences') - plt.xticks(np.arange(1,len(sort_df),1), (sort_df)) - plt.xlabel('Range of frequency differences [Hz]') - plt.ylabel('Mean firing rate of the cell') - +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 = 24) +plt.xticks(np.arange(1,len(sort_df),1), (sort_df)) +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() #Adaption der Zellen: #wie viel Prozent der Anfangsrate macht die Adaption von Zellen aus? - adapt = adaptation_df(sort_df, dct_rate) - plt.figure() - plt.boxplot(adapt) - plt.title('Adaptation of cell firing rate during a trial') - plt.xlabel('Cell') - plt.ylabel('Adaptation size [Hz]') +adapt = adaptation_df(sort_df, dct_rate) +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 = 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() +''' #Vatriablen speichern, die man für die Übersicht aller Zellen braucht - name = str(dataset.strip('invivo-1')) + name = str(dataset.replace('-invivo-1', '')) f = open('../results/Nyquist/Ny_' + name + '.txt' , 'w') f.write(str(sort_df)) f.write(str(df_map)) @@ -91,3 +103,4 @@ for dataset in data_chirps: f.write(str(over_r)) f.write(str(adapt)) f.close() +''' diff --git a/code/eod_chirp_beat.py b/code/eod_chirp_beat.py new file mode 100644 index 0000000..072cca7 --- /dev/null +++ b/code/eod_chirp_beat.py @@ -0,0 +1,96 @@ +import numpy as np +import matplotlib.pyplot as plt +import thunderfish.peakdetection as pd + + +def create_chirp(eodf): + stimulusrate = eodf # the eod frequency of the fake fish + currentchirptimes = [0.0] + chirpwidth = 0.014 # ms + chirpsize = 100. + chirpampl = 0.02 + chirpkurtosis = 1. + p = 0. + stepsize = 0.00001 + + + time = np.arange(-0.05, 0.05, stepsize) + signal = np.zeros(time.shape) + ampl = np.ones(time.shape) + freq = np.ones(time.shape) + + ck = 0 + csig = 0.5 * chirpwidth / np.power(2.0*np.log(10.0), 0.5/chirpkurtosis) + + for k, t in enumerate(time): + a = 1. + f = stimulusrate + if ck < len(currentchirptimes): + if np.abs(t - currentchirptimes[ck]) < 2.0 * chirpwidth: + x = t - currentchirptimes[ck] + g = np.exp(-0.5 * (x/csig)**2) + f = chirpsize * g + stimulusrate + a *= 1.0 - chirpampl * g + elif t > currentchirptimes[ck] + 2.0 * chirpwidth: + ck += 1 + freq[k] = f + ampl[k] = a + p += f * stepsize + signal[k] = a * np.sin(6.28318530717959 * p) + + return time, signal + + +def plot_chirp(eodf, eodf1, phase, axis): + time, chirp_eod = create_chirp(eodf) + eod = np.sin(time * 2 * np.pi * eodf1 + phase) + + y = chirp_eod * 0.4 + eod + p, t = pd.detect_peaks(y, 0.1) + axis.plot(time*1000, y, color = 'royalblue') + axis.plot(time[p]*1000, (y)[p], lw=2, color='k') + axis.plot(time[t]*1000, (y)[t], lw=2, color='k') + axis.spines["top"].set_visible(False) + axis.spines["right"].set_visible(False) + + + +inch_factor = 2.54 + +fig = plt.figure(figsize=(20 / inch_factor, 10 / inch_factor)) +ax1 = fig.add_subplot(221) +ax2 = fig.add_subplot(222) +ax3 = fig.add_subplot(223) +ax4 = fig.add_subplot(224) + +plot_chirp(600, 650, 0, ax2) +plot_chirp(600, 650, np.pi, ax4) + +plot_chirp(600, 620, 0, ax1) +plot_chirp(600, 620, np.pi, ax3) + +ax1.set_ylabel('EOD [mV]', fontsize=22) +ax1.set_title('$\Delta$f = 20 Hz', fontsize = 18) +ax1.yaxis.set_tick_params(labelsize=18) +ax1.set_xticklabels([]) + +ax2.set_title('$\Delta$f = 50 Hz', fontsize = 18) +ax2.set_xticklabels([]) +ax2.set_yticklabels([]) +ax3.set_ylabel('EOD [mV]', fontsize=22) +ax3.xaxis.set_tick_params(labelsize=18) +ax3.yaxis.set_tick_params(labelsize=18) + +ax3.set_xlabel('time [ms]', fontsize=22) +ax4.set_xlabel('time [ms]', fontsize=22) +ax4.xaxis.set_tick_params(labelsize=18) +ax4.set_yticklabels([]) + + + + + + +fig.tight_layout() +#plt.show() +plt.savefig('chirps_while_beat.png') diff --git a/code/eod_freq_beat.py b/code/eod_freq_beat.py index 8dab8e7..04ed6e7 100644 --- a/code/eod_freq_beat.py +++ b/code/eod_freq_beat.py @@ -2,38 +2,43 @@ from read_baseline_data import * from IPython import embed import matplotlib.pyplot as plt import numpy as np +import thunderfish.peakdetection as pd +from IPython import embed ## beat data_dir = '../data' dataset = '2018-11-09-ad-invivo-1' +inch_factor = 2.54 time, eod = read_baseline_eod(os.path.join(data_dir, dataset)) eod_norm = eod - np.mean(eod) +eod_norm = eod_norm[10000:20000] + +x = np.arange(0., len(eod_norm)) -# calculate eod times and indices by zero crossings -threshold = 0 -shift_eod = np.roll(eod_norm, 1) -eod_times = time[(eod_norm >= threshold) & (shift_eod < threshold)] +y = np.sin(time[10000:20000]*2*np.pi*600)*0.5 +ampl = eod_norm + y +p, t = pd.detect_peaks(ampl, 0.1) -#x = eod_times*40000 -x = np.arange(0., len(eod_times)-1) -y = np.sin(x*2*np.pi*600) -eod_freq_beat = 1/(np.diff(eod_times) + y) +time_axis = np.arange(len(ampl)) -# glätten -kernel = np.ones(7)/7 -smooth_eod_freq_beat = np.convolve(eod_freq_beat, kernel, mode = 'valid') -time_axis = np.arange(len(smooth_eod_freq_beat)) -plt.plot(time_axis, smooth_eod_freq_beat) +fig, ax = plt.subplots(figsize=(20/inch_factor, 10/inch_factor)) +plt.plot(time[10000:20000], ampl) +plt.plot(time[10000:20000][p], ampl[p], lw=2, color='k') +plt.plot(time[10000:20000][t], ampl[t], lw=2, color='k') +ax.set_xlabel("time [ms]", fontsize = 22) +plt.xticks(fontsize = 18) +ax.set_ylabel("eod amplitude [mV]", fontsize = 22) +plt.yticks(fontsize = 18) +ax.spines["top"].set_visible(False) +ax.spines["right"].set_visible(False) +fig.tight_layout() + plt.show() +#plt.savefig('beat.png') + -#eod_freq_beat = eod_freq_normal + y -#smooth_eod_freq_beat = np.convolve(eod_freq_beat, kernel, mode = 'valid') -#fig = plt.plot(time_axis,smooth_eod_freq_beat) -#plt.xlabel("time [ms]") -#plt.ylabel("eod frequency [mV]") -#plt.show() diff --git a/code/func_chirp.py b/code/func_chirp.py index 836b52d..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,25 +21,29 @@ 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))*3, 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))*3, 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))*3, 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))*3, color = 'green', s= 22) - - fig.suptitle('EOD for chirps', fontsize = 16) - axs[0,0].set_ylabel('Amplitude [mV]') - axs[0,1].set_xlabel('Amplitude [mV]') - axs[1,0].set_xlabel('Time [ms]') - axs[1,1].set_xlabel('Time [ms]') - plt.close() + 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() + + @@ -60,16 +66,18 @@ def cut_chirps(freq, eod, times): ls_beat.extend(beat_cut) beat_mod = np.std(ls_beat) #Std vom Bereich vor dem Chirp - plt.figure() - plt.scatter(np.arange(0,len(ls_mod),1), ls_mod) - plt.scatter(np.arange(0,len(ls_mod),1), np.ones(len(ls_mod))*beat_mod, color = 'violet') - plt.close() + #plt.figure() + #plt.scatter(np.arange(0,len(ls_mod),1), ls_mod) + #plt.scatter(np.arange(0,len(ls_mod),1), np.ones(len(ls_mod))*beat_mod, color = 'violet') return(ls_mod, beat_mod) -def plot_std_chirp(sort_df, df_map, chirp_spikes, ls_mod): - plt.figure() + + +def plot_std_chirp(sort_df, df_map, chirp_spikes, chirp_mods): + + 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) @@ -80,8 +88,15 @@ def plot_std_chirp(sort_df, df_map, chirp_spikes, ls_mod): for k in freq: for phase in chirp_spikes[k]: dct_phase[i].append(phase[1]) - - plt.scatter(dct_phase[i], ls_mod[i]) - plt.title('Change of std depending on the phase where the chirp occured') - plt.close() + + for i in sort_df: + 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 d512986..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,9 +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') - plt.xlabel('# of trials') - plt.ylabel('Instant firing rate of the cell') + 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 new file mode 100644 index 0000000..1657ef3 --- /dev/null +++ b/code/order_eff.py @@ -0,0 +1,48 @@ +from read_chirp_data import * +from func_spike import * +import matplotlib.pyplot as plt +import numpy as np +from IPython import embed #Funktionen importieren + + + + +data_dir = "../data" +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") + + + +inch_factor = 2.54 +data_rate_dict = {} +for dataset in data_chirps: + + data_rate_dict[dataset] = [] + chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) + times = read_chirp_times(os.path.join(data_dir, dataset)) + eod = read_chirp_eod(os.path.join(data_dir, dataset)) + df_map = map_keys(chirp_spikes) + + + for i in df_map.keys(): + freq = list(df_map[i]) + k = freq[0] + phase = list(chirp_spikes[k].keys())[0] + + spikes = chirp_spikes[k][phase] + 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 = 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_eodform_spikehist.py b/code/plot_eodform_spikehist.py index 9512ca6..e7a5cdd 100644 --- a/code/plot_eodform_spikehist.py +++ b/code/plot_eodform_spikehist.py @@ -14,6 +14,9 @@ dataset = '2018-11-14-ad-invivo-1' # read eod and time of baseline time, eod = read_baseline_eod(os.path.join(data_dir, dataset)) + + + eod_norm = eod - np.mean(eod) # calculate eod times and indices by zero crossings @@ -23,6 +26,10 @@ eod_times = time[(eod_norm >= threshold) & (shift_eod < threshold)] eod_duration = eod_times[2]- eod_times[1] #time in s +eod_duration = eod_times[2]- eod_times[1] + + + # read spikes during baseline activity spikes = read_baseline_spikes(os.path.join(data_dir, dataset)) #spikes in s # calculate interpike intervals and plot them @@ -37,9 +44,8 @@ plt.yticks(fontsize = 18) ax.spines["top"].set_visible(False) ax.spines["right"].set_visible(False) fig.tight_layout() -plt.show() +#plt.show() #plt.savefig('isis.pdf') -exit() plt.savefig('isis.png') @@ -93,10 +99,10 @@ plt.yticks(fontsize=18) ax1.spines['top'].set_visible(False) ax2 = ax1.twinx() -ax2.fill_between(time_axis, mu_eod+std_eod, mu_eod-std_eod, color='navy', alpha=0.5) +ax2.fill_between(time_axis, mu_eod+std_eod, mu_eod-std_eod, color='royalblue', alpha=0.5) ax2.plot(time_axis, mu_eod, color='black', lw=2) ax2.set_ylabel('voltage [mV]', fontsize=22) -ax2.tick_params(axis='y', labelcolor='navy') +ax2.tick_params(axis='y', labelcolor='royalblue') ax2.spines['top'].set_visible(False) plt.yticks(fontsize=18) 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) diff --git a/code/stimulus_chirp.py b/code/stimulus_chirp.py index f45827c..0e01894 100644 --- a/code/stimulus_chirp.py +++ b/code/stimulus_chirp.py @@ -35,13 +35,13 @@ for k, t in enumerate(time): p += f * stepsize signal[k] = a * np.sin(6.28318530717959 * p) -fig = plt.figure(figsize = (20/inch_factor, 15/inch_factor)) +fig = plt.figure(figsize = (20/inch_factor, 12/inch_factor)) ax1 = fig.add_subplot(211) plt.yticks(fontsize=18) ax2 = fig.add_subplot(212, sharex=ax1) plt.setp(ax1.get_xticklabels(), visible=False) -ax1.plot(time*1000, signal, color = 'midnightblue', lw = 1) -ax2.plot(time*1000, freq, color = 'midnightblue', lw = 3) +ax1.plot(time*1000, signal, color = 'royalblue', lw = 1) +ax2.plot(time*1000, freq, color = 'royalblue', lw = 3) ax1.set_ylabel("field [mV]", fontsize = 22) @@ -53,5 +53,5 @@ ax2.yaxis.set_label_coords(-0.1, 0.5) plt.xticks(fontsize=18) plt.yticks(fontsize=18) fig.tight_layout() -#plt.show() -plt.savefig('stimulus_chirp.png') +plt.show() +#plt.savefig('stimulus_chirp.png') diff --git a/code/vector_phase.py b/code/vector_phase.py deleted file mode 100644 index c902e6e..0000000 --- a/code/vector_phase.py +++ /dev/null @@ -1,25 +0,0 @@ -from read_baseline_data import * -from utility import * -#import nix_helpers as nh -import matplotlib.pyplot as plt -import numpy as np -from IPython import embed #Funktionen importieren - - -#Zeitpunkte einer EOD über Zero-crossings finden, die in einer Steigung liegen -data_dir = "../data" -dataset = "2018-11-09-ad-invivo-1" -time,eod = read_baseline_eod(os.path.join(data_dir, dataset)) -spike_times = read_baseline_spikes(os.path.join(data_dir, dataset)) -print(len(spike_times)) - -eod_times = zero_crossing(eod,time) -eod_durations = np.diff(eod_times) -print(len(spike_times)) -print(len(eod_durations)) - -#for st in spike_times: - #et = eod_times[eod_times < st] - #dt = st - et - -#vs = vector_strength(spike_times, eod_durations)