benda/susceptibility1
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

neues bild hinzugefügt

Browse Source
This commit is contained in:
saschuta 2024-06-13 12:10:14 +02:00
parent 08388ae10e
commit 75106e9095
6 changed files with 80 additions and 93 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

BIN
model_full.pdf
View File

Binary file not shown.

BIN
model_full.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 178 KiB

After

Width:  |  Height:  |  Size: 178 KiB

BIN
nonlin_regime.pdf
View File

Binary file not shown.

BIN
nonlin_regime.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 114 KiB

After

Width:  |  Height:  |  Size: 56 KiB

150
nonlin_regime.py
View File

@ -40,7 +40,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
#reshuffled = 'reshuffled' # ,
# standard combination with intruder small
a_f2s = [0.1]
min_amps = '_minamps_'
dev_name = ['05']
@ -60,7 +60,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
position_diff = 0
plot_style()
default_figsize(column=2, length=7.5)
default_figsize(column=2, length=5.5)
########################################
@ -75,7 +75,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
DF2_frmult[d] = recalc_fr_to_DF1(DF2_frmult, d, frame_cell.fr_data.iloc[0])
DF1_frmult[d] = recalc_fr_to_DF1(DF1_frmult, d, frame_cell.fr_data.iloc[0])
##DF2_frmult[d] = recalc_fr_to_DF1(DF2_frmult, d, frame_cell.fr_data.iloc[0])
freqs = [(DF1_frmult[2], DF2_frmult[2])]
freqs = [(DF1_frmult[3], DF2_frmult[3])]
@ -88,6 +88,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
nfft = 32768
cells_here = ['2012-07-03-ak-invivo-1']
cells_here = ["2013-01-08-aa-invivo-1"]
for cell_here in cells_here:
###########################################
@ -100,24 +101,26 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
full_names = ['calc_model_amp_freqs-_F1_0.22833333333333333Fr_F2_1Fr_af_coupled__C1Len_50_FirstC1_0.0001_LastC1_1.0_FRrelavtiv__start_0.0001_end_1_StimLen_100_nfft_20000_trialsnr_1_mult_minimum_1_power_1_minamps__dev_original_05_point_1temporal',
'calc_model_amp_freqs-_F1_0.22833333333333333Fr_F2_1Fr_af_coupled__C1Len_50_FirstC1_0.0001_LastC1_1.0_FRrelavtiv__start_0.0001_end_1_StimLen_100_nfft_20000_trialsnr_1_mult_minimum_1_power_1_minamps__dev_original_05_point_1_old_fit_temporal']
full_names = ['calc_model_amp_freqs-_F1_0.22833333333333333Fr_F2_1Fr_af_coupled__C1Len_50_FirstC1_0.0001_LastC1_1.0_FRrelavtiv__start_0.0001_end_1_StimLen_100_nfft_20000_trialsnr_1_mult_minimum_1_power_1_minamps__dev_original_05_point_1_fft2_temporal']
c_grouped = ['c1'] # , 'c2']
c_nrs_orig = [0.02, 0.2] # 0.0002, 0.05, 0.5
c_nrs_orig = [0.01,0.025, 0.04, 0.1] # 0.0002, 0.05, 0.5
trials_nr = 20 # 20
redo = False # True
log = 'log' # 'log'
grid0 = gridspec.GridSpec(1, 1, bottom=0.08, top=0.93, left=0.11,
grid0 = gridspec.GridSpec(1, 1, bottom=0.15, top=0.88, left=0.11,
right=0.95, wspace=0.04) #
grid00 = gridspec.GridSpecFromSubplotSpec(2, 1,
wspace=0.04, hspace=0.1,
wspace=0.04, hspace=0.45,
subplot_spec=grid0[0], height_ratios=[1.5,1],) # height_ratios=[2,1],
grid_up = gridspec.GridSpecFromSubplotSpec(len(c_nrs_orig) + 1, 3,
grid_up = gridspec.GridSpecFromSubplotSpec(1, len(c_nrs_orig),
hspace=0.75,
wspace=0.1, height_ratios=[1, 1, 0.7],
wspace=0.1,
subplot_spec=grid00[
0]) # 1.2hspace=0.4,wspace=0.2,len(chirps)
grid_down = gridspec.GridSpecFromSubplotSpec(1, 2,
0]) # height_ratios=[1, 1, 0.7], 1.2hspace=0.4,wspace=0.2,len(chirps)
grid_down = gridspec.GridSpecFromSubplotSpec(1, 1,
hspace=0.75,
wspace=0.1,
subplot_spec=grid00[1]) # 1.2hspace=0.4,wspace=0.2,len(chirps)
@ -179,16 +182,17 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
c_dist_recalc=c_dist_recalc)
if not c_dist_recalc:
c_nrs = np.array(c_nrs) * 100
#embed()
letters = ['A', 'B']
letters = ['A', 'B', 'C', 'D']
indexes = [[0, 1, 2, 3]]
#scores_all = [scores]
#embed()
#
scores = ['amp_B1_012_mean_original', 'amp_B2_012_mean_original', 'amp_B1+B2_012_mean_original', 'amp_B1-B2_012_mean_original']
scores = ['c_B1_012_original', 'c_B2_012_original', 'c_B1+B2_012_original', 'c_B1-B2_012_original']
color01, color012, color01_2, color02, color0_burst, color0 = colors_suscept_paper_dots()
colors = [color01, color02, color012, color01_2]
linestyles = ['-','-','-','-']
@ -201,10 +205,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
except:
print('grid search problem4')
embed()
if 'old_fit' in full_name:
ax_u1.set_title('old fit (before 2021)')
else:
ax_u1.set_title('new fit (after 2021)')
plt_single_trace([], ax_u1, frame_cell_orig, freq1, freq2,
scores=np.array(scores)[index], labels=np.array(labels)[index],
colors=np.array(colors)[index],
@ -217,14 +218,14 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
frame_cell = frame_cell_orig[
(frame_cell_orig.df1 == freq1) & (frame_cell_orig.df2 == freq2)]
c1 = c_dist_recalc_here(c_dist_recalc, frame_cell)
ax_u1.set_xlim(0, 50)
ax_u1.set_xlim(0, 25)
if i != 0:
ax_u1.set_ylabel('')
remove_yticks(ax_u1)
if i < 2:
ax_u1.fill_between(c1, frame_cell[np.array(scores)[index][0]],
frame_cell[np.array(scores)[index][1]], color='grey',
alpha=0.1)
#if i < 2:
# ax_u1.fill_between(c1, frame_cell[np.array(scores)[index][0]],
# frame_cell[np.array(scores)[index][1]], color='grey',
# alpha=0.1)
ax_u1.scatter(c_nrs, (np.array(yposs[i]) - 35) * np.ones(len(c_nrs)), color='black',
marker='v',
@ -232,7 +233,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
#
for c_nn, c_nr in enumerate(c_nrs):
ax_u1.text(c_nr, yposs[i][c_nn] + 50, letters[c_nn], color='black', ha='center',
ax_u1.text(c_nr, yposs[i][c_nn] + 15, letters[c_nn], color='black', ha='center',
va='top')
# ax_u1.plot([c_nr, c_nr], [0, 435], color='black', linewidth=0.8, clip_on=False)
@ -251,7 +252,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
axes = []
p_arrays_all = []
model_fit = '_old_fit_' # ''#'_old_fit_'#''#'_old_fit_'#''#'_old_fit_'#''###'_old_fit_'
model_fit = ''#'_old_fit_' # ''#'_old_fit_'#''#'_old_fit_'#''#'_old_fit_'#''###'_old_fit_'
model_cells, reshuffled = chose_old_vs_new_model(model_fit=model_fit)
for c_nn, c_nr in enumerate(c_nrs):
@ -279,9 +280,15 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
nfft_here = 2 ** 14
stimulus_length_here = 0.5
stimulus_length_here = values_stimuluslength_model_full()
trials_nr = 1
nfft_here = values_nfft_full_model()
# #
#
a_f2s = [c_nrs_orig[c_nn]]
_, arrays_spikes, arrays_stim, results_diff, position_diff, auci_wo, auci_w, arrays, names, p_arrays_p, ff_p = calc_roc_amp_core_cocktail(
[freq1 + eodf], [freq2 + eodf], datapoints, auci_wo, auci_w, results_diff,
a_f2s,
@ -294,7 +301,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
n=n,
reshuffled=reshuffled, min_amps=min_amps)
# ff_p, arrays, names, p_arrays_p, arrays_spikes, arrays_stim,
p_arrays_here = p_arrays_p[1::]
p_arrays_here = [p_arrays_p[3]]
xlimp = (0, 300)
for p in range(len(p_arrays_here)):
p_arrays_here[p][0] = p_arrays_here[p][0][ff_p < xlimp[1]]
@ -305,12 +312,11 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
time = time * 1000
# plot the first array
arrays_time = arrays[1::] # [v_mems[1],v_mems[3]]#[1,2]#[1::]
arrays_here = arrays[1::] # [arrays[1],arrays[3]]#arrays[1::]#
arrays_time = [arrays[3]] # [v_mems[1],v_mems[3]]#[1,2]#[1::]
arrays_here = [arrays[3]] # [arrays[1],arrays[3]]#arrays[1::]#
arrays_st = arrays_stim[1::] # [arrays_stim[1],arrays_stim[3]]#
arrays_sp = arrays_spikes[
1::] # [arrays_spikes[1],arrays_spikes[3]]#arrays_spikes[1::]
arrays_st = [arrays_stim[3]] #1:: [arrays_stim[1],arrays_stim[3]]#
arrays_sp = [arrays_spikes[3]] # [arrays_spikes[1],arrays_spikes[3]]#arrays_spikes[1::]
colors_array_here = ['grey', 'grey', 'grey'] # colors_array[1::]
p_arrays_all.append(p_arrays_here)
@ -323,13 +329,14 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
freqs = [np.abs(freq2)]
else:
freqs = [np.abs(freq1), np.abs(freq2)]
grid_pt = gridspec.GridSpecFromSubplotSpec(5, 1,
grid_pt = gridspec.GridSpecFromSubplotSpec(2, 1,
hspace=0.3,
wspace=0.2,
subplot_spec=grid_up[c_nn, a],
height_ratios=[1, 0.7, 1, 0.25,
2.5]) # hspace=0.4,wspace=0.2,len(chirps)
subplot_spec=grid_up[a, c_nn],
) # hspace=0.4,wspace=0.2,len(chirps)
stim = False
if stim:
axe = plt.subplot(grid_pt[0])
axes.append(axe)
plt_stim_saturation(a, [], arrays_st, axe, colors_array_here, f,
@ -342,83 +349,58 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
if not c_dist_recalc:
a_f2_cm = np.array(a_f2_cm) * 100
if a == 2: # if (a_f1s[0] != 0) & (a_f2s[0] != 0):
fish = 'Three fish: $' + label_f_eod_name_core_rm() + '$\,\&\,' + f_vary_name() + '\,\&\,' + f_stable_name() # + '$'#' $\Delta '$\Delta$
beat_here = twobeat_cond(big=True, double=True,
cond=False) + '\,' + f_vary_name(
freq=int(freq1), delta=True) + ',\,$c_{1}=%s$' % (
int(np.round(c_nrs[c_nn]))) + '$\%$' + '\n' + f_stable_name(
freq=int(freq2), delta=True) + ',\,$c_{2}=%s$' % (
#if a == 2: # if (a_f1s[0] != 0) & (a_f2s[0] != 0):
#fish = '3 fish: $' + label_f_eod_name_core_rm() + '$\,\&\,' + f_vary_name() + '\,\&\,' + f_stable_name() # + '$'#' $\Delta '$\Delta$
beat_here = '$c_{1}=%s$' % (
int(np.round(c_nrs[c_nn]))) + '$\%$' + ',\,$c_{2}=%s$' % (
int(np.round(a_f2_cm[0]))) + '$\%$' # +'$'
title_name = fish + '\n' + beat_here # +c1+c2
elif a == 0: # elif (a_f1s[0] != 0):
beat_here = ' ' + onebeat_cond(big=True, double=True,
cond=False) + '\,' + f_vary_name(
freq=int(freq1), delta=True) # +'$' + ' $\Delta '
fish = 'Two fish: $' + label_f_eod_name_core_rm() + '$\,\&\,' + f_vary_name() # +'$'
c1 = ',\,$c_{1}=%s$' % (int(np.round(c_nrs[c_nn]))) + '$\%$ \n '
title_name = fish + '\n' + beat_here + c1 # +'cm'+'cm'+'cm'
elif a == 1: # elif (a_f2s[0] != 0):
beat_here = ' ' + onebeat_cond(big=True, double=True,
cond=False) + '\,' + f_stable_name(
freq=int(freq2),
delta=True) # +'$'
fish = '\n Two fish: $' + label_f_eod_name_core_rm() + '$\,\&\,' + f_stable_name() # +'$'
c1 = ',\,$c_{2}=%s$' % (int(np.round(a_f2_cm[0]))) + '$\%$ \n'
title_name = fish + '\n' + beat_here + c1 # +'cm'
axe.text(1, 1.1, title_name, va='bottom', ha='right',
transform=axe.transAxes)
plt.suptitle(f_vary_name(freq=int(freq1), delta=True)+f_stable_name(freq=int(freq2), delta=True))
title_name = beat_here # fish + '\n' + +c1+c2#twobeat_cond(big=True, double=True,cond=False)
#############################
spikes = False
if spikes:
axs = plt.subplot(grid_pt[1])
plt_spikes_ROC(axs, 'grey', np.array(arrays_sp[a]) * 1000, xlim, lw=1)
#############################
axt = plt.subplot(grid_pt[2])
axt = plt.subplot(grid_pt[0])
axts.append(axt)
plt_vmem_saturation(a, arrays_sp, arrays_time, axt, colors_array_here, f,
time, xlim=xlim)
axt.text(1, 1.1, title_name, va='bottom', ha='right',
transform=axt.transAxes)
axp = plt.subplot(grid_pt[-1])
axps.append(axp)
if a == 0:
if c_nn == 0:
axt.show_spines('')
axt.xscalebar(0.1, -0.1, 10, 'ms', va='right', ha='bottom')
axt.yscalebar(-0.02, 0.35, 600, 'Hz', va='left', ha='top')
axt.xscalebar(0.1, -0.1, 30, 'ms', va='right', ha='bottom')
axt.yscalebar(-0.02, 0.35, 200, 'Hz', va='left', ha='top')
f_counter += 1
if (not os.path.exists(name_psd)) | (redo == True):
''' if ((not os.path.exists(name_psd)) | (redo == True)) & do:
np.save(name_psd, p_arrays_all)
np.save(name_psd_f, ff_p)
else:
ff_p = np.load(name_psd_f) # p_arrays_p
p_arrays_all = np.load(name_psd) # p_arrays_p
p_arrays_all = np.load(name_psd) # p_arrays_p'''
for c_nn, c_nr in enumerate(c_nrs):
for a in range(len(arrays_here)):
axps_here = [[axps[0], axps[1], axps[2]], [axps[3], axps[4], axps[5]]]
axp = axps_here[c_nn][a]
axps_here = [[axps[0], axps[1], axps[2], axps[3]]]#[axps[3], axps[4], axps[5]
axp = axps_here[a][c_nn]
pp = log_calc_psd(log, p_arrays_all[c_nn][a][0],
np.nanmax(p_arrays_all))
markeredgecolors = []
if a == 0:
colors_peaks = [color01] # , 'red']
freqs = [np.abs(freq1)] # ], np.abs(freq2)],
elif a == 1:
colors_peaks = [color02] # , 'red']
freqs = [np.abs(freq2)]
else:
colors_peaks = [color01_012, color02_012] # , 'red']
freqs = [np.abs(freq1), np.abs(freq2)]
markeredgecolors = [color01, color02]
colors_peaks = [color01, color02, color012, color01_2]
markeredgecolors = [color01, color02, color012, color01_2]
freqs = [np.abs(freq1), np.abs(freq2), np.abs(freq1)+np.abs(freq2), np.abs(np.abs(freq1)-np.abs(freq2))]
@ -430,20 +412,20 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
scalebar = False
if scalebar:
axp.show_spines('b')
if a == 0:
if c_nn == 0:
axp.yscalebar(-0.05, 0.5, 20, 'dB', va='center', ha='left')
axp.set_ylim(-33, 5)
else:
axp.show_spines('lb')
if a == 0:
if c_nn == 0:
axp.set_ylabel('dB') # , va='center', ha='left'
else:
remove_yticks(axp)
axp.set_ylim(-39, 5)
else:
axp.show_spines('lb')
if a != 0:
if c_nn != 0:
remove_yticks(axp)
else:
axp.set_ylabel(power_spectrum_name())
@ -470,8 +452,8 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
#join_y(ax_us)
fig = plt.gcf()
fig.tag([[axes[0], axes[1], axes[2]]], xoffs=0, yoffs=3.7)
fig.tag([[axes[3], axes[4], axes[5]]], xoffs=0, yoffs=3.7)
#fig.tag([[axes[0], axes[1], axes[2]]], xoffs=0, yoffs=3.7)
#fig.tag([[axes[3], axes[4], axes[5]]], xoffs=0, yoffs=3.7)
#fig.tag([ax_us[0], ax_us[1], ax_us[2]], xoffs=-2.3, yoffs=1.4)
@ -481,4 +463,4 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
if __name__ == '__main__':
#embed()
sys.excepthook = info
nonlin_regime(yposs = [[170,180], [200,200]],)#, [430,470]
nonlin_regime(yposs = [[270,270, 270, 270]],)#, [430,470], [200,200], [200,200, 200, 200], [200,200, 200, 200], [200,200, 200, 200]

5
susceptibility1.tex
View File

@ -439,6 +439,11 @@ Theoretical work on leaky integrate-and-fire and conductance-based models sugges
\end{figure*}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{nonlin_regime}
\caption{\label{fig:nonlin_regime} }
\end{figure*}
Without any external stimulation, a P-unit fires action potentials at a spontaneous baseline rate \fbase{} to the fish's own EOD of frequency \feod{}. Accordingly, a peak at \fbase{} is present in the power spectrum of this baseline activity (\subfigrefb{fig:motivation}{A}). Superposition of the receiver's EOD with an EOD of another fish with frequency $f_1$ results in a beat, a periodic amplitude modulation of the receiver's EOD. The frequency of the beat is given by the difference frequency $\Delta f_1 = f_1 - \feod$ between the two fish. P-units encode this beat in their firing rate \citep{Bastian1981a,Barayeu2023} and consequently a peak at this beat frequency appears the the power spectrum of the response (\subfigrefb{fig:motivation}{B}). A second peak at the first harmonic indicates a nonlinear response that here is easily identified by the clipping of the P-unit's firing rate at zero. Pairing the fish with another fish with a higher beat frequency $\Delta f_2 = f_2 - \feod$ results in a weaker response with a single peak in the response power spectrum (\subfigrefb{fig:motivation}{C}). Note that $\Delta f_2$ has been choosen to match the P-unit's baseline firing rate.
When stimulating the fish with both frequencies, additional peaks appear in the response power spectrum at the sum \bsum{} and the difference frequency \bdiff{} (\subfigrefb{fig:motivation}{D}). Thus, the response is not equal to sum of the responses to the two beats presented separately. These peaks at the sum and the difference of the two stimulus frequencies are a hallmark of nonlinear interactions that by definition are absent in linear systems.