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

updating nonlin figure

Browse Source
This commit is contained in:
saschuta 2024-06-14 11:15:05 +02:00
parent 29e74f868c
commit 83ab41c4ac
16 changed files with 91 additions and 91 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

BIN
data/nonlin_regime_afe_0.01_nfft_400000_len_100_psd.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.01_nfft_400000_len_100_psdf.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.04_nfft_400000_len_100_psd.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.04_nfft_400000_len_100_psdf.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.1_nfft_400000_len_100_psd.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.1_nfft_400000_len_100_psdf.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.2_nfft_400000_len_100_psd.npy Normal file
View File

Binary file not shown.

BIN
data/nonlin_regime_afe_0.2_nfft_400000_len_100_psdf.npy Normal file
View File

Binary file not shown.

14
data_overview_mod.py
View File

@ -287,14 +287,12 @@ def printing_values_data_overview(cell_types, frame_load_sp, scores, species, x_
if len(frame_file) > 0:
# ja der fisch hatte halt eine ziemlich niedriege EODf
frame_here = frame_file[frame_file.cell != '2022-01-05-af-invivo-1']
try:
print(start_name(cell_type_here,
species) + ' eodf_min ' + str(
int(np.round(np.nanmin(frame_here.eod_fr)))) + ' eodf_max ' + str(
int(np.round(np.nanmax(frame_here.eod_fr)))) + ' n ' + str(len(frame_here)))
except:
print('eodf thing')
embed()
print(start_name(cell_type_here,
species) + ' eodf_min ' + str(
int(np.round(np.nanmin(frame_here.eod_fr)))) + ' eodf_max ' + str(
int(np.round(np.nanmax(frame_here.eod_fr)))) + ' n ' + str(len(frame_here)))
my_list = np.unique(frame_here.cell)
new_list = [item[0:10] for item in my_list]
fish.extend(new_list)

6
model_and_data.py
View File

@ -241,11 +241,7 @@ def model_and_data2(eod_metrice=False, width=0.005, nffts=['whole'], powers=[1],
ax_model = []
for s, sav_name in enumerate(save_names):
try:
ax_external = plt.subplot(grid[nrs_s[s]])
except:
print('vers something')
embed()
ax_external = plt.subplot(grid[nrs_s[s]])
ax_model.append(ax_external)
save_name = find_folder_name('calc_model') + '/' + sav_name

5
model_full.py
View File

@ -439,10 +439,7 @@ def plt_data_full_model(c1, chose_score, detections, devs, dfs, end, grid, mult_
concat=True)
# mt_sorted['m1, m2']
# except:
# print('grouped thing')
# embed()
###################
# groups sorted by repro tag
# todo: evnetuell die tuples gleich hier umspeichern vom csv ''

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: 46 KiB

After

Width:  |  Height:  |  Size: 57 KiB

129
nonlin_regime.py
View File

@ -8,7 +8,7 @@ from matplotlib import gridspec, pyplot as plt
from plotstyle import plot_style
from threefish.calc_time import extract_am
from threefish.core import find_folder_name, info
from threefish.core import find_folder_name, find_project_data, info
from threefish.defaults import default_figsize, default_ticks_talks
from threefish.load import load_savedir, save_visualization
from threefish.plot.limits import join_x, join_y, set_same_ylim
@ -31,11 +31,11 @@ from threefish.twobeat.subplots import plt_psd_saturation, plt_single_trace, plt
from threefish.values import values_nfft_full_model, values_stimuluslength_model_full
def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False, beat='',
nfft_for_morph=4096 * 4,
gain=1,
cells_here=["2013-01-08-aa-invivo-1"], fish_jammer='Alepto', us_name='',
show=True):
def nonlin_regime(yposs=[450, 450, 450], printing=False, beat='',
nfft_for_morph=4096 * 4,
gain=1,
fish_jammer='Alepto', us_name='',
show=True):
runs = 1
n = 1
@ -113,7 +113,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
trials_nr = 20 # 20
redo = False # True
log = 'log' # 'log'
grid0 = gridspec.GridSpec(1, 1, bottom=0.13, top=0.88, left=0.11,
grid0 = gridspec.GridSpec(1, 1, bottom=0.13, top=0.94, left=0.11,
right=0.95, wspace=0.04) #
grid00 = gridspec.GridSpecFromSubplotSpec(2, 1,
wspace=0.04, hspace=0.75,
@ -133,11 +133,6 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
frame_cell_orig = frame[(frame.cell == cell_here)]
if len(frame_cell_orig) > 0:
try:
pass
except:
print('min thing')
embed()
get_frame_cell_params(c_grouped, cell_here, frame, frame_cell_orig)
@ -200,18 +195,12 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
color01, color012, color01_2, color02, color0_burst, color0 = colors_suscept_paper_dots()
colors = [color01, color02, color012, color01_2]
linestyles = ['-','-','-','-']
#frame_cell_orig['amp_B1+B2_012_mean_original']
#frame_cell_orig['amp_B1-B2_012_mean_original']
#for i, index in enumerate(indexes):
index = [0, 1, 2, 3]
try:
ax_u1 = plt.subplot(grid_down[0, i])
except:
print('grid search problem4')
embed()
ax_u1 = plt.subplot(grid_down[0, i])
labels = [label_deltaf1(), label_deltaf2(),
label_sum(), label_diff(), label_fbasename_small()]
#ax_u1.legend(ncol = 4, loc = (0,1.2))
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],
@ -225,14 +214,10 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
(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, 25)
ax_u1.legend(ncol = 2, loc = (0,1.2))
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)
ax_u1.scatter(c_nrs, (np.array(yposs[i]) - 0) * np.ones(len(c_nrs)), color='black',
marker='v',
@ -245,6 +230,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
# ax_u1.plot([c_nr, c_nr], [0, 435], color='black', linewidth=0.8, clip_on=False)
ylim = ax_u1.get_ylim()
ax_u1.set_ylim(0, ylim[-1])
ax_u1.set_xlabel('Contrast [$\%$]')
start = 200 # 1000
mults_period = 3
@ -263,6 +249,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
p_arrays_all = []
model_fit = ''#'_old_fit_' # ''#'_old_fit_'#''#'_old_fit_'#''#'_old_fit_'#''###'_old_fit_'
model_cells, reshuffled = chose_old_vs_new_model(model_fit=model_fit)
axe_all = []
for c_nn, c_nr in enumerate(c_nrs):
@ -272,18 +259,19 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
trials_nr = 1
nfft_here = values_nfft_full_model()
stimulus_length_here = 100 # values_stimuluslength_model_full()
stimulus_length_here = 500 # values_stimuluslength_model_full()
nfft_here = 20 * 20000
a_f2s = [c_nrs_orig[c_nn]]
a_f2s = [c_nrs_orig[c_nn]]
save_dir = load_savedir(level=0).split('/')[0] +'_afe_'+str(a_f2s[0])+'_nfft_'+str(nfft_here)+'_len_' + str(stimulus_length_here)
name_psd = save_dir + '_psd.npy'
name_psd_f = save_dir + '_psdf.npy'
name_psd = find_project_data() + save_dir + '_psd.npy'
name_psd_f = find_project_data() + save_dir + '_psdf.npy'
do = True
if ((not os.path.exists(name_psd)) | (redo == True)) | (do == True):
#embed()
xlimp = (0, 300)
if ((not os.path.exists(name_psd)) | (redo == True)):# | (do == True)
# psd generierung
@ -298,11 +286,19 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
a_f1s=[c_nrs_orig[c_nn]],
n=n,
reshuffled=reshuffled, min_amps=min_amps, mean_choice='first')
np.save(name_psd, p_arrays_p)
p_arrays_here = [p_arrays_p[3]]
# embed()
for p in range(len(p_arrays_here)):
p_arrays_here[p][0] = p_arrays_here[p][0][ff_p < xlimp[1]]
ff_p = ff_p[ff_p < xlimp[1]]
np.save(name_psd, p_arrays_here)
np.save(name_psd_f, ff_p)
else:
ff_p = np.load(name_psd_f) # p_arrays_p
p_arrays_p = np.load(name_psd) # p_arrays_p
p_arrays_here = np.load(name_psd) # p_arrays_p
# #
#
@ -313,7 +309,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
stimulus_length_here = 1 # values_stimuluslength_model_full()
nfft_here = stimulus_length_here * 20000
a_f2s = [c_nrs_orig[c_nn]]
_, arrays_spikes, _, _, _, _, _, _, _, _, _ = calc_roc_amp_core_cocktail(
_, arrays_spikes, arrays_stim, _, _, _, _, arrays, names, _, _ = calc_roc_amp_core_cocktail(
[freq1 + eodf], [freq2 + eodf], datapoints, auci_wo, auci_w, results_diff,
a_f2s,
fish_jammer, trials_nr, nfft_here, us_name, gain, runs, a_fr, nfft_for_morph,
@ -329,13 +325,9 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
##################################################################
# ff_p, arrays, names, p_arrays_p, arrays_spikes, arrays_stim,
p_arrays_here = [p_arrays_p[3]]
xlimp = (0, 300)
#embed()
for p in range(len(p_arrays_here)):
p_arrays_here[p][0] = p_arrays_here[p][0][ff_p < xlimp[1]]
ff_p = ff_p[ff_p < xlimp[1]]
time = np.arange(0, len(arrays[a][0]) / sampling, 1 / sampling)
time = time * 1000
@ -359,19 +351,24 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
else:
freqs = [np.abs(freq1), np.abs(freq2)]
grid_pt = gridspec.GridSpecFromSubplotSpec(3, 1,
hspace=0.2,
hspace=0.3,
wspace=0.2,
subplot_spec=grid_up[a, c_nn], height_ratios = [1, 1,2]
subplot_spec=grid_up[a, c_nn], height_ratios = [1.2, 1,2]
) # hspace=0.4,wspace=0.2,len(chirps)
stim = False
#if stim:
axe = plt.subplot(grid_pt[0])
axes.append(axe)
am, time_am = extract_am(arrays_st, time, extract='',norm = False)
plt_stim_saturation(a, [], am, axe, colors_array_here, f,
f_counter, names, time,
xlim=xlim) # np.array(arrays_sp)*1000
am, time_am = extract_am(arrays_st[0], time, extract='',norm = False)
am = am-np.mean(am)
if c_nn == 0:
spines = 'l'
else:
spines = ''
plt_stim_saturation(0, [], [am*100], axe, colors_array_here, f,
c_nn, ['Contrast [$\%$]'], time,
xlim=xlim, lw = 1, spines = spines) # np.array(arrays_sp)*1000
a_f2_cm = c_dist_recalc_func(frame_cell, c_nrs=[a_f2s[0]], cell=cell_here,
c_dist_recalc=c_dist_recalc)
@ -381,23 +378,30 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
#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]))) + '$\%$' # +'$'
beat_here = 'Contrast$\,=%s$' % (int(np.round(c_nrs[c_nn]))) + '$\%$' # +'$'
#'$c_{1}=%s$' % (int(np.round(c_nrs[c_nn]))) + '$\%$' + ',\,$c_{2}=%s$' % (int(np.round(a_f2_cm[0]))) + '$\%$' # +'$'
plt.suptitle(f_vary_name(freq=int(freq1), delta=True)+', '+f_stable_name(freq=int(freq2), delta=True))
#plt.suptitle(f_vary_name(freq=int(freq1), delta=True)+', '+f_stable_name(freq=int(freq2), delta=True))
print(f_vary_name(freq=int(freq1), delta=True)+', '+f_stable_name(freq=int(freq2)))
title_name = beat_here # fish + '\n' + +c1+c2#twobeat_cond(big=True, double=True,cond=False)
axe.text(1, 1.1, title_name, va='bottom', ha='right',
transform=axe.transAxes)
axe_all.append(axe)
#embed()
#embed()
#if c_nn == 0:
#axe.xscalebar(0.1, 0, 30, 'ms', va='right', ha='bottom')#0.1
#axe.yscalebar(-0.02, 0.35, 5, 'Contrast [$\%$]', va='left', ha='top')
#############################
axs = plt.subplot(grid_pt[1])
plt_spikes_ROC(axs, 'grey', np.array(arrays_sp[a]) * 1000, xlim, lw=1)
if c_nn == 0:
axs.xscalebar(0.15,-0.1, 30, 'ms', va='right', ha='bottom')
axs.text(1, 1.1, title_name, va='bottom', ha='right',
transform=axs.transAxes)
axss.append(axs)
spikes = False
@ -410,7 +414,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
if c_nn == 0:
axt.show_spines('')
axt.xscalebar(0.1, -0.1, 30, 'ms', va='right', ha='bottom')
axt.xscalebar(0.3, -0.1, 30, 'ms', va='right', ha='bottom')
axt.yscalebar(-0.02, 0.35, 200, 'Hz', va='left', ha='top')
axp = plt.subplot(grid_pt[-1])
@ -437,13 +441,17 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
np.nanmax(p_arrays_all))
pps.append(pp)
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))]
colors_peaks = [color01, color02, color012, color01_2, color0]
markeredgecolors = [color01, color02, color012, color01_2, color0]
fr = frame_cell_orig.fr.unique()[0]
labels = [label_deltaf1(), label_deltaf2(),
label_sum(), label_diff(), label_fbasename_small(), label_fbasename_small()]
freqs = [np.abs(freq1), np.abs(freq2), np.abs(freq1)+np.abs(freq2), np.abs(np.abs(freq1)-np.abs(freq2)), fr]
plt_psd_saturation(pp, ff_p, a, axp, colors_array_here, freqs=freqs,
colors_peaks=colors_peaks, xlim=xlimp,
markeredgecolor=markeredgecolors, )
markeredgecolor=markeredgecolors, labels = labels)
if c_nn == 0:
axp.legend(ncol=5, loc=(-0, -1.1))
if log:
scalebar = False
if scalebar:
@ -452,7 +460,6 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
axp.yscalebar(-0.05, 0.5, 20, 'dB', va='center', ha='left')
axp.set_ylim(-33, 5)
else:
axp.show_spines('lb')
if c_nn == 0:
axp.set_ylabel('dB') # , va='center', ha='left'
@ -482,6 +489,8 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
#axts[0].get_shared_x_axes().join(*axts)
axps[0].get_shared_y_axes().join(*axps)
axps[0].get_shared_x_axes().join(*axps)
set_same_ylim(axe_all)
#join_y(axe_all)
#join_y(axts)
#set_same_ylim(axts)
@ -495,7 +504,7 @@ def nonlin_regime(yposs=[450, 450, 450], freqs=[(39.5, -210.5)], printing=False,
#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([axss[0], axss[1], axss[2], axss[3], ax_u1], xoffs=-2.3, yoffs=1.4)
fig.tag([axes[0], axes[1], axes[2], axes[3], ax_u1], xoffs=-2.3, yoffs=1.4)
save_visualization(cell_here, show)

2
susceptibility1.tex
View File

@ -442,7 +442,7 @@ Theoretical work on leaky integrate-and-fire and conductance-based models sugges
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{nonlin_regime.pdf}
\caption{\label{fig:nonlin_regime} The model used has the identifier 2013-01-08-aa.}
\caption{\label{fig:nonlin_regime} The model used has the identifier 2013-01-08-aa.\fone{} is 30\,Hz and \ftwo{} is 130\, Hz, that is eqaul to \fbase{}.}
\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.

26
susceptibility1.tex.bak
View File

@ -405,7 +405,7 @@ Neuronal processing is inherently nonlinear --- spiking thresholds or rectificat
\section{Introduction}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{plot_chi2}
\includegraphics[width=\columnwidth]{plot_chi2.pdf}
\caption{\label{fig:lifresponse} First- (linear) and second-order response functions of the leaky integrate-and-fire model. \figitem{A} Magnitude of the first-order response function $|\chi_1(f_1)|$, also known as ``gain'' quantifies the response amplitude relative to the stimulus amplitude, both measured at the stimulus frequency. \figitem{B} Magnitude of the second-order response function $|\chi_2(f_1, f_2)|$ quantifies the response at the sum of two stimulus frequencies. For linear systems, the second-order response functions is zero, because linear systems do not create new frequencies and thus there is no response at the summed frequencies. The plots show the analytical solutions from \citep{Lindner2001} and \citep{Voronenko2017} with $\mu = 1.1$ and $D = 0.001$.}
\end{figure*}
@ -434,15 +434,15 @@ Theoretical work on leaky integrate-and-fire and conductance-based models sugges
\subsection{Nonlinear responses in P-units stimulated with two beat frequencies}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{motivation}
\includegraphics[width=\columnwidth]{motivation.pdf}
\caption{\label{fig:motivation} Nonlinearity in an electrophysiologically recorded P-unit of \lepto{} in a three-fish setting (cell identifier ``2021-08-03-ac"). Receiver with EOD frequency $\feod{} =664$\,Hz encounters fish with EOD frequencies $f_{1}=631$\,Hz and $f_{2}=797$\,Hz. Both encountered fish lead to a beat contrast of 10\,\%. Top: Sheme of a nonlinear system. Second row: Interference of the receiver EOD with the EODs of other fish. Third row: Spike trains of the P-unit. Forth row: Firing rate, retrieved as the convolution of the spike trains with a Gaussian kernel ($\sigma = 1$\,ms). Bottom row: Power spectrum of the firing rate. \figitem{A} Baseline condition: Only the receiver is present. The baseline firing rate \fbase{} dominates the power spectrum of the firing rate. \figitem{B} The receiver and the fish with EOD frequency $f_{1}=631$\,Hz are present. \figitem{C} The receiver and the fish with EOD frequency $f_{2}=797$\,Hz are present. \figitem{D} All three fish with the EOD frequencies \feod{}, $f_{1}$ and $f_{2}$ are present. Nonlinear peaks occur at the sum and difference of the two beat frequencies in the power spectrum of the firing rate.
}
\end{figure*}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{nonlin_regime}
\caption{\label{fig:nonlin_regime} The model used has the identifier 2013-01-08-aa.}
\includegraphics[width=\columnwidth]{nonlin_regime.pdf}
\caption{\label{fig:nonlin_regime} The model used has the identifier 2013-01-08-aa.\fone{} is 30\,Hz and \ftwo{} is 130\, Hz, that is eqaul to \fbase{}. Both contrasts increase equally in strength.}
\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.
@ -456,7 +456,7 @@ Weakly nonlinear responses are expected in cells with sufficiently low intrinsic
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{cells_suscept}
\includegraphics[width=\columnwidth]{cells_suscept.pdf}
\caption{\label{fig:cells_suscept} Linear and nonlinear stimulus encoding in a low-CV P-unit (cell identifier ``2010-06-21-ai"). \figitem{A} Interspike interval (ISI) distribution of the cell's baseline activity, i.e. the cell is driven only by the unperturbed own electric field. The low CV of the ISIs indicates quite regular firing. \figitem{B} Power spectral density of the baseline response with peaks at the cell's baseline firing rate \fbase{} and the fish's EOD frequency \feod{}. \figitem{C} Random amplitude modulation stimulus (top, with cutoff frequency of 300\,Hz) and evoked responses (spike raster, bottom) of the same P-unit. The stimulus contrast (right) measures the strength of the AM. \figitem{D} Gain of the transfer function (first-order susceptibility), \Eqnref{linearencoding_methods}, computed from the responses to 10\,\% (light purple) and 20\,\% contrast (dark purple) RAM stimulation of 10\,s duration. \figitem{E, F} Absolute value of the second-order susceptibility, \Eqnref{eq:susceptibility}, for both the low and high stimulus contrast. Pink triangles mark vertical, horizontal, and diagonal lines where \fone, \ftwo{} or \fsum{} are equal to \fbase{}. \figitem{G} Second-order susceptibilities projected onto the diagonal (means of all anti-diagonals of the matrices shown in \panel{E, F}). Dots mark \fbase{}, horizontal dashed lines mark medians of the projected susceptibilities.}
\end{figure*}
@ -472,7 +472,7 @@ In contrast, a high-CV P-unit (CV$_{\text{base}}=0.4$) does not exhibit pronounc
\subsection{Ampullary afferents exhibit strong nonlinear interactions}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{ampullary}
\includegraphics[width=\columnwidth]{ampullary.pdf}
\caption{\label{fig:ampullary} Linear and nonlinear stimulus encoding in an ampullary afferent (cell identifier ``2012-04-26-ae"). \figitem{A} Interspike interval (ISI) distribution of the cell's baseline activity. The very low CV of the ISIs indicates almost perfect periodic spiking. \figitem{B} Power spectral density of baseline activity with peaks at the cell's baseline firing rate and its harmonics. \figitem{C} Bad-limited white noise stimulus (top, with cutoff frequency of 150\,Hz) added to the fish's self-generated electric field and spike raster of the evoked responses (bottom) for two stimulus contrasts as indicated (right). \figitem{D} Gain of the transfer function, \Eqnref{linearencoding_methods}, of the responses to stimulation with 2\,\% (light green) and 20\,\% contrast (dark green) of 10\,s duration. \figitem{E, F} Absolute value of the second-order susceptibility, \Eqnref{eq:susceptibility}, for both stimulus contrasts as indicated. Pink triangles indicate baseline firing rate. \figitem{G} Projections of the second-order susceptibilities in \panel{E, F} onto the diagonal. }
\end{figure*}
@ -484,8 +484,8 @@ Irrespective of the CV, neither of the two example P-units shows the complete ex
In the example recordings shown above (\figsrefb{fig:cells_suscept} and \fref{fig:ampullary}), we only observe nonlinear responses on the anti-diagnal of the second-order susceptibility, where the sum of the two stimulus frequencies matches the neuron's baseline firing rate, which is in line with theoretical expectations \citep{Voronenko2017}. However, a pronounced nonlinear response at frequencies \foneb{} or \ftwob{}, although predicted by theory, cannot be observed. Here we investigate how these discrepancies can be understood.
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{model_and_data}
\caption{\label{model_and_data} Estimation of second-order susceptibilities in the limit of weak stimuli. \figitem{A} \suscept{} estimated from $N=11$ 0.5\,slong trials of an electrophysiological recording of another low-CV P-unit (cell 2012-07-03-ak, $\fbase=120$\,Hz, CV=0.20) driven with a weak RAM stimulus with contrast 2.5\,\%. Pink edges mark baseline firing rate where enhanced nonlinear responses are expected. \figitem[i]{B} \textit{Standard condition} of model simulations with intrinsic noise (bottom) and a RAM stimulus (top). \notejb{Since the model overestimated the sensitivity of the real P-unit, we adjusted the RAM contrast to 0.009\,\%, such that the resulting spike trains had the same CV as the electrophysiolgical recorded P-unit during the 2.5\,\% contrast stimulation (see table~\ref{modelparams} for model parameters).} \figitem[ii]{B} \suscept{} estimated from simulations of the cell's LIF model counterpart (cell 2012-07-03-ak, table~\ref{modelparams}) based on the same number of trials as in the electrophysiological recording. \figitem[iii]{B} Same as \panel[ii]{B} but using $10^6$ stimulus repetitions. \figitem[i-iii]{C} Same as in \panel[i-iii]{B} but in the \textit{noise split} condition: there is no external RAM signal driving the model. Instead, a large part (90\,\%) of the total intrinsic noise is treated as signal and is presented as an equivalent amplitude modulation (\signalnoise, center), while the intrinsic noise is reduced to 10\,\% of its original strength (see methods for details). In addition to one million trials, this reveals the full expected structure of the second-order susceptibility.}
\includegraphics[width=\columnwidth]{model_and_data.pdf}
\caption{\label{model_and_data} Estimation of second-order susceptibilities in the limit of weak stimuli. \figitem{A} \suscept{} estimated from $N=11$ 0.5\,slong trials of an electrophysiological recording of another low-CV P-unit (cell 2012-07-03-ak, $\fbase=120$\,Hz, CV=0.20) driven with a weak RAM stimulus with contrast 2.5\,\%. Pink edges mark baseline firing rate where enhanced nonlinear responses are expected. \figitem[i]{B} \textit{Standard condition} of model simulations with intrinsic noise (bottom) and a RAM stimulus (top). \notejb{Since the model overestimated the sensitivity of the real P-unit, we adjusted the RAM contrast to 0.9\,\%, such that the resulting spike trains had the same CV as the electrophysiolgical recorded P-unit during the 2.5\,\% contrast stimulation (see table~\ref{modelparams} for model parameters).} \notejb{chi2 scale is higher than in real cell} \figitem[ii]{B} \suscept{} estimated from simulations of the cell's LIF model counterpart (cell 2012-07-03-ak, table~\ref{modelparams}) based on the same number of trials as in the electrophysiological recording. \figitem[iii]{B} Same as \panel[ii]{B} but using $10^6$ stimulus repetitions. \figitem[i-iii]{C} Same as in \panel[i-iii]{B} but in the \textit{noise split} condition: there is no external RAM signal driving the model. Instead, a large part (90\,\%) of the total intrinsic noise is treated as signal and is presented as an equivalent amplitude modulation (\signalnoise, center), while the intrinsic noise is reduced to 10\,\% of its original strength (see methods for details). In addition to one million trials, this reveals the full expected structure of the second-order susceptibility.}
\end{figure*}
One reason could be simply too little data for a good estimate of the second-order susceptibility. Electrophysiological recordings are limited in time, and therefore responses to only a limited number of trials, i.e. repetitions of the same RAM stimulus, are available. As a consequence, the cross-spectra, \Eqnref{eq:crosshigh}, are insufficiently averaged and the full structure of the second-order susceptibility might be hidden in finite-data noise. This experimental limitation can be overcome by using a computational model for the P-unit, a stochastic leaky integrate-and-fire model with adaptation current and dendritic preprocessing, and parameters fitted to the experimentally recorded P-unit (\figrefb{flowchart}) \citep{Barayeu2023}. The model faithfully reproduces the second-order susceptibility of a low-CV cell estimated from the same low number of trials as in the experiment ($\n{}=11$, compare \panel{A} and \panel[ii]{B} in \figrefb{model_and_data}).
@ -503,7 +503,7 @@ With high levels of intrinsic noise, we would not expect the nonlinear response
We estimated the second-order susceptibility of P-unit responses using RAM stimuli. In particular, we found pronounced nonlinear responses in the limit of weak stimulus amplitudes. How do these findings relate to the situation of two pure sinewave stimuli with finite amplitudes that approximates the interference of EODs of real animals? For the P-units the relevant signals are the beat frequencies \bone{} and \btwo{} that arise from the interference of either of the two foreign EODs with the receiving fish's own EOD (\figref{fig:motivation}). In the introductory example, the response power spectrum showed peaks from nonlinear interactions at the sum of the two beat frequencies (orange marker, \subfigrefb{fig:motivation}{D}) and at the difference between the two beat frequencies (red marker, \subfigrefb{fig:motivation}{D}). In this example, $\Delta f_{2}$ was similar to \fbase{}, corresponding to the horizontal line of the second-order susceptibility estimated for a vanishing external RAM stimulus (\subfigrefb{model_and_data}\,\panel[iii]{C}). In the three-fish example, there was a second nonlinearity at the difference between the two beat frequencies (red dot, \subfigrefb{fig:motivation}{D}), that is not covered by the so-far shown part of the second-order susceptibility (\subfigrefb{model_and_data}\,\panel[iii]{C}), in which only the response at the sum of the two stimulus frequencies is addressed. % less prominent,
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{model_full}
\includegraphics[width=\columnwidth]{model_full.pdf}
\caption{\label{fig:model_full} Using second-order susceptibility to predict responses to sine-wave stimuli. \figitem[]{A} Absolute value of the second-order susceptibility, \Eqnref{eq:susceptibility}, for both positive and negative frequencies. \susceptf{} was estimated from $N=10^6$ trials of model simulations in the noise-split condition (cell 2013-01-08-aa, see table~\ref{modelparams} for model parameters). White lines indicate zero frequencies. Nonlinear responses at \fsum{} are quantified in the upper right and lower left quadrants. Nonlinear responses at \fdiff{} are quantified in the upper left and lower right quadrants. Baseline firing rate of this cell was at $\fbase=120$\,Hz. The position of the orange/red letters correspond to the beat frequencies used for the stimulation with pure sine-waves in the subsequent panels and indicate the sum/difference of those beat frequencies. \figitem[]{B} Absolute value of the first-order susceptibility. \figitem{C--F} Black line -- power spectral density of model simulations in response to stimulation with two pure sine waves, \fone{} and \ftwo, in addition to the receiving fish's own EOD (three fish scenario). The contrasts of beat beats is 0.0065. Colored circles highlight the height of selected peaks in the power spectrum. Black circles highlight the peak height that can be predicted from \panel{A, B}. Grey line -- power spectral density of model in the baseline condition. \figitem{B} The sum of the two beat frequencies match \fbase{}. \figitem{C} The difference of \fone{} and \ftwo{} match \fbase{}. \figitem{D} Only the first beat frequency matches \fbase{}. \figitem{C} None of the two beat frequencies matches \fbase{}.}
\end{figure*}
@ -511,7 +511,7 @@ However, the second-order susceptibility \Eqnref{eq:susceptibility} is a spectra
Is it possible based on the second-order susceptibility estimated by means of RAM stimuli (\subfigrefb{fig:model_full}{A}) to predict nonlinear responses in a three-fish setting? We can test this by stimulating the same model with two beats with weak amplitudes (\subfigrefb{fig:model_full}{C--F}). If we choose a frequency combination where the sum of the two beat frequencies is equal to the model's baseline firing rate \fbase{}, a peak at the sum of the two beat frequencies appears in the power spectrum of the response (\subfigrefb{fig:model_full}{C}), as expected from \suscept. If instead we choose two beat frequencies that differ by \fbase{}, a peak is present at the difference frequency (\subfigrefb{fig:model_full}{D}). If only one beat frequency is equal to \fbase{}, both a peak at the sum and at the difference frequency is present in the P-unit response (\subfigrefb{fig:model_full}{E}). And if none of these conditions are met, neither a peak at the sum nor at the difference of the two beat frequencies appears (\subfigrefb{fig:model_full}{F}).
\begin{figure*}[tp]
\includegraphics[width=\columnwidth]{data_overview_mod}
\includegraphics[width=\columnwidth]{data_overview_mod.pdf}
\caption{\label{fig:data_overview} Nonlinear responses in P-units and ampullary cells. The second-order susceptibility is condensed into the peakedness of the nonlinearity, \nli{} \Eqnref{eq:nli_equation}, that relates the amplitude of the projected susceptibility at a cell's baseline firing rate to its median (see \subfigrefb{fig:cells_suscept}{G}). Each of the recorded neurons contributes at maximum with two stimulus contrasts. Black squares and circles highlight recordings conducted in a single cell. Squares in \panel{A, C, E} correspond to the cell in \figrefb{fig:cells_suscept} and circles to the cell in \figrefb{fig:cells_suscept_high_CV}. Squares in \panel{B, D, F} correspond to the cell in \figrefb{fig:ampullary}. \figitem{A, B} There is a negative correlation between the CV during baseline and \nli. \figitem{C, D} There is a negative correlation between the CV during stimulation and \nli. \figitem{E, F} \nli{} is plotted against the response modulation, (see methods), an indicator of the subjective stimulus strength for a cell. There is a negative correlation between response modulation and \nli. Restricting the analysis to the weakest stimulus that was presented to each unique neuron, does not change the results. The number of unique neurons is 221 for P-units and 45 for ampullary cells.
% The two example P-units shown before are highlighted with dark markers in \subfigrefb{fig:data_overview}{A, C, E} (squares -- \figrefb{fig:cells_suscept}, circles -- \figrefb{fig:cells_suscept_high_CV}).
% Several of the recorded neurons contribute with two samples to the population analysis as their responses have been recorded to two different contrast of the same RAM stimulus. Higher stimulus contrasts lead to a stronger drive and thus stronger response modulations (see color code bar in \subfigref{fig:data_overview}{A}, see methods).
@ -766,7 +766,7 @@ Whenever the membrane voltage $V_m(t)$ crossed the spiking threshold $\theta=1$
%, \figrefb{eif}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{flowchart}
\includegraphics[width=\columnwidth]{flowchart.pdf}
\caption{\label{flowchart}
Components of the P-unit model. The main steps of the model are illustrated in the left column. The three other columns show the corresponding signals in three different settings: (i) the baseline situation, no external stimulus, only the fish's self-generated EOD (i.e. the carrier) is present. (ii) RAM stimulation, the carrier is amplitude modulated with a weak (2\,\% contrast) nad-limited white-noise stimulus. (iii) Noise split condition in which 90\,\% of the internal noise is used as a driving RAM stimulus scaled with the correction factor $\rho$ (see text). Note that the mean firing rate and the CV of the ISI distribution is the same in this and the baseline condition. As an example, simulations of the model for cell ``2012-07-03-ak'' are shown (see table~\ref{modelparams} for model parameters). \figitem{A} Thresholding: a simple linear threshold was applied to the EOD carrier, \Eqnref{eq:eod}. \figitem{B} Subsequent dendritic low-pass filtering attenuates the carrier and carves out the AM signal. \figitem{C} Gaussian white-noise is added to the signal in \panel{B}. Note the reduced internal noise amplitude in the noise split (iii) condition. \figitem{D} Spiking output of the LIF model in response to the sum of B and C. \figitem{E} Power spectra of the LIF neuron's spiking activity. Under the baseline condition (\panel[i]{E}) there are several peaks, from left to right, at the baseline firing rate $\fbase{}$, $f_{EOD} - \fbase{}$, $f_{EOD}$, and $f_{EOD} + \fbase{}$. In the stimulus driven regime (\panel[ii]{E}), there is only a peak at \feod, while under the noise split condition (\panel[iii]{E}) the same peaks as in the baseline condition are present.}
\end{figure*}
@ -878,13 +878,13 @@ CVs in P-units can range up to 1.5 \citep{Grewe2017, Hladnik2023}. We show the s
\label{S1:highcvpunit}
\begin{figure*}[!ht]
\includegraphics[width=\columnwidth]{cells_suscept_high_CV}
\includegraphics[width=\columnwidth]{cells_suscept_high_CV.pdf}
\caption{\label{fig:cells_suscept_high_CV} Response of experimentally measured noisy P-units (cell identifier ``2018-08-24-af") with a relatively high CV of 0.34 to RAM stimuli with two different contrasts. \figitem{A} Interspike intervals (ISI) distribution during baseline. \figitem{B} Baseline power spectrum. \figitem{C} Top: EOD carrier (gray) with RAM (red). Center: Spike trains in response to the 5\,\% RAM contrast. Bottom: Spike trains in response to the 10\,\% RAM contrast. \figitem{D} First-order susceptibility (\Eqnref{linearencoding_methods}). \figitem{E} Absolute value $|\chi_2(f_1, f_2)|$ of the second-order susceptibility, \Eqnref{eq:susceptibility}, for the 5\,\% RAM contrast. Pink lines -- edges of the structure when \fone, \ftwo{} or \fsum{} are equal to \fbase{}. \figitem{F} $|\chi_2(f_1, f_2)|$ for the 10\,\% RAM contrast. \figitem{G} Projected diagonals, calculated as the mean of the anti-diagonals of the matrices in \panel{E--F}. Gray dots: \fbase{}. Dashed lines: Medians of the projected diagonals.}
\end{figure*}
\begin{figure*}[t]
\includegraphics[width=\columnwidth]{trialnr}
\includegraphics[width=\columnwidth]{trialnr.pdf}
\caption{\label{fig:trialnr} Dependence of the estimate of the second-order susceptibility on the number of trials $\n$. While the estimate of the noise floor (10th and 90th percentile) of the $|\chi_2(f_1, f_2)|$ matrix does not saturate yet, the estimates of the high values in the matrix that make up the characteristic ridges saturate for $N>10^6$. The model used has the identifier 2013-01-08-aa.
}
\end{figure*}