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

updated captions to new figures

Browse Source
This commit is contained in:
Jan Benda 2025-05-26 11:36:02 +02:00
parent d043a8c47c
commit 04578fd77e
3 changed files with 80 additions and 69 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

68
dataoverview.py
View File

@ -29,55 +29,57 @@ ampul_examples = (ampul_example, [], ampul_examples)
def plot_corr(ax, data, xcol, ycol, zcol, zmin, zmax, xpdfmax, cmap, color,
si_thresh, example=[], split_example=[], examples=[]):
xdata = data[xcol]
ydata = data[ycol]
ax.axhline(si_thresh, color='k', ls=':', lw=0.5)
xmax = ax.get_xlim()[1]
ymax = ax.get_ylim()[1]
mask = (data[xcol] < xmax) & (data[ycol] < ymax)
mask = (xdata < xmax) & (ydata < ymax)
if 'stimindex' in data:
for cell, run in example + split_example + examples:
mask &= ~((data['cell'] == cell) & (data['stimindex'] == run))
else: # simulations
for cell, alpha in example + split_example + examples:
mask &= ~((data['cell'] == cell) & (data['contrast'] == alpha))
sc = ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
sc = ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=4, marker='o', linewidth=0, edgecolors='none',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax, zorder=20)
elw = 0.3
if 'stimindex' in data:
for cell, run in example:
mask = (data['cell'] == cell) & (data['stimindex'] == run)
ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=6, marker='^', linewidth=elw, edgecolors='black',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax,
zorder=20)
for cell, run in split_example:
mask = (data['cell'] == cell) & (data['stimindex'] == run)
ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=6, marker='s', linewidth=elw, edgecolors='black',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax,
zorder=20)
for cell, run in examples:
mask = (data['cell'] == cell) & (data['stimindex'] == run)
ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=6, marker='o', linewidth=elw, edgecolors='black',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax,
zorder=20)
else: # simulations
for cell, alpha in example:
mask = (data['cell'] == cell) & (data['contrast'] == alpha)
ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=6, marker='^', linewidth=elw, edgecolors='black',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax,
zorder=20)
for cell, alpha in split_example:
mask = (data['cell'] == cell) & (data['contrast'] == alpha)
ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=6, marker='s', linewidth=elw, edgecolors='black',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax,
zorder=20)
for cell, alpha in examples:
mask = (data['cell'] == cell) & (data['contrast'] == alpha)
ax.scatter(data[mask, xcol], data[mask, ycol], c=data[mask, zcol],
ax.scatter(xdata[mask], ydata[mask], c=data[mask, zcol],
s=6, marker='o', linewidth=elw, edgecolors='black',
clip_on=False, cmap=cmap, vmin=zmin, vmax=zmax,
zorder=20)
@ -89,7 +91,7 @@ def plot_corr(ax, data, xcol, ycol, zcol, zmin, zmax, xpdfmax, cmap, color,
cb.outline.set_color('none')
cb.outline.set_linewidth(0)
# pdf x-axis:
kde = gaussian_kde(data[xcol], 0.02*xmax/np.std(data[xcol], ddof=1))
kde = gaussian_kde(xdata, 0.02*xmax/np.std(xdata, ddof=1))
xx = np.linspace(0, ax.get_xlim()[1], 400)
pdf = kde(xx)
xax = ax.inset_axes([0, 1.05, 1, 0.2])
@ -99,7 +101,7 @@ def plot_corr(ax, data, xcol, ycol, zcol, zmin, zmax, xpdfmax, cmap, color,
xax.set_ylim(bottom=0)
xax.set_ylim(0, xpdfmax)
# pdf y-axis:
kde = gaussian_kde(data[ycol], 0.02*ymax/np.std(data[ycol], ddof=1))
kde = gaussian_kde(ydata, 0.02*ymax/np.std(ydata, ddof=1))
xx = np.linspace(0, ax.get_ylim()[1], 400)
pdf = kde(xx)
yax = ax.inset_axes([1.05, 0, 0.2, 1])
@ -109,12 +111,12 @@ def plot_corr(ax, data, xcol, ycol, zcol, zmin, zmax, xpdfmax, cmap, color,
yax.set_xlim(left=0)
# threshold:
if 'cvbase' in xcol:
ax.text(xmax, 0.4*ymax, f'{100*np.sum(data[ycol] > si_thresh)/len(data):.0f}\\%',
ax.text(xmax, 0.4*ymax, f'{100*np.sum(ydata > si_thresh)/len(data):.0f}\\%',
ha='right', va='bottom', fontsize='small')
ax.text(xmax, 0.3, f'{100*np.sum(data[ycol] < si_thresh)/len(data):.0f}\\%',
ax.text(xmax, 0.3, f'{100*np.sum(ydata < si_thresh)/len(data):.0f}\\%',
ha='right', va='center', fontsize='small')
# statistics:
r, p = pearsonr(data[xcol], data[ycol])
r, p = pearsonr(xdata, ydata)
ax.text(1, 0.9, f'$R={r:.2f}$ **', ha='right',
transform=ax.transAxes, fontsize='small')
#ax.text(1, 0.77, f'{significance_str(p)}', ha='right',
@ -128,16 +130,17 @@ def plot_corr(ax, data, xcol, ycol, zcol, zmin, zmax, xpdfmax, cmap, color,
def si_stats(title, data, sicol, si_thresh, nsegscol):
print(title)
sidata = data[sicol]
cells = np.unique(data['cell'])
ncells = len(cells)
nrecs = len(data)
print(f' cells: {ncells}')
print(f' recordings: {nrecs}')
print(f' SI threshold: {si_thresh:.1f}')
hcells = np.unique(data[data(sicol) > si_thresh, 'cell'])
hcells = np.unique(data[sidata > si_thresh, 'cell'])
print(f' high SI cells: n={len(hcells):3d}, {100*len(hcells)/ncells:4.1f}%')
print(f' high SI recordings: n={np.sum(data(sicol) > si_thresh):3d}, '
f'{100*np.sum(data(sicol) > si_thresh)/nrecs:4.1f}%')
print(f' high SI recordings: n={np.sum(sidata > si_thresh):3d}, '
f'{100*np.sum(sidata > si_thresh)/nrecs:4.1f}%')
nsegs = data[nsegscol]
print(f' number of segments: {np.min(nsegs):4.0f} - {np.max(nsegs):4.0f}, median={np.median(nsegs):4.0f}, mean={np.mean(nsegs):4.0f}, std={np.std(nsegs):4.0f}')
nrecs = []
@ -155,9 +158,10 @@ def si_stats(title, data, sicol, si_thresh, nsegscol):
def plot_cvbase_si_punit(ax, data, ycol, si_thresh, color):
ax.set_xlabel('CV$_{\\rm base}$')
ax.set_ylabel('SI($r$)')
ax.set_xlim(0, 1.5)
ax.set_ylim(0, 7.2)
ax.set_xticks_delta(0.5)
ax.set_ylabel('SI($r$)')
ax.set_ylim(0, 6.5)
ax.set_yticks_delta(2)
examples = punit_examples if 'stimindex' in data else model_examples
cax = plot_corr(ax, data, 'cvbase', ycol, 'respmod2', 0, 250, 3,
@ -167,10 +171,10 @@ def plot_cvbase_si_punit(ax, data, ycol, si_thresh, color):
def plot_cvstim_si_punit(ax, data, ycol, si_thresh, color):
ax.set_xlabel('CV$_{\\rm stim}$')
ax.set_ylabel('SI($r$)')
ax.set_xlim(0, 1.6)
ax.set_ylim(0, 7.2)
ax.set_xticks_delta(0.5)
ax.set_ylabel('SI($r$)')
ax.set_ylim(0, 6.5)
ax.set_yticks_delta(2)
examples = punit_examples if 'stimindex' in data else model_examples
#cax = plot_corr(ax, data, 'cvstim', ycol, 'respmod2', 0, 250, 2,
@ -189,9 +193,10 @@ def plot_cvstim_si_punit(ax, data, ycol, si_thresh, color):
def plot_rmod_si_punit(ax, data, ycol, si_thresh, color):
ax.set_xlabel('Response modulation', 'Hz')
ax.set_ylabel('SI($r$)')
ax.set_xlim(0, 250)
ax.set_ylim(0, 7.2)
ax.set_xticks_delta(100)
ax.set_ylabel('SI($r$)')
ax.set_ylim(0, 6.5)
ax.set_yticks_delta(2)
examples = punit_examples if 'stimindex' in data else model_examples
cax = plot_corr(ax, data, 'respmod2', ycol, 'cvbase', 0, 1.5, 0.016,
@ -201,11 +206,11 @@ def plot_rmod_si_punit(ax, data, ycol, si_thresh, color):
def plot_cvbase_si_ampul(ax, data, ycol, si_thresh, color):
ax.set_xlabel('CV$_{\\rm base}$')
ax.set_ylabel('SI($r$)')
ax.set_xlim(0, 0.2)
ax.set_ylim(0, 15)
ax.set_xticks_delta(0.1)
ax.set_yticks_delta(5)
ax.set_ylabel('SI($r$)')
ax.set_ylim(0, 10)
ax.set_yticks_delta(2)
cax = plot_corr(ax, data, 'cvbase', ycol, 'respmod2', 0, 80, 20,
'coolwarm', color, si_thresh, *ampul_examples)
cax.set_ylabel('Response mod.', 'Hz')
@ -213,11 +218,11 @@ def plot_cvbase_si_ampul(ax, data, ycol, si_thresh, color):
def plot_cvstim_si_ampul(ax, data, ycol, si_thresh, color):
ax.set_xlabel('CV$_{\\rm stim}$')
ax.set_ylabel('SI($r$)')
ax.set_xlim(0, 0.85)
ax.set_ylim(0, 15)
ax.set_xticks_delta(0.2)
ax.set_yticks_delta(5)
ax.set_ylabel('SI($r$)')
ax.set_ylim(0, 10)
ax.set_yticks_delta(2)
#cax = plot_corr(ax, data, 'cvstim', ycol, 'respmod2', 0, 80, 6,
# 'coolwarm', color, si_thresh, *ampul_examples)
#cax.set_ylabel('Response mod.', 'Hz')
@ -233,11 +238,11 @@ def plot_cvstim_si_ampul(ax, data, ycol, si_thresh, color):
def plot_rmod_si_ampul(ax, data, ycol, si_thresh, color):
ax.set_xlabel('Response modulation', 'Hz')
ax.set_ylabel('SI($r$)')
ax.set_xlim(0, 80)
ax.set_ylim(0, 15)
ax.set_xticks_delta(20)
ax.set_yticks_delta(5)
ax.set_ylabel('SI($r$)')
ax.set_ylim(0, 10)
ax.set_yticks_delta(2)
cax = plot_corr(ax, data, 'respmod2', ycol, 'cvbase', 0, 0.2, 0.06,
'coolwarm', color, si_thresh, *ampul_examples)
cax.set_ylabel('CV$_{\\rm base}$')
@ -255,6 +260,7 @@ if __name__ == '__main__':
ampul_data = TableData(data_path /
'Apteronotus_leptorhynchus-Ampullary-data.csv',
sep=';')
#si = ''
si = '_nmax'
si_thresh = 1.8

73
noisesplit.py
View File

@ -65,14 +65,43 @@ def plot_overn(ax, s, files, nmax=1e6):
ax.set_xlabel('segments')
ax.set_ylabel(r'$|\chi_2|$', r'Hz/\%$^2$')
def plot_chi2_data(ax, s, cell_name, run):
vmax = 15
data_file = data_path / f'{cell_name}-baseline.npz'
data = np.load(data_file)
eodf = float(data['eodf'])
ratebase = float(data['ratebase/Hz'])
cvbase = float(data['cvbase'])
data_file = data_path / f'{cell_name}-spectral-100-s{run:02d}.npz'
data = np.load(data_file)
nsegs = data['nsegs']
fcutoff = data['fcutoff']
nfft = data['nfft']
deltat = data['deltat']
alpha = data['contrast']
freqs = data['freqs']
pss = data['pss']
prss = data['prss']
chi2 = np.abs(prss)*0.5/(pss.reshape(1, -1)*pss.reshape(-1, 1))
print(f'Measured cell {"-".join(data_file.name.split("-")[:-3])} at {100*alpha:4.1f}% contrast:')
print(f' r={ratebase:3.0f}Hz, CV={cvbase:4.2f}, dt={1000*deltat:4.2f}ms, nfft={nfft}, win={1000*deltat*nfft:6.1f}ms, nsegs={nsegs}')
print()
ax.text(1, 1.1, f'$N={nsegs}$',
ha='right', transform=ax.transAxes)
plot_chi2(ax, s, freqs, chi2, fcutoff, ratebase, vmax)
return alpha, ratebase, eodf
def plot_chi2_contrast(ax1, ax2, s, files, nums, nsmall, nlarge, rate):
vmax = {0.05: {nsmall: 15, nlarge: 1.2},
0.01: {nsmall: 400, nlarge: 4}}
for ax, n in zip([ax1, ax2], [nsmall, nlarge]):
i = nums.index(n)
data = np.load(files[i])
nsegs = data['nsegs']
fcutoff = data['fcutoff']
alpha = data['contrast']
fcutoff = float(data['fcutoff'])
alpha = float(data['contrast'])
freqs = data['freqs']
pss = data['pss']
prss = data['prss']
@ -83,19 +112,20 @@ def plot_chi2_contrast(ax1, ax2, s, files, nums, nsmall, nlarge, rate):
ax.text(1, 1.1, f'$N=10^{{{np.log10(nsegs):.0f}}}$',
ha='right', transform=ax.transAxes)
chi2 = np.abs(prss)*0.5/(pss.reshape(1, -1)*pss.reshape(-1, 1))
cax = plot_chi2(ax, s, freqs, chi2, fcutoff, rate)
cax = plot_chi2(ax, s, freqs, chi2, fcutoff, rate, vmax[alpha][n])
cax.set_ylabel('')
print(f'Modeled cell {"-".join(files[i].name.split("-")[2:-2])} at {100*alpha:4.1f}% contrast: noise_frac={1:3.1f}, nsegs={n}')
print(f'Modeled cell {"-".join(files[i].name.split("-")[:-4])} at {100*alpha:4.1f}% contrast: noise_frac={1:3.1f}, nsegs={n}')
print()
def plot_chi2_split(ax1, ax2, s, files, nums, nsmall, nlarge, rate):
vmax = {nsmall: 4, nlarge: 1.2}
for ax, n in zip([ax1, ax2], [nsmall, nlarge]):
i = nums.index(n)
data = np.load(files[i])
nsegs = data['nsegs']
fcutoff = data['fcutoff']
alpha = data['contrast']
fcutoff = float(data['fcutoff'])
alpha = float(data['contrast'])
noise_frac = data['noise_frac']
freqs = data['freqs']
pss = data['pss']
@ -107,37 +137,12 @@ def plot_chi2_split(ax1, ax2, s, files, nums, nsmall, nlarge, rate):
else:
ax.text(1, 1.1, f'$N=10^{{{np.log10(nsegs):.0f}}}$',
ha='right', transform=ax.transAxes)
cax = plot_chi2(ax, s, freqs, chi2, fcutoff, rate)
cax = plot_chi2(ax, s, freqs, chi2, fcutoff, rate, vmax[n])
cax.set_ylabel('')
print(f'Modeled cell {"-".join(files[i].name.split("-")[2:-1])} at {100*alpha:4.1f}% contrast: noise_frac={noise_frac:3.1f}, nsegs={n}')
print(f'Modeled cell {"-".join(files[i].name.split("-")[:-3])} at {100*alpha:4.1f}% contrast: noise_frac={noise_frac:3.1f}, nsegs={n}')
print()
return alpha, noise_frac
def plot_chi2_data(ax, s, cell_name, run):
data_file = data_path / f'{cell_name}-baseline.npz'
data = np.load(data_file)
eodf = float(data['eodf'])
ratebase = float(data['ratebase/Hz'])
cvbase = float(data['cvbase'])
data_file = data_path / f'{cell_name}-spectral-100-s{run:02d}.npz'
data = np.load(data_file)
nsegs = data['nsegs']
fcutoff = data['fcutoff']
nfft = data['nfft']
deltat = data['deltat']
alpha = data['contrast']
freqs = data['freqs']
pss = data['pss']
prss = data['prss']
chi2 = np.abs(prss)*0.5/(pss.reshape(1, -1)*pss.reshape(-1, 1))
print(f'Measured cell {"-".join(data_file.name.split("-")[:-2])} at {100*alpha:4.1f}% contrast: r={ratebase:3.0f}Hz, CV={cvbase:4.2f}, dt={1000*deltat:4.2f}ms, nfft={nfft}, win={1000*deltat*nfft:6.1f}ms, nsegs={nsegs}')
print()
ax.text(1, 1.1, f'$N={nsegs}$',
ha='right', transform=ax.transAxes)
plot_chi2(ax, s, freqs, chi2, fcutoff, ratebase, 15)
return alpha, ratebase, eodf
def plot_ram(ax, contrast, eodf, wtime, wnoise):
tmax = 50
@ -217,7 +222,7 @@ if __name__ == '__main__':
s = plot_style()
fig, axs = plt.subplots(4, 4, cmsize=(s.plot_width, 0.85*s.plot_width),
width_ratios=[1, 0, 1, 1, 0.15, 0.9])
width_ratios=[1, 0, 1, 1, 0.2, 0.85])
fig.subplots_adjust(leftm=8, rightm=1.5, topm=3.5, bottomm=4,
wspace=0.25, hspace=0.6)
axs[0, 2].set_visible(False)

8
nonlinearbaseline.tex
View File

@ -471,7 +471,7 @@ Weakly nonlinear responses are expected in cells with sufficiently low intrinsic
\begin{figure*}[p]
\includegraphics[width=\columnwidth]{punitexamplecell}
\caption{\label{fig:punit} Linear and nonlinear stimulus encoding in example P-units. \figitem{A} Interspike interval (ISI) distribution of a cell's baseline activity, i.e. the cell is driven only by the unperturbed own electric field (cell identifier ``2020-10-27-ag''). This cell has a rather high baseline firing rate $r$ and an intermediate CV$_{\text{base}}$ of its interspike intervals, as indicated. \figitem{B} Power spectral density of the cell's baseline response with peaks at the cell's baseline firing rate $r$ and the fish's EOD frequency $f_{\text{EOD}}$. \figitem{C} Random amplitude modulation stimulus (top, red, with cutoff frequency of 300\,Hz) and evoked responses (spike raster, bottom) of the same P-unit for two different stimulus contrasts (right). The stimulus contrast quantifies the standard deviation of the AM relative to the fish's EOD amplitude. \figitem{D} Gain of the transfer function (first-order susceptibility), \eqnref{linearencoding_methods}, computed from the responses to 10\,\% (light blue) and 20\,\% contrast (dark blue) RAM stimulation of 5\,s duration. \figitem{E} Absolute value of the second-order susceptibility, \eqnref{eq:susceptibility}, for both the low and high stimulus contrast. At the lower stimulus contrast an anti-diagonal where the sum of the two stimulus frequencies equals the neuron's baseline frequency clearly sticks out of the noise floor. \figitem{F} At the higher contrast, the anti-diagonal is weaker. \figitem{G} Second-order susceptibilities projected onto the diagonal (means of all anti-diagonals of the matrices shown in \panel{E, F}). The anti-diagonals from \panel{E} and \panel{F} show up as a peak close to the cell's baseline firing rate $r$. The susceptibility index, SI($r$), quantifies the height of this peak relative to the values in the vicinity \notejb{See equation XXX}. \figitem{H} ISI distributions (top) and second-order susceptibilities (bottom) of four more example P-units (``2021-06-18-ae'', ``2012-03-30-ah'', ``2018-08-24-ak'', ``2018-08-14-ac'') covering the range of baseline firing rates and CV$_{\text{base}}$s as indicated. The first two cells show an anti-diagonal, but not the full expected triangular structure. The second-order susceptibilities of the other two cells are mostly flat and consequently the SI($r$) values are close to one.}
\caption{\label{fig:punit} Linear and nonlinear stimulus encoding in example P-units. \figitem{A} Interspike interval (ISI) distribution of a cell's baseline activity, i.e. the cell is driven only by the unperturbed own electric field (cell identifier ``2020-10-27-ag''). This cell has a rather high baseline firing rate $r=405$\,Hz and an intermediate CV$_{\text{base}}=0.49$ of its interspike intervals. \figitem{B} Power spectral density of the cell's baseline response with marked peaks at the cell's baseline firing rate $r$ and the fish's EOD frequency $f_{\text{EOD}}$. \figitem{C} Random amplitude modulation (RAM) stimulus (top, red, with cutoff frequency of 300\,Hz) and evoked responses (spike raster, bottom) of the same P-unit for two different stimulus contrasts (right). The stimulus contrast quantifies the standard deviation of the RAM relative to the fish's EOD amplitude. \figitem{D} Gain of the transfer function (first-order susceptibility), \eqnref{linearencoding_methods}, computed from the responses to 10\,\% (light blue) and 20\,\% contrast (dark blue) RAM stimulation of 5\,s duration. \figitem{E} Absolute value of the second-order susceptibility, \eqnref{eq:susceptibility}, for both the low and high stimulus contrast. At the lower stimulus contrast an anti-diagonal where the sum of the two stimulus frequencies equals the neuron's baseline frequency clearly sticks out of the noise floor. \figitem{F} At the higher contrast, the anti-diagonal is much weaker. \figitem{G} Second-order susceptibilities projected onto the diagonal (averages over all anti-diagonals of the matrices shown in \panel{E, F}). The anti-diagonals from \panel{E} and \panel{F} show up as a peak close to the cell's baseline firing rate $r$. The susceptibility index, SI($r$), quantifies the height of this peak relative to the values in the vicinity \notejb{See equation XXX}. \figitem{H} ISI distributions (top) and second-order susceptibilities (bottom) of four more example P-units (``2021-06-18-ae'', ``2017-07-18-ai'', ``2018-08-24-ak'', ``2018-08-14-ac'') covering the range of baseline firing rates and CV$_{\text{base}}$s as indicated. The first two cells show an anti-diagonal, but not the full expected triangular structure. The second-order susceptibilities of the other two cells are flat and consequently the SI($r$) values are close to one.}
\end{figure*}
Noise stimuli, here random amplitude modulations (RAM) of the EOD (\subfigref{fig:punit}{C}, top trace, red line), are commonly used to characterize stimulus-driven responses of sensory neurons using transfer functions (first-order susceptibility), spike-triggered averages, or stimulus-response coherences. Here, we additionally estimate the second-order susceptibility to quantify nonlinear encoding. P-unit spikes align more or less clearly with fluctuations in the RAM stimulus. A higher stimulus intensity, here a higher contrast of the RAM relative to the EOD amplitude (see methods), entrains the P-unit response more clearly (light and dark purple for low and high contrast stimuli, \subfigrefb{fig:punit}{C}). Linear encoding, quantified by the transfer function, \eqnref{linearencoding_methods}, is similar for the two RAM contrasts in this low-CV P-unit (\subfigrefb{fig:punit}{D}), as expected for a linear system. The first-order susceptibility is low for low frequencies, peaks in the range below 100\,Hz and then falls off again \citep{Benda2005}.
@ -487,7 +487,7 @@ In contrast, a high-CV P-unit (CV$_{\text{base}}=0.4$) does not exhibit pronounc
\begin{figure*}[p]
\includegraphics[width=\columnwidth]{ampullaryexamplecell}
\caption{\label{fig:ampullary} Linear and nonlinear stimulus encoding in example ampullary afferents. \figitem{A} Interspike interval (ISI) distribution of the cell's baseline activity (cell identifier ``2012-04-26-ae''). 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. Ampullary afferents do not respond to the fish's EOD frequency, $f_{\text{EOD}}$. \figitem{C} Band-limited white noise stimulus (top, red, with a cutoff frequency of 150\,Hz) added to the fish's self-generated electric field (no amplitude modulation) 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 5\,\% (light green) and 10\,\% 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. Both show a strong anti-diagonal where the two stimulus frequencies add up to the afferent's baseline firing rate. \figitem{G} Projections of the second-order susceptibilities in \panel{E, F} onto the diagonal. \figitem{H} ISI distributions (top) and second-order susceptibilites (bottom) of four more example afferents (``2010-11-26-an'', ``2011-10-25-ac'', ``2011-02-18-ab'', and ``2014-01-16-aj'').}
\caption{\label{fig:ampullary} Linear and nonlinear stimulus encoding in example ampullary afferents. \figitem{A} Interspike interval (ISI) distribution of the cell's baseline activity (cell identifier ``2012-05-15-ac''). 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. Ampullary afferents do not respond to the fish's EOD frequency, $f_{\text{EOD}}$ --- a sharp peak at $f_{\text{EOD}}$ is missing. \figitem{C} Band-limited white noise stimulus (top, red, with a cutoff frequency of 150\,Hz) added to the fish's self-generated electric field (no amplitude modulation!) 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 5\,\% (light green) and 10\,\% 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. Both show a clear anti-diagonal where the two stimulus frequencies add up to the afferent's baseline firing rate. \figitem{G} Projections of the second-order susceptibilities in \panel{E, F} onto the diagonal. \figitem{H} ISI distributions (top) and second-order susceptibilites (bottom) of four more example afferents (``2010-11-26-an'', ``2010-11-08-aa'', ``2011-02-18-ab'', and ``2014-01-16-aj'').}
\end{figure*}
Electric fish possess an additional electrosensory system, the passive or ampullary electrosensory system, that responds to low-frequency exogenous electric stimuli. The population of ampullary afferents is much less heterogeneous, and known for the much lower CVs of their baseline ISIs (CV$_{\text{base}}=0.06$ to $0.22$, \citealp{Grewe2017}). Ampullary cells do not phase-lock to the EOD and the ISIs are unimodally distributed (\subfigrefb{fig:ampullary}{A}). As a consequence of the high regularity of their baseline spiking activity, the corresponding power spectrum shows distinct peaks at the baseline firing rate \fbase{} and its harmonics. Since the cells do not respond to the self-generated EOD, there is no peak at \feod{} (\subfigrefb{fig:ampullary}{B}). When driven by a low-contrast noise stimulus (note: this is no longer an AM stimulus, \subfigref{fig:ampullary}{C}), ampullary cells exhibit very pronounced bands in the second-order susceptibility, where \fsum{} is equal to \fbase{} or its harmonic (yellow diagonals in \subfigrefb{fig:ampullary}{E}), implying strong nonlinear response components at these frequency combinations (\subfigrefb{fig:ampullary}{G}, top). With higher stimulus contrasts these bands disappear (\subfigrefb{fig:ampullary}{F}), the projection onto the diagonal loses its distinct peak at \fsum{} and its overall level is reduced (\subfigrefb{fig:ampullary}{G}, bottom).
@ -498,7 +498,7 @@ In the example recordings shown above (\figsrefb{fig:punit} and \fref{fig:ampull
\begin{figure*}[p]
\includegraphics[width=\columnwidth]{noisesplit}
\caption{\label{fig:noisesplit} Estimation of second-order susceptibilities. \figitem{A} \suscept{} (right) estimated from $N=198$ 256\,ms long FFT segments of an electrophysiological recording of another P-unit (cell ``2017-07-18-ai'', $r=78$\,Hz, CV$_{\text{base}}=0.22$) driven with a RAM stimulus with contrast 5\,\% (left). \figitem[i]{B} \textit{Standard condition} of model simulations with intrinsic noise (bottom) and a RAM stimulus (top). \figitem[ii]{B} \suscept{} estimated from simulations of the cell's LIF model counterpart (cell ``2017-07-18-ai'', table~\ref{modelparams}) based on a similar number of $N=100$ FFT segments. As in the electrophysiological recording only a weak anti-diagonal is visible. \figitem[iii]{B} Same as \panel[ii]{B} but using $10^6$ FFT segments. Now, the expected triangular structure is revealed. \figitem[iv]{B} Convergence of the \suscept{} estimate as a function of FFT segments. \figitem{C} At a lower stimulus contrast of 1\,\% the estimate did not converge yet even for $10^6$ FFT segments. \figitem[i]{D} Same as in \panel[i]{B} but in the \textit{noise split} condition: there is no external RAM signal (red) driving the model. Instead, a large part (90\,\%) of the total intrinsic noise is treated as a signal and is presented as an equivalent amplitude modulation ($s_{\xi}(t)$, orange), while the intrinsic noise is reduced to 10\,\% of its original strength (bottom, see methods for details). \figitem[i]{D} 100 FFT segments are still not sufficient for estimating \suscept{}. \figitem[iii]{D} Simulating one million segments reveals the full expected trangular structure of the second-order susceptibility. \figitem[iv]{D} In the noise-split condition, the \suscept{} estimate converges already at about $10^{4}$ FFT segments.}
\caption{\label{fig:noisesplit} Estimation of second-order susceptibilities. \figitem{A} \suscept{} (right) estimated from $N=198$ 256\,ms long FFT segments of an electrophysiological recording of another P-unit (cell ``2017-07-18-ai'', $r=78$\,Hz, CV$_{\text{base}}=0.22$) driven with a RAM stimulus with contrast 5\,\% (left). \figitem[i]{B} \textit{Standard condition} of model simulations with intrinsic noise (bottom) and a RAM stimulus (top). \figitem[ii]{B} \suscept{} estimated from simulations of the cell's LIF model counterpart (cell ``2017-07-18-ai'', table~\ref{modelparams}) based on the same RAM contrast and number of $N=100$ FFT segments. As in the electrophysiological recording only a weak anti-diagonal is visible. \figitem[iii]{B} Same as \panel[ii]{B} but using $10^6$ FFT segments. Now, the expected triangular structure is revealed. \figitem[iv]{B} Convergence of the \suscept{} estimate as a function of FFT segments. \figitem{C} At a lower stimulus contrast of 1\,\% the estimate did not converge yet even for $10^6$ FFT segments. The triangular structure is not revealed yet. \figitem[i]{D} Same as in \panel[i]{B} but in the \textit{noise split} condition: there is no external RAM signal (red) driving the model. Instead, a large part (90\,\%) of the total intrinsic noise is treated as a signal and is presented as an equivalent amplitude modulation ($s_{\xi}(t)$, orange, 10.6\,\% contrast), while the intrinsic noise is reduced to 10\,\% of its original strength (bottom, see methods for details). \figitem[i]{D} 100 FFT segments are still not sufficient for estimating \suscept{}. \figitem[iii]{D} Simulating one million segments reveals the full expected trangular structure of the second-order susceptibility. \figitem[iv]{D} In the noise-split condition, the \suscept{} estimate converges already at about $10^{4}$ FFT segments.}
\end{figure*}
%\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 electrophysiological 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}
@ -562,7 +562,7 @@ Overall, observing \nli{} values greater than at least 1.6, even for a number of
experimental data, the SI($r$) was estimated based on 100 FFT
segments \notejb{dt and nfft}. The SI($r$) is plotted against the
cells' CV of its baseline interspike intervals (left column), the
response modulation (the standard deviation of firing rate evoked
response modulation (standard deviation of firing rate evoked
by the band-limited white-noise stimulus) --- a measure of
effective stimulus strength (center column), and the CV of the
interspike intervals during stimulation with the white-noise