From 454ac6984921d4f10a93b047add8035f7fcede2f Mon Sep 17 00:00:00 2001 From: Ramona Date: Mon, 19 Nov 2018 16:06:57 +0100 Subject: [PATCH 01/10] phases --- code/spikes_analysis.py | 45 +++++++++++++++++++++++++++-------------- 1 file changed, 30 insertions(+), 15 deletions(-) diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index 5ae3da4..e059c38 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -9,7 +9,7 @@ dataset = "2018-11-09-ad-invivo-1" spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) eod = read_chirp_eod(os.path.join(data_dir, dataset)) -times = read_chirp_times(os.path.join(data_dir, dataset)) +chirp_times = read_chirp_times(os.path.join(data_dir, dataset)) df_map = {} for k in spikes.keys(): @@ -20,31 +20,46 @@ for k in spikes.keys(): df_map[df] = [k] # make phases together, 12 phases -spikes_mat = {} +phase_vec = np.arange(0, 1+1/12, 1/12) +phase_mat_df = {} + for deltaf in df_map.keys(): + phase_df = {} for rep in df_map[deltaf]: for phase in spikes[rep]: #print(phase) - spikes_one_chirp = spikes[rep][phase] - if deltaf == '-50Hz' and phase == (9, 0.54): - spikes_mat[deltaf, rep, phase] = spikes_one_chirp + for idx in range(len(phase_vec)-1): + if phase[1] > phase_vec[idx] and phase[1] < phase_vec[idx+1]: + if phase_vec[idx] in phase_df.keys(): + phase_df[phase_vec[idx]].append(phase[1]) + else: + phase_df[phase_vec[idx]] = [phase[1]] -plot_spikes = spikes[(0, '-50Hz', '20%', '100Hz')][(0, 0.789)] + #spikes_one_chirp = spikes[rep][phase] -mu = 1 -sigma = 1 -time_gauss = np.arange(-4, 4, 1) -gauss = gaussian(time_gauss, mu, sigma) -# spikes during time vec (00010000001)? -smoothed_spikes = np.convolve(plot_spikes, gauss, 'same') -window = np.mean(np.diff(plot_spikes)) -time_vec = np.arange(plot_spikes[0], plot_spikes[-1]+window, window) + + phase_mat_df[deltaf] = phase_df + +embed() +exit() +plot_spikes = spikes[(0, '-50Hz', '20%', '100Hz')][(0, 0.789)] fig, ax = plt.subplots() ax.scatter(plot_spikes, np.ones(len(plot_spikes))*10, marker='|', color='k') -ax.plot(time_vec, smoothed_spikes) plt.show() +#mu = 1 +#sigma = 1 +#time_gauss = np.arange(-4, 4, 1) +#gauss = gaussian(time_gauss, mu, sigma) +# spikes during time vec (00010000001)? + +#smoothed_spikes = np.convolve(plot_spikes, gauss, 'same') +#window = np.mean(np.diff(plot_spikes)) +#time_vec = np.arange(plot_spikes[0], plot_spikes[-1]+window, window) + +#ax.plot(time_vec, smoothed_spikes) + #embed() #exit() #hist_data = plt.hist(plot_spikes, bins=np.arange(-200, 400, 20)) From ae6e78f9a749481a2ba56035c0bfe17b900549c4 Mon Sep 17 00:00:00 2001 From: Ramona Date: Mon, 19 Nov 2018 17:12:09 +0100 Subject: [PATCH 02/10] still problems with access --- code/spikes_analysis.py | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index e059c38..2d72447 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -31,21 +31,25 @@ for deltaf in df_map.keys(): for idx in range(len(phase_vec)-1): if phase[1] > phase_vec[idx] and phase[1] < phase_vec[idx+1]: if phase_vec[idx] in phase_df.keys(): - phase_df[phase_vec[idx]].append(phase[1]) + phase_df[phase_vec[idx]].append(phase) else: - phase_df[phase_vec[idx]] = [phase[1]] + phase_df[phase_vec[idx]] = [phase] #spikes_one_chirp = spikes[rep][phase] phase_mat_df[deltaf] = phase_df -embed() -exit() -plot_spikes = spikes[(0, '-50Hz', '20%', '100Hz')][(0, 0.789)] +#df_spikes = spikes[df_map['-50Hz'][0]] +#phase_spikes = df_spikes[phase_mat_df['-50Hz'][0.0][0]] +trial_spikes = np.asarray(spikes[(0, '-50Hz', '20%', '100Hz')][(0, 0.789)]) +plot_spikes = trial_spikes[(trial_spikes < 50.0) & (trial_spikes > -50.0)] + +hist_data = plt.hist(plot_spikes, color='w') fig, ax = plt.subplots() ax.scatter(plot_spikes, np.ones(len(plot_spikes))*10, marker='|', color='k') +ax.plot(hist_data[1][:-1], hist_data[0]) plt.show() #mu = 1 From 1daf5c244ed736f8113df00a20e8be9550393060 Mon Sep 17 00:00:00 2001 From: Ramona Date: Tue, 20 Nov 2018 17:40:37 +0100 Subject: [PATCH 03/10] data is ready for plotting --- code/spikes_analysis.py | 42 +++++++++++++++++++++++++---------------- 1 file changed, 26 insertions(+), 16 deletions(-) diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index 2d72447..a665a80 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -4,6 +4,7 @@ from read_chirp_data import * from utility import * from IPython import embed +sampling_rate = 40 #kHz data_dir = "../data" dataset = "2018-11-09-ad-invivo-1" @@ -21,37 +22,46 @@ for k in spikes.keys(): # make phases together, 12 phases phase_vec = np.arange(0, 1+1/12, 1/12) -phase_mat_df = {} +cut_range = np.arange(-50*sampling_rate, 50*sampling_rate, 1) + +df_phase_time = {} +df_phase_binary = {} for deltaf in df_map.keys(): - phase_df = {} + df_phase_time[deltaf] = {} + df_phase_binary[deltaf] = {} for rep in df_map[deltaf]: for phase in spikes[rep]: #print(phase) for idx in range(len(phase_vec)-1): if phase[1] > phase_vec[idx] and phase[1] < phase_vec[idx+1]: - if phase_vec[idx] in phase_df.keys(): - phase_df[phase_vec[idx]].append(phase) - else: - phase_df[phase_vec[idx]] = [phase] - - #spikes_one_chirp = spikes[rep][phase] + spikes_to_cut = np.asarray(spikes[rep][phase]) + spikes_cut = spikes_to_cut[(spikes_to_cut > -50) & (spikes_to_cut < 50)] + spikes_idx = np.round(spikes_cut*sampling_rate) + binary_spikes = np.isin(cut_range, spikes_idx)*1 - phase_mat_df[deltaf] = phase_df + if phase_vec[idx] in df_phase_time[deltaf].keys(): + df_phase_time[deltaf][phase_vec[idx]].append(spikes[rep][phase]) + df_phase_binary[deltaf][phase_vec[idx]] = np.vstack((df_phase_binary[deltaf][phase_vec[idx]], binary_spikes)) + else: + df_phase_time[deltaf][phase_vec[idx]] = [spikes[rep][phase]] + df_phase_binary[deltaf][phase_vec[idx]] = binary_spikes -#df_spikes = spikes[df_map['-50Hz'][0]] -#phase_spikes = df_spikes[phase_mat_df['-50Hz'][0.0][0]] -trial_spikes = np.asarray(spikes[(0, '-50Hz', '20%', '100Hz')][(0, 0.789)]) -plot_spikes = trial_spikes[(trial_spikes < 50.0) & (trial_spikes > -50.0)] +plot_trials = df_phase_binary['-50Hz'][0.0] +#hist_data = plt.hist(plot_trials) +#ax.plot(hist_data[1][:-1], hist_data[0]) -hist_data = plt.hist(plot_spikes, color='w') fig, ax = plt.subplots() -ax.scatter(plot_spikes, np.ones(len(plot_spikes))*10, marker='|', color='k') -ax.plot(hist_data[1][:-1], hist_data[0]) +for i, trial in enumerate(plot_trials): + embed() + exit() + trial[trial == 0] = np.nan + ax.scatter(np.ones(len(trial)), trial, marker='|', color='k', size=12) plt.show() + #mu = 1 #sigma = 1 #time_gauss = np.arange(-4, 4, 1) From 9d971b44e5d9dac044ac7780cfc112346f67af52 Mon Sep 17 00:00:00 2001 From: efish Date: Wed, 21 Nov 2018 10:01:05 +0100 Subject: [PATCH 04/10] Mo --- code/base_chirps.py | 48 ++++++++++++++++++++++++++++++++++----------- code/base_spikes.py | 40 ++++++++++++++++++++++++++++++++----- 2 files changed, 72 insertions(+), 16 deletions(-) diff --git a/code/base_chirps.py b/code/base_chirps.py index dfdfe37..80c14bb 100644 --- a/code/base_chirps.py +++ b/code/base_chirps.py @@ -1,5 +1,5 @@ from read_chirp_data import * -#import nix_helpers as nh +import nix_helpers as nh import matplotlib.pyplot as plt import numpy as np from IPython import embed @@ -27,14 +27,14 @@ for k in eod.keys(): else: df_map[df] = [k] -print(ch) #die Chirphöhe wird ausgegeben, um zu bestimmen, ob Chirps oder Chirps large benutzt wurde +#print(ch) #die Chirphöhe wird ausgegeben, um zu bestimmen, ob Chirps oder Chirps large benutzt wurde #die äußere Schleife geht für alle Keys durch und somit durch alle dfs #die innnere Schleife bildet die 16 Wiederholungen einer Frequenz in 4 Subplots ab -for idx in df_map.keys(): - freq = list(df_map[idx]) +for i in df_map.keys(): + freq = list(df_map[i]) fig,axs = plt.subplots(2, 2, sharex = True, sharey = True) for idx, k in enumerate(freq): @@ -58,18 +58,44 @@ for idx in df_map.keys(): fig.suptitle('EOD for chirps', fontsize = 16) -plt.show() +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]') +#axs kann nur einzelne Label erzeugen, nicht generell möglich wie beim Titel -#Problem: axs hat keine label-Funktion, also müsste axes nochmal definiert werden. Momentan erscheint Schrift nur auf einem der Subplots +for i in df_map.keys(): + freq = list(df_map[i]) -#ax = plt.gca() -#ax.set_ylabel('Time [ms]') -#ax.set_xlabel('Amplitude [mV]') -#ax.label_outer() + ct = times[freq[1]] + ct1 = ct[1] + e1 = eod[k] + zeit = np.asarray(e1[0]) + ampl = np.asarray(e1[1]) + + time_cut = zeit[(zeit > ct1-25) & (zeit < ct1+25)] + eod_cut = ampl[(zeit > ct1-25) & (zeit < ct1+25)] + + change = ampl[int(ct1)] + + plt.figure(12) + plt.plot(time_cut, eod_cut) + plt.scatter(ct1, 3, color = 'green', s= 30) + plt.title('Chirp reaction Ampl.') + plt.xlabel('Time [ms]') + plt.ylabel('Amplitude[mV]') + #plt.show() + + +#4. Chirps einer Phase zuordnen - zusammen plotten? + + + +#next Step: EOD-Amplitudenmodulation für beat aber OHNE Chirps plotten +#allerdings: in der Aufnahme sind nur kurze Zeitfenster ohne Chirps zu finden! -#next Step: relative Amplitudenmodulation berechnen, Max und Min der Amplitude bestimmen, EOD und Chirps zuordnen, Unterschied berechnen diff --git a/code/base_spikes.py b/code/base_spikes.py index 7ab065f..8e6793b 100644 --- a/code/base_spikes.py +++ b/code/base_spikes.py @@ -1,16 +1,18 @@ from read_baseline_data import * -#import nix_helpers as nh +from read_chirp_data import * +import nix_helpers as nh import matplotlib.pyplot as plt import numpy as np from IPython import embed #Funktionen importieren + data_dir = "../data" -dataset = "2018-11-09-aa-invivo-1" +dataset = "2018-11-09-ad-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") -spike_times = read_baseline_spikes(os.path.join(data_dir, dataset)) -#spike_frequency = len(spike_times) / spike_times[-1] + +spike_times = read_baseline_spikes(os.path.join(data_dir, dataset)) #inst_frequency = 1. / np.diff(spike_times) spike_rate = np.diff(spike_times) @@ -21,7 +23,6 @@ plt.hist(spike_rate,x) mu = np.mean(spike_rate) sigma = np.std(spike_rate) cv = sigma/mu -print(cv) plt.title('A.lepto ISI Histogramm', fontsize = 14) plt.xlabel('duration ISI[ms]', fontsize = 12) @@ -32,3 +33,32 @@ plt.yticks(fontsize = 12) plt.show() + + +#Nyquist-Theorem Plot: + +chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) + + +df_map = {} +#Keys werden nach df sortiert ausgegeben +for k in chirp_spikes.keys(): + df = k[1] + ch = k[3] + if df in df_map.keys(): + df_map[df].append(k) + else: + df_map[df] = [k] + + +for i in df_map.keys(): + freq = list(df_map[i]) + for k in freq: + + spikes = chirp_spikes[k] + trial = spikes[1] + print(trial) +# +# plt.plot(spikes, rate) +# plt.show() + From a151ebd886704ddde4f278fcd6ee01eaf3269f8f Mon Sep 17 00:00:00 2001 From: efish Date: Wed, 21 Nov 2018 10:10:19 +0100 Subject: [PATCH 05/10] nFkt --- code/utility.py | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/code/utility.py b/code/utility.py index 68102d3..0f3ce19 100644 --- a/code/utility.py +++ b/code/utility.py @@ -23,3 +23,14 @@ def gaussian(x, mu, sig): y = np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.))) return y +def map_keys(input): + df_map = {} + for k in input.keys(): + df = k[1] + ch = k[3] + if df in df_map.keys(): + df_map[df].append(k) + else: + df_map[df] = [k] + return df_map + print(ch) From c713cbc1ffd4e9c890af2e505c0fdebbf533a8b4 Mon Sep 17 00:00:00 2001 From: Ramona Date: Wed, 21 Nov 2018 11:01:45 +0100 Subject: [PATCH 06/10] finally :D --- code/spikes_analysis.py | 30 ++++++++++++++---------------- 1 file changed, 14 insertions(+), 16 deletions(-) diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index a665a80..41451a4 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -42,33 +42,31 @@ for deltaf in df_map.keys(): binary_spikes = np.isin(cut_range, spikes_idx)*1 if phase_vec[idx] in df_phase_time[deltaf].keys(): - df_phase_time[deltaf][phase_vec[idx]].append(spikes[rep][phase]) + df_phase_time[deltaf][phase_vec[idx]].append(spikes_cut) df_phase_binary[deltaf][phase_vec[idx]] = np.vstack((df_phase_binary[deltaf][phase_vec[idx]], binary_spikes)) else: - df_phase_time[deltaf][phase_vec[idx]] = [spikes[rep][phase]] + df_phase_time[deltaf][phase_vec[idx]] = [spikes_cut] df_phase_binary[deltaf][phase_vec[idx]] = binary_spikes -plot_trials = df_phase_binary['-50Hz'][0.0] -#hist_data = plt.hist(plot_trials) -#ax.plot(hist_data[1][:-1], hist_data[0]) + +plot_trials = df_phase_time['-50Hz'][0.0] +plot_trials_binary = np.mean(df_phase_binary['-50Hz'][0.0], axis=0) + +mu = 1 +sigma = 100 +time_gauss = np.arange(-4*sigma, 4*sigma, 1) +gauss = gaussian(time_gauss, mu, sigma) +smoothed_spikes = np.convolve(plot_trials_binary, gauss, 'same') +time_axis = np.arange(-50, 50, 1/sampling_rate) fig, ax = plt.subplots() for i, trial in enumerate(plot_trials): - embed() - exit() - trial[trial == 0] = np.nan - ax.scatter(np.ones(len(trial)), trial, marker='|', color='k', size=12) + ax.scatter(trial, np.ones(len(trial))+i, marker='|', color='k') +ax.plot(time_axis, smoothed_spikes) plt.show() -#mu = 1 -#sigma = 1 -#time_gauss = np.arange(-4, 4, 1) -#gauss = gaussian(time_gauss, mu, sigma) -# spikes during time vec (00010000001)? - -#smoothed_spikes = np.convolve(plot_spikes, gauss, 'same') #window = np.mean(np.diff(plot_spikes)) #time_vec = np.arange(plot_spikes[0], plot_spikes[-1]+window, window) From 15e7fbc51fdae0e15e0f97475eac56eec5dd30b0 Mon Sep 17 00:00:00 2001 From: Ramona Date: Wed, 21 Nov 2018 11:06:00 +0100 Subject: [PATCH 07/10] smoothing function --- code/spikes_analysis.py | 7 ++----- code/utility.py | 10 ++++++++++ 2 files changed, 12 insertions(+), 5 deletions(-) diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index 41451a4..e369e8e 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -53,11 +53,8 @@ for deltaf in df_map.keys(): plot_trials = df_phase_time['-50Hz'][0.0] plot_trials_binary = np.mean(df_phase_binary['-50Hz'][0.0], axis=0) -mu = 1 -sigma = 100 -time_gauss = np.arange(-4*sigma, 4*sigma, 1) -gauss = gaussian(time_gauss, mu, sigma) -smoothed_spikes = np.convolve(plot_trials_binary, gauss, 'same') +window = 100 +smoothed_spikes = smooth(plot_trials_binary, window) time_axis = np.arange(-50, 50, 1/sampling_rate) fig, ax = plt.subplots() diff --git a/code/utility.py b/code/utility.py index 0f3ce19..3bfdb30 100644 --- a/code/utility.py +++ b/code/utility.py @@ -23,6 +23,16 @@ def gaussian(x, mu, sig): y = np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.))) return y + +def smooth(data, window): + mu = 1 + sigma = window + time_gauss = np.arange(-4 * sigma, 4 * sigma, 1) + gauss = gaussian(time_gauss, mu, sigma) + smoothed_data = np.convolve(data, gauss, 'same') + return smoothed_data + + def map_keys(input): df_map = {} for k in input.keys(): From a7aec1b12f33148cc1d98c1c120093eea06bdefd Mon Sep 17 00:00:00 2001 From: Jan Grewe Date: Wed, 21 Nov 2018 12:45:00 +0100 Subject: [PATCH 08/10] [chirp data] read the Stimulus from file --- code/read_chirp_data.py | 39 ++++++++++++++++++++++++++++++++++----- 1 file changed, 34 insertions(+), 5 deletions(-) diff --git a/code/read_chirp_data.py b/code/read_chirp_data.py index c3afe76..4491779 100644 --- a/code/read_chirp_data.py +++ b/code/read_chirp_data.py @@ -1,5 +1,7 @@ import numpy as np import os +import nixio as nix +from IPython import embed def read_chirp_spikes(dataset): @@ -85,10 +87,37 @@ def read_chirp_times(dataset): return chirp_times +def read_chirp_stimulus(dataset): + base = dataset.split(os.path.sep)[-1] + ".nix" + nix_file = nix.File.open(os.path.join(dataset, base), nix.FileMode.ReadOnly) + b = nix_file.blocks[0] + data = {} + for t in b.tags: + if "Chirps" in t.name: + stims = [] + index = int(t.name.split("_")[-1]) + df = t.metadata["RePro-Info"]["settings"]["deltaf"] + cs = t.metadata["RePro-Info"]["settings"]["chirpsize"] + stim_da = t.references["GlobalEFieldStimulus"] + si = stim_da.dimensions[0].sampling_interval + for mt in b.multi_tags: + if mt.positions[0] >= t.position[0] and \ + mt.positions[0] < (t.position[0] + t.extent[0]): + break + for i in range(len(mt.positions)): + start_index = int(mt.positions[i] / si) + end_index = int((mt.positions[i] + mt.extents[i]) / si) - 1 + stim = stim_da[start_index:end_index] + time = stim_da.dimensions[0].axis(len(stim)) + mt.positions[i] + stims.append((time, stim)) + data[(index, df, cs)] = stims + nix_file.close() + return data + if __name__ == "__main__": data_dir = "../data" - dataset = "2018-11-09-ad-invivo-1" - spikes = load_chirp_spikes(os.path.join(data_dir, dataset)) - chirp_times = load_chirp_times(os.path.join(data_dir, dataset)) - chirp_eod = load_chirp_eod(os.path.join(data_dir, dataset)) - + dataset = "2018-11-20-ad-invivo-1" + #spikes = load_chirp_spikes(os.path.join(data_dir, dataset)) + #chirp_times = load_chirp_times(os.path.join(data_dir, dataset)) + #chirp_eod = load_chirp_eod(os.path.join(data_dir, dataset)) + stim = read_chirp_stimulus(os.path.join(data_dir, dataset)) From 54a78ec29687def199adf67dee6ef62945d2d551 Mon Sep 17 00:00:00 2001 From: Ramona Date: Wed, 21 Nov 2018 13:49:22 +0100 Subject: [PATCH 09/10] still mistakes when normalizing --- code/spikes_analysis.py | 66 ++++++++++++++++++++++------------------- code/utility.py | 4 ++- 2 files changed, 39 insertions(+), 31 deletions(-) diff --git a/code/spikes_analysis.py b/code/spikes_analysis.py index e369e8e..2d4a329 100644 --- a/code/spikes_analysis.py +++ b/code/spikes_analysis.py @@ -4,14 +4,21 @@ 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-09-ad-invivo-1" +# parameters for binning, smoothing and plotting +num_bin = 12 +window = sampling_rate +time_axis = np.arange(-50, 50, 1/sampling_rate) +# read data from files spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) eod = read_chirp_eod(os.path.join(data_dir, dataset)) chirp_times = read_chirp_times(os.path.join(data_dir, dataset)) +# make a delta f map for the quite more complicated keys df_map = {} for k in spikes.keys(): df = k[1] @@ -20,56 +27,55 @@ for k in spikes.keys(): else: df_map[df] = [k] -# make phases together, 12 phases -phase_vec = np.arange(0, 1+1/12, 1/12) +# differentiate between phases +phase_vec = np.arange(0, 1+1/num_bin, 1/num_bin) cut_range = np.arange(-50*sampling_rate, 50*sampling_rate, 1) +# make dictionaries for spiketimes df_phase_time = {} df_phase_binary = {} +# iterate over delta f, repetition, phases and a single chirp for deltaf in df_map.keys(): df_phase_time[deltaf] = {} df_phase_binary[deltaf] = {} for rep in df_map[deltaf]: for phase in spikes[rep]: - #print(phase) - for idx in range(len(phase_vec)-1): + for idx in np.arange(num_bin): + # check the phase if phase[1] > phase_vec[idx] and phase[1] < phase_vec[idx+1]: + # get spikes between 50 ms befor and after the chirp spikes_to_cut = np.asarray(spikes[rep][phase]) spikes_cut = spikes_to_cut[(spikes_to_cut > -50) & (spikes_to_cut < 50)] 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 - if phase_vec[idx] in df_phase_time[deltaf].keys(): - df_phase_time[deltaf][phase_vec[idx]].append(spikes_cut) - df_phase_binary[deltaf][phase_vec[idx]] = np.vstack((df_phase_binary[deltaf][phase_vec[idx]], binary_spikes)) + # add the spikes to the dictionaries with the correct df and phase + if idx in df_phase_time[deltaf].keys(): + df_phase_time[deltaf][idx].append(spikes_cut) + df_phase_binary[deltaf][idx] = np.vstack((df_phase_binary[deltaf][idx], binary_spikes)) else: - df_phase_time[deltaf][phase_vec[idx]] = [spikes_cut] - df_phase_binary[deltaf][phase_vec[idx]] = binary_spikes - - - -plot_trials = df_phase_time['-50Hz'][0.0] -plot_trials_binary = np.mean(df_phase_binary['-50Hz'][0.0], axis=0) - -window = 100 -smoothed_spikes = smooth(plot_trials_binary, window) -time_axis = np.arange(-50, 50, 1/sampling_rate) + df_phase_time[deltaf][idx] = [spikes_cut] + df_phase_binary[deltaf][idx] = binary_spikes -fig, ax = plt.subplots() -for i, trial in enumerate(plot_trials): - ax.scatter(trial, np.ones(len(trial))+i, marker='|', color='k') -ax.plot(time_axis, smoothed_spikes) -plt.show() +# for plotting iterate over delta f and phases +for df in df_phase_time.keys(): + for phase in df_phase_time[df].keys(): + plot_trials = df_phase_time[df][phase] + plot_trials_binary = np.mean(df_phase_binary[df][phase], axis=0) + smoothed_spikes = smooth(plot_trials_binary, window) -#window = np.mean(np.diff(plot_spikes)) -#time_vec = np.arange(plot_spikes[0], plot_spikes[-1]+window, window) + fig, ax = plt.subplots(2, 1) + for i, trial in enumerate(plot_trials): + ax[0].scatter(trial, np.ones(len(trial))+i, marker='|', color='k') + ax[1].plot(time_axis, smoothed_spikes) -#ax.plot(time_vec, smoothed_spikes) + ax[0].set_title(df) + ax[0].set_ylabel('repetition', fontsize=12) -#embed() -#exit() -#hist_data = plt.hist(plot_spikes, bins=np.arange(-200, 400, 20)) -#ax.plot(hist_data[1][:-1], hist_data[0]) \ No newline at end of file + ax[1].set_xlabel('time [ms]', fontsize=12) + ax[1].set_ylabel('firing rate [?]', fontsize=12) + plt.show() diff --git a/code/utility.py b/code/utility.py index 3bfdb30..41f75c9 100644 --- a/code/utility.py +++ b/code/utility.py @@ -1,4 +1,5 @@ import numpy as np +from IPython import embed def zero_crossing(eod, time): @@ -29,7 +30,8 @@ def smooth(data, window): sigma = window time_gauss = np.arange(-4 * sigma, 4 * sigma, 1) gauss = gaussian(time_gauss, mu, sigma) - smoothed_data = np.convolve(data, gauss, 'same') + gauss_norm = gauss/(np.sum(gauss)/len(gauss)) + smoothed_data = np.convolve(data, gauss_norm, 'same') return smoothed_data From 470240df1d362f5112b808cff31d3043c7a5a37b Mon Sep 17 00:00:00 2001 From: efish Date: Wed, 21 Nov 2018 13:53:52 +0100 Subject: [PATCH 10/10] oops --- code/base_chirps.py | 60 ++++++++++++++++----------------------------- code/base_spikes.py | 23 ++++++----------- code/utility.py | 4 +-- 3 files changed, 31 insertions(+), 56 deletions(-) diff --git a/code/base_chirps.py b/code/base_chirps.py index 80c14bb..1112d78 100644 --- a/code/base_chirps.py +++ b/code/base_chirps.py @@ -1,5 +1,6 @@ from read_chirp_data import * -import nix_helpers as nh +from utility import * +#import nix_helpers as nh import matplotlib.pyplot as plt import numpy as np from IPython import embed @@ -15,24 +16,12 @@ data = ("2018-11-09-ad-invivo-1", "2018-11-09-ae-invivo-1", "2018-11-09-ag-inviv #for dataset in data: eod = read_chirp_eod(os.path.join(data_dir, dataset)) times = read_chirp_times(os.path.join(data_dir, dataset)) - - - -df_map = {} #Keys werden nach df sortiert ausgegeben -for k in eod.keys(): - df = k[1] - ch = k[3] - if df in df_map.keys(): - df_map[df].append(k) - else: - df_map[df] = [k] - -#print(ch) #die Chirphöhe wird ausgegeben, um zu bestimmen, ob Chirps oder Chirps large benutzt wurde +df_map = map_keys(eod) #die äußere Schleife geht für alle Keys durch und somit durch alle dfs -#die innnere Schleife bildet die 16 Wiederholungen einer Frequenz in 4 Subplots ab +#die innnere Schleife bildet die 16 Wiederholungen einer Frequenz ab for i in df_map.keys(): freq = list(df_map[i]) fig,axs = plt.subplots(2, 2, sharex = True, sharey = True) @@ -62,40 +51,33 @@ 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]') -#axs kann nur einzelne Label erzeugen, nicht generell möglich wie beim Titel - -for i in df_map.keys(): - freq = list(df_map[i]) - ct = times[freq[1]] - ct1 = ct[1] +#for i in df_map.keys(): +freq = list(df_map['-50Hz']) +ls_mod = [] +beat_mods = [] +for k in freq: e1 = eod[k] zeit = np.asarray(e1[0]) ampl = np.asarray(e1[1]) - time_cut = zeit[(zeit > ct1-25) & (zeit < ct1+25)] - eod_cut = ampl[(zeit > ct1-25) & (zeit < ct1+25)] - - change = ampl[int(ct1)] - - plt.figure(12) - plt.plot(time_cut, eod_cut) - plt.scatter(ct1, 3, color = 'green', s= 30) - plt.title('Chirp reaction Ampl.') - plt.xlabel('Time [ms]') - plt.ylabel('Amplitude[mV]') - #plt.show() - - -#4. Chirps einer Phase zuordnen - zusammen plotten? - + ct = times[k] + for chirp in ct: + time_cut = zeit[(zeit > chirp-10) & (zeit < chirp+10)] + eods_cut = ampl[(zeit > chirp-10) & (zeit < chirp+10)] + beat_cut = ampl[(zeit > chirp-55) & (zeit < chirp-10)] + chirp_mod = np.std(eods_cut) #Std vom Bereich um den Chirp + beat_mod = np.std(beat_cut) #Std vom Bereich vor dem Chirp + ls_mod.append(chirp_mod) + beat_mods.append(beat_mod) -#next Step: EOD-Amplitudenmodulation für beat aber OHNE Chirps plotten -#allerdings: in der Aufnahme sind nur kurze Zeitfenster ohne Chirps zu finden! +#Länge des Mods ist 160, 16 Wiederholungen mal 10 Chirps pro Trial +#Verwendung der Std für die Amplitudenmodulation? +#Chirps einer Phase zuordnen - zusammen plotten? diff --git a/code/base_spikes.py b/code/base_spikes.py index 8e6793b..2fd0e9f 100644 --- a/code/base_spikes.py +++ b/code/base_spikes.py @@ -1,6 +1,7 @@ from read_baseline_data import * from read_chirp_data import * -import nix_helpers as nh +from utility import * +#import nix_helpers as nh import matplotlib.pyplot as plt import numpy as np from IPython import embed #Funktionen importieren @@ -38,26 +39,18 @@ plt.show() #Nyquist-Theorem Plot: chirp_spikes = read_chirp_spikes(os.path.join(data_dir, dataset)) - - -df_map = {} -#Keys werden nach df sortiert ausgegeben -for k in chirp_spikes.keys(): - df = k[1] - ch = k[3] - if df in df_map.keys(): - df_map[df].append(k) - else: - df_map[df] = [k] +df_map = map_keys(chirp_spikes) for i in df_map.keys(): freq = list(df_map[i]) for k in freq: - spikes = chirp_spikes[k] - trial = spikes[1] - print(trial) + phase_map = map_keys(spikes) + for p in phase_map: + spike_rate = 1./ np.diff(p) + +print(spike_rate) # # plt.plot(spikes, rate) # plt.show() diff --git a/code/utility.py b/code/utility.py index 0f3ce19..2a3f083 100644 --- a/code/utility.py +++ b/code/utility.py @@ -27,10 +27,10 @@ def map_keys(input): df_map = {} for k in input.keys(): df = k[1] - ch = k[3] + #ch = k[3] if df in df_map.keys(): df_map[df].append(k) else: df_map[df] = [k] return df_map - print(ch) + #print(ch)