From 04578fd77e34afe7f6bf320b6d9b59ed8601cfc7 Mon Sep 17 00:00:00 2001 From: Jan Benda Date: Mon, 26 May 2025 11:36:02 +0200 Subject: [PATCH] updated captions to new figures --- dataoverview.py | 68 ++++++++++++++++++++++------------------ noisesplit.py | 73 +++++++++++++++++++++++-------------------- nonlinearbaseline.tex | 8 ++--- 3 files changed, 80 insertions(+), 69 deletions(-) diff --git a/dataoverview.py b/dataoverview.py index 4847908..284cf47 100644 --- a/dataoverview.py +++ b/dataoverview.py @@ -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 diff --git a/noisesplit.py b/noisesplit.py index af0e797..7e1c853 100644 --- a/noisesplit.py +++ b/noisesplit.py @@ -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) diff --git a/nonlinearbaseline.tex b/nonlinearbaseline.tex index 7d79bc0..eee0217 100644 --- a/nonlinearbaseline.tex +++ b/nonlinearbaseline.tex @@ -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