From eacf6ee04c2caee1c6628d1e9d96944606ce5cbb Mon Sep 17 00:00:00 2001 From: Jan Benda Date: Wed, 24 Sep 2025 16:01:36 +0200 Subject: [PATCH] submitted to eNeuro --- Makefile | 2 +- cover2.tex | 76 ++++++++ nonlinearbaseline.tex | 417 ++++++++++++++++++++++-------------------- 3 files changed, 297 insertions(+), 198 deletions(-) create mode 100644 cover2.tex diff --git a/Makefile b/Makefile index 6cec36f..2172ed2 100644 --- a/Makefile +++ b/Makefile @@ -1,7 +1,7 @@ TEXBASE=nonlinearbaseline BIBFILE=references.bib REBUTTALBASE= -COVERBASE=cover1 +COVERBASE=cover2 TEXFILE=$(TEXBASE).tex PDFFILE=$(TEXBASE).pdf diff --git a/cover2.tex b/cover2.tex new file mode 100644 index 0000000..07d114b --- /dev/null +++ b/cover2.tex @@ -0,0 +1,76 @@ +\documentclass[11pt]{article} + +\usepackage[utf8]{inputenc} +\usepackage{textcomp} +\usepackage{xcolor} +\usepackage{graphicx} + +\usepackage[ngerman,english]{babel} + +\usepackage[left=25mm, right=25mm, top=20mm, bottom=20mm]{geometry} +\setlength{\parskip}{2ex} + +\usepackage[mediumqspace,Gray,squaren]{SIunits} % \ohm, \micro + +\usepackage{natbib} +%\bibliographystyle{jneurosci} + +\usepackage[breaklinks=true,bookmarks=true,bookmarksopen=true,pdfpagemode=UseNone,pdfstartview=FitH,colorlinks=true,citecolor=blue,urlcolor=blue]{hyperref} + +\setlength{\parindent}{0em} + + +\begin{document} +\hspace*{\fill} September 24, 2025 + +\noindent Dear Arvind Kumar, + +we would like to submit our manuscript ``Spike generation in +electroreceptor afferents introduces additional spectral response +components by weakly nonlinear interactions'' for publication in +eNeuro. It is the result of a collaborative effort on the theoretical +side, Benjamin Linder, HU Berlin, and the experimental and numerical +modeling side, Jan Grewe and myself at the University of Tuebingen, +within the DFG priority program 2205 ``Evolutionary optimization of +neuronal processing''. + +Computational neuroscientists have a strong interest in characterizing +non-linearities and study their functional consequences. However, +experimental backing of these theoretical findings are scarce. For +example, encoding of dynamic stimuli in spike trains often can be well +approximated by linear response theory, in particular when driving the +neuron in the supra-threshold regime at low signal-to-noise +ratios. This has been exploited in numerous theoretical studies. At +less noise or stronger stimulus amplitudes non-linear effects become +more prominent. In the weakly non-linear regime, the second term of +the Volterra series becomes relevant. Benjamin Lindner developed +analytical solutions of this term for the leaky integrate-and-fire +neuron, which predicts non-linear interactions whenever one or the sum +of two stimulus frequencies matches the neuron's baseline firing +rate. Until now, however, these fundamental non-linearities arising +from the core mechanism of spike generation have not been reported in +any experimental data. + +With our work we set out to fill this gap. We scan a large set of +electrophysiological data measured in two types of electrosensory +neurons of the electric fish \textit{Apteronotus leptorhynchus} for +signatures of these non-linearities. In ampullary cells, non-linear +interaction between two stimulus frequencies are prominent, whereas in +P-units they are harder to find. Estimating the second-order +susceptibilites from real data turns out to be a hard problem as +limited data leads to poor estimates. Comparison with models that have +been fitted to individual P-units, we are able to deduce the presence +of non-linear interactions also in some of the P-units. Finally, we +discuss our findings and the relevance of non-linear interactions in +the neuroethological context. + +We believe that this analysis of electrophysiological data close to +expectations from theoretical work is of strong interest to many +readers of eNeuro and can inspire future research in other sensory +systems. + +Best regards,\\%[-2ex] +%\hspace*{0.17\textwidth}\includegraphics[width=0.3\textwidth]{JanBenda-Signature2020}\\ +Prof. Dr. Jan Benda, on behalf of all authors + +\end{document} diff --git a/nonlinearbaseline.tex b/nonlinearbaseline.tex index b6ff7f6..426c4b6 100644 --- a/nonlinearbaseline.tex +++ b/nonlinearbaseline.tex @@ -11,10 +11,10 @@ Jan Benda\textsuperscript{1, 4, $\dagger$}} \date{\normalsize - \textsuperscript{1} Institute for Neurobiology, Eberhard Karls Universit\"at T\"ubingen, Germany\\ - \textsuperscript{2} Bernstein Center for Computational Neuroscience Berlin, Berlin, Germany\\ - \textsuperscript{3} Department of Physics, Humboldt University Berlin, Berlin, Germany\\ - \textsuperscript{4} Bernstein Center for Computational Neuroscience Tübingen, Germany\\ + \textsuperscript{1} Institute for Neurobiology, Eberhard Karls Universit\"at T\"ubingen,\\ Auf der Morgenstelle 28, 72076 T\"ubingen, Germany\\ + \textsuperscript{2} Bernstein Center for Computational Neuroscience Berlin, 10115 Berlin, Germany\\ + \textsuperscript{3} Department of Physics, Humboldt University Berlin, 12489 Berlin, Germany\\ + \textsuperscript{4} Bernstein Center for Computational Neuroscience Tübingen, 72076 T\"ubingen, Germany\\ \textsuperscript{$\dagger$} corresponding author: \url{jan.benda@uni-tuebingen.de}} \newcommand{\runningtitle}{} @@ -22,7 +22,7 @@ %%%%% overall style %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % number the lines: -\newcommand{\setlinenumbers}{false} +\newcommand{\setlinenumbers}{true} % use double spacing: \newcommand{\setdoublespacing}{false} @@ -236,10 +236,27 @@ \maketitle -%\paragraph{Short title:} -%\paragraph{Corresponding author:}Jan Grewe, E-mail: jan.grewe@uni-tuebingen.de -%\paragraph{Conflict of interest:}The authors declare no conflict of interest. +\thispagestyle{empty} + +\noindent +\begin{tabular}{@{}lr} + Number of figures & 10 \\ + Number of tables & 0 \\ + Number of multimedia & 0 \\ + Number of words in abstract & 194 \\ + Number of words in introduction & 659 \\ + Number of words in discussion & 1891 +\end{tabular} + +\paragraph{Acknowledgements} +We thank Tim Hladnik, Henriette Walz, Franziska Kuempfbeck, Fabian Sinz, Laura Seidler, Eva Vennemann, and Ibrahim Tunc for the data they recorded over the years in our lab. + +\paragraph{Conflict of interest:}The authors declare no conflict of interest. + %\paragraph{Author Contributions:} All authors designed the study and discussed the results. AB performed the data analyses and modeling, AB and JG drafted the paper, all authors discussed and revised the manuscript. + +\paragraph{Founding sources} +Supported by SPP 2205 ``Evolutionary optimisation of neuronal processing'' by the DFG, project number 430157666. \begin{keywords} \item Volterra series @@ -249,18 +266,25 @@ \end{keywords} -% Please keep the abstract below 300 words +\newpage +\setcounter{page}{1} + +\begin{center} + \sffamily\bfseries\LARGE Spike generation in electroreceptor afferents + introduces\\ additional spectral response components by\\ weakly + nonlinear interactions +\end{center} + +% 250 words \section{Abstract} -%Neuronal processing is inherently nonlinear --- spiking thresholds or rectification at synapses are essential to neuronal computations.} Spiking thresholds in neurons or rectification at synapses are essential for neuronal computations rendering neuronal processing inherently nonlinear. Nevertheless, linear response theory has been instrumental for understanding, for example, the impact of noise or neuronal synchrony on signal transmission, or the emergence of oscillatory activity, but is valid only at low stimulus amplitudes or large levels of intrinsic noise. At higher signal-to-noise ratios, however, nonlinear response components become relevant. Theoretical results for leaky integrate-and-fire neurons in the weakly nonlinear regime suggest strong responses at the sum of two input frequencies if these frequencies or their sum match the neuron's baseline firing rate. We here analyze nonlinear responses in two types of primary electroreceptor afferents, the P-units of the active and the ampullary cells of the passive electrosensory system of the wave-type electric fish \textit{Apteronotus leptorhynchus}. In our combined experimental and modeling approach we identify these predicted nonlinear responses in low-noise P-units and, much stronger, in ampullary cells. Our results provide experimental evidence for nonlinear responses of spike generators in the weakly nonlinear regime. We conclude that such nonlinear responses occur in any sensory neuron that operates in similar regimes particularly at near-threshold stimulus conditions. -% Please keep the Author Summary between 150 and 200 words -% Use the first person. -\section{Author summary} -The generation of action potentials involves a strong threshold nonlinearity. Nevertheless, the encoding of stimuli with small amplitudes by neurons with sufficient intrinsic noise can be well described as a linear system. As the stimulus amplitude is increased, new spectral components start to appear in the so called weakly nonlinear regime. This regime has been well characterized by theoretical analysis based on simple neuron models such as the leaky integrate-and-fire model. These findings predict nonlinear interactions whenever one or the sum of two stimulus frequencies matches the neuron's baseline firing rate. We set out to find these signatures in a large set of electrophysiological recordings from electroreceptive neurons of a weakly electric fish. In ampullary cells these interactions are prominent, whereas in P-units they are harder to find. Estimating nonlinear response kernels from limited real data turns out to be a hard problem. By comparison with models that have been fitted to individual P-units we are then able to interpret the poor estimates. The nonlinear response components could boost sensory responses to weak signals emitted, for example, by distant conspecifics. +% max 120 words +\section{Significance statement} +The generation of action potentials involves a strong threshold nonlinearity. Nevertheless, the encoding of stimuli with small amplitudes by neurons with sufficient intrinsic noise can be well described as a linear system. As the stimulus amplitude is increased, new spectral components start to appear in the so called weakly nonlinear regime. Theory predicts nonlinear interactions whenever one or the sum of two stimulus frequencies matches the neuron's baseline firing rate. Indeed, we find these interactions in a large set of electrophysiological recordings from primary electroreceptive afferents of a weakly electric fish. The nonlinear response components could boost sensory responses to weak signals emitted, for example, by distant conspecifics. + - \section{Introduction} \begin{figure*}[t] @@ -276,6 +300,187 @@ Noise in nonlinear systems, however, linearizes the system's response properties Here we search for such weakly nonlinear responses in electroreceptors of the two electrosensory systems of the wave-type electric fish \textit{Apteronotus leptorhynchus}, i.e. the tuberous (active) and the ampullary (passive) electrosensory system. The p-type electroreceptor afferents of the active system (P-units) are driven by the fish's high-frequency, quasi-sinusoidal electric organ discharges (EOD) and encode disturbances of it \citep{Bastian1981a}. The electroreceptors of the passive system are tuned to lower-frequency exogeneous electric fields such as caused by muscle activity of prey \citep{Kalmijn1974}. As different animals have different EOD-frequencies, being exposed to stimuli of multiple distinct frequencies is part of the animal's everyday life \citep{Benda2020,Henninger2020} and weakly nonlinear interactions may occur in the electrosensory periphery. In communication contexts \citep{Walz2014, Henninger2018} the EODs of interacting fish superimpose and lead to periodic amplitude modulations (AMs or beats) of the receiver's EOD. Nonlinear mechanisms in P-units, enable encoding of AMs in their time-dependent firing rates \citep{Bastian1981a, Walz2014, Middleton2006, Barayeu2023}. When multiple animals interact, the EOD interferences induce second-order amplitude modulations referred to as envelopes \citep{Yu2005, Fotowat2013, Stamper2012Envelope} and saturation nonlinearities allow also for the encoding of these in the electrosensory periphery \citep{Savard2011}. Field observations have shown that courting males were able to react to the extremely weak signals of distant intruding males despite the strong foreground EOD of the nearby female \citep{Henninger2018}. Weakly nonlinear interactions at particular combinations of signals can be of immediate relevance in such settings as they could boost detectability of the faint signals \citep{Schlungbaum2023}. + +\section{Materials and Methods} + +\subsection{Experimental subjects and procedures} + +Within this project, we re-evaluated datasets that were recorded between 2010 and 2023 at the Ludwig Maximilian University (LMU) M\"unchen and the Eberhard-Karls University T\"ubingen. All experimental protocols complied with national and European law and were approved by the respective Ethics Committees of the Ludwig-Maximilians Universität München (permit no. 55.2-1-54-2531-135-09) and the Eberhard-Karls Unversität Tübingen (permit no. ZP 1/13 and ZP 1/16). +The final sample consisted of 172 P-units and 30 ampullary electroreceptor afferents recorded in 80 weakly electric fish of both sexes of the species \lepto{}. Fish were obtained from a commercial supplier for tropical fish (Aquarium Glaser GmbH, Rodgau, Germany) and kept in tanks with a water temperature of $25\pm1\,^\circ$C and a conductivity of around $270\,\micro\siemens\per\centi\meter$ under a 12\,h:12\,h light-dark cycle. + +Before surgery, the animals were deeply anesthetized via bath application of a solution of MS222 (120\,mg/l, PharmaQ, Fordingbridge, UK) buffered with Sodium Bicarbonate (120\,mg/l). The posterior anterior lateral line nerve (pALLN) was exposed by making a small cut into the skin covering the nerve. The cut was placed dorsal of the operculum just before the nerve descends towards the anterior lateral line ganglion (ALLNG). Those parts of the skin that were to be cut were locally anesthetized by cutaneous application of liquid lidocaine hydrochloride (20\,mg/ml, bela-pharm GmbH). During the surgery, water supply was ensured by a mouthpiece to maintain anesthesia with a solution of MS222 (100\,mg/l) buffered with Sodium Bicarbonate (100\,mg/l). After surgery, fish were immobilized by intramuscular injection of from 25\,$\micro$l to 50\,$\micro$l of tubocurarine (5\,mg/ml dissolved in fish saline; Sigma-Aldrich). +Respiration was then switched to normal tank water and the fish was transferred to the experimental tank. + +\subsection{Electrophysiological recordings} +For the recordings fish were positioned centrally in the experimental tank, with the major parts of their body submerged into the water. Those body parts that were above the water surface were covered with paper tissue to avoid drying of the skin. Local analgesia was refreshed in intervals of two hours by cutaneous application of Lidocaine (2\,\%; bela-pharm, Vechta, Germany) around the surgical wounds. Electrodes (borosilicate; 1.5\,mm outer diameter; GB150F-8P; Science Products, Hofheim, Germany) were pulled to a resistance of 50--100\,\mega\ohm{} (model P-97; Sutter Instrument, Novato, CA) and filled with 1\,M KCl solution. Electrodes were fixed in a microdrive (Luigs-Neumann, Ratingen, Germany) and lowered into the nerve. Recordings of electroreceptor afferents were amplified and lowpass filtered at 10\,kHz (SEC-05, npi-electronics, Tamm, Germany, operated in bridge mode). All signals, neuronal recordings, recorded EOD, and the generated stimulus, were digitized with sampling rates of 20 or 40\,kHz (PCI-6229, National Instruments, Austin, TX). RELACS (\url{https://github.com/relacs/relacs}) running on a Linux computer was used for online spike and EOD detection, stimulus generation, and calibration. Recorded data was then stored on the hard drive for offline analysis. + +\subsection{Identification of P-units and ampullary cells} +Recordings were classified as P-units if baseline action potentials phase locked to the EOD with vectors strengths between 0.7 and 0.95, a baseline firing rate larger than 30\,Hz, a serial correlation of subsequent interspike intervals below zero, a coefficient of variation of baseline interspike intervals below 1.5 und during stimulation below 2. As ampullary cells we classified recordings with vector strengths below 0.15, baseline firing rate above 10\,Hz, baseline CV below 0.18, CV during stimulation below 1.0, and a response modulation during stimulation below 80\,Hz \citep{Grewe2017}. We here selected only those cells of which the neuron's baseline activity as well as the responses to band-limited white noise stimuli were recorded. + +\subsection{Electric field recordings} +For monitoring the EOD without the stimulus, two vertical carbon rods ($11\,\centi\meter$ long, 8\,mm diameter) in a head-tail configuration were placed isopotential to the stimulus. Their signal was differentially amplified with a gain factor between 100 and 500 (depending on the recorded animal) and band-pass filtered (3 to 1500\,Hz pass-band, DPA2-FX; npi electronics, Tamm, Germany). For an estimate of the transdermal potential that drives the electroreceptors, two silver wires spaced by 1\,cm were located next to the left gill of the fish and orthogonal to the fish's longitudinal body axis (amplification 100 to 500 times, band-pass filtered with 3 to 1\,500\,Hz pass-band, DPA2-FX; npi-electronics, Tamm, Germany). This local EOD measurement recorded the combination of the fish's own EOD and the applied stimulus. + + +\subsection{Stimulation}\label{rammethods} +Electric stimuli were attenuated (ATN-01M, npi-electronics, Tamm, Germany), isolated from ground (ISO-02V, npi-electronics, Tamm, Germany) and delivered via two horizontal carbon rods (30 cm length, 8 mm diameter) located $15\,\centi\meter$ laterally to the fish. +The fish were stimulated with band-limited white noise stimuli with a cut-off frequency of 150\,Hz (ampullary), 300 or 400\,Hz (P-units). The stimulus intensity is given as a contrast, i.e. the standard deviation of the noise stimulus relative to the fish's EOD amplitude. The contrast varied between 1 and 20\,\% (median 5\,\%) for P-units and 2.5 and 20\,\% (median 5\,\%) for ampullary cells. Only noise stimuli with a duration of at least 2\,s (maximum of 50\,s, median 10\,s) and enough repetitions to results in at least 100 FFT segments (see below, P-units: 100--1520, median 313, ampullary cells: 105 -- 3648, median 722) were included into the analysis. When ampullary cells were recorded, the white noise was directly applied as the stimulus. To create random amplitude modulations (RAM) for P-unit recordings, the noise stimulus was first multiplied with the EOD of the fish (MXS-01M; npi electronics). + +\subsection{Data analysis} Data analysis was done in Python 3 using the packages matplotlib \citep{Hunter2007}, numpy \citep{Walt2011}, scipy \citep{scipy2020}, pandas \citep{Mckinney2010}, nixio \citep{Stoewer2014}, and thunderlab (\url{https://github.com/bendalab/thunderlab}). + +\paragraph{Code accessibility} +The P-unit model parameters and spectral analysis algorithms are available at \url{https://github.com/bendalab/punitmodel/tree/v1}. + +\paragraph{Baseline analysis}\label{baselinemethods} +The baseline firing rate $r$ was calculated as the number of spikes divided by the duration of the baseline recording (median 32\,s). The coefficient of variation (CV) of the interspike intervals (ISI) is their standard deviation relative to their mean: $\rm{CV}_{\rm base} = \sqrt{\langle (ISI- \langle ISI \rangle) ^2 \rangle} / \langle ISI \rangle$. If the baseline was recorded several times in a recording, the measures from the longest recording were taken. + +\paragraph{White noise analysis} \label{response_modulation} +When stimulated with band-limited white noise stimuli, neuronal activity is modulated around the average firing rate that is similar to the baseline firing rate and in that way encodes the time-course of the stimulus. For an estimate of the time-dependent firing rate $r(t)$ we convolved each spike train with normalized Gaussian kernels with standard deviation of 2\,ms, if not mentioned otherwise, and averaged the resulting single-trail firing rates over trials. The response modulation quantifies the variation of $r(t)$ computed as the standard deviation in time $\sigma_{s} = \sqrt{\langle (r(t)-\langle r(t) \rangle_t )^2\rangle_t}$, where $\langle \cdot \rangle_t$ denotes averaging over time. + +\paragraph{Spectral analysis}\label{susceptibility_methods} +To characterize the relation between the spiking response evoked by white-noise stimuli, we estimated the first- and second-order susceptibilities in the frequency domain. For this we converted spike times into binary vectors $x_k$ with $\Delta t = 0.5$\,ms wide bins that are set to 2\,kHz where a spike occurred and zero otherwise. Fast Fourier transforms (FFT) $S(\omega)$ and $X(\omega)$ of the stimulus $s_k$ (also down-sampled to a sampling rate of 2\,kHz) and $x_k$, respectively, were computed numerically according to +\begin{equation} + \label{fourier} + X(\omega) = \sum_{k=0}^{N-1} \, x_k e^{- i \omega k} +\end{equation} +for $N=512$ long segments of $T=N \Delta t = 256$\,ms duration with no overlap, resulting in a spectral resolution of about 4\,Hz. Note, that for a real Fourier integral a factor $\Delta t$ is missing. For simplicity we use angular frequencies $\omega=2\pi f$ instead of frequencies $f$. + +In the experimental data the duration of the noise stimuli varied and they were presented once or repeatedly (frozen noise). For the analysis we discarded the responses within the initial 200\,ms of stimulation in each trial. To make the recordings comparable we always used the first 100 segments from as many trials as needed for the following analysis. + +In the simulations we generated for each trial a new realization of the noise stimulus. We discarded the first 500\,ms of the response and used the following 10 FFT segments for the analysis. This was repeated until $10^6$ or $10^7$ FFT segments were collected. + +The power spectrum of the stimulus $s(t)$ was estimated as +\begin{equation} + \label{powerss} + S_{ss}(\omega) = \frac{\Delta t}{N} \langle S(\omega) S^*(\omega) \rangle +\end{equation} +with $S^*$ being the complex conjugate of $S$ and $\langle \cdot \rangle$ denoting the average over the FFT segments. The factor in front is $\Delta t^2$ from the missing integration factors of the two Fourier transforms, \eqnref{fourier}, divided by $T=N \Delta t$, needed to make this a proper power spectral density. + +The power spectrum of the spike trains $S_{xx}(\omega)$ was estimated accordingly. Likewise, the cross-spectrum $S_{xs}(\omega)$ between stimulus and evoked spike trains was estimated according to +\begin{equation} + \label{crossxs} + S_{xs}(\omega) = \frac{\Delta t}{N} \langle X(\omega) S^*(\omega) \rangle +\end{equation} + +The first-order susceptibility (transfer function) +\begin{equation} + \label{linearencoding_methods} + \chi_{1}(\omega) = \frac{S_{xs}(\omega) }{S_{ss}(\omega) } +\end{equation} +was then computed from $S_{xs}(\omega)$ and $ S_{ss}(\omega)$. We report $\chi_{1}(\omega)$ in Hz/\%, i.e. firing rate per percent stimulus contrast. Multiplying $\chi_{1}(\omega)$ with the contrast of a sinusoidal stimulus in percent results in the amplitude of the evoked firing rate modulation in Hertz. + +The second-order cross-spectrum +\begin{equation} + \label{crossxss} + S_{xss}(\omega_{1},\omega_{2}) = \frac{\Delta t^2}{N} \langle X(\omega_{1} + \omega_{2}) S^*(\omega_{1}) S^*(\omega_{2}) +\end{equation} +quantifies nonlinear interactions that generate responses at the sum and difference (for negative $\omega_1$ or $\omega_2$) evoked by two stimulus frequencies $\omega_1$ and $\omega_2$. +The second-order susceptibility +\begin{equation} + \label{chi2} + \chi_{2}(\omega_{1}, \omega_{2}) = \frac{S_{xss} (\omega_{1},\omega_{2})}{2 S_{ss} (\omega_{1}) S_{ss} (\omega_{2})} +\end{equation} +normalizes the second-order cross-spectrum by the spectral power at the two stimulus frequencies. We report $\chi_{2}(\omega_{1}, \omega_{2})$ in Hz/\%$^2$, i.e. firing rate per percent stimulus contrast squared. +Throughout the manuscript we only show the absolute values of the complex-valued second-order susceptibility matrix and ignore the corresponding phases. + +\paragraph{Susceptibility index}\label{projected_method} +We expected to see a sharp ridge in the second-order susceptibility at $\omega_1 + \omega_2 = r$ \citep{Voronenko2017,Franzen2023}. To characterize this in a single number we computed a susceptibility index. First, we projected the absolute values of the second-order susceptibility matrix onto the diagonal by averaging over anti-diagonal elements. In this projection $D(f)$ we took the position of the maximum +\begin{equation} + \label{dmax} + f_{\rm peak} = {\rm argmax} D(r - 50\,\rm{Hz} \leq f \leq r + 50\,\rm{Hz}) +\end{equation} +within $\pm 50$\,Hz of the neuron's baseline firing rate $r$ as the position of the expected peak. For an estimate of the noise-floor surrounding this peak we averaged over 10\,Hz wide windows 10\,Hz to the left and right of the peak: +\begin{equation} + \label{dref} + D_{\rm ref} = \frac{1}{2}\left(\langle D(f_{\rm peak}-20\,\rm{Hz} \leq f \leq f_{\rm peak}-10\,\rm{Hz})\rangle_f + \langle D(f_{\rm peak}+10\,\rm{Hz} \leq f \leq f_{\rm peak}+20\,\rm{Hz}))\rangle_f \right) +\end{equation} +The size of the peak relative to this reference is then the susceptibility index +\begin{equation} + \label{siindex} + SI = \frac{D(f_{\rm peak})}{D_{\rm ref}} +\end{equation} +Values larger than one indicate a sharp ridge in the susceptibility matrix close to where the stimulus frequencies add up to the baseline firing rate. + + +\begin{figure*}[t] + \includegraphics[width=\columnwidth]{flowchart.pdf} + \caption{\label{flowchart} + Architecture of the P-unit model. Each row illustrates subsequent processing steps for three different stimulation regimes: (i) baseline activity without external stimulus, only the fish's self-generated EOD (the carrier, \eqnref{eq:eod}) is present. (ii) RAM stimulation, \eqnref{eq:ram_equation}. The amplitude of the EOD carrier is modulated with a weak (2\,\% contrast) band-limited white-noise stimulus. (iii) Noise split, \eqnsref{eq:ram_split}--\eqref{eq:Noise_split_intrinsic}, where 90\,\% of the intrinsic noise is replaced by a RAM stimulus, whose amplitude is scaled to maintain the mean firing rate and the CV of the ISIs of the model's baseline activity. As an example, simulations of the model for cell ``2012-07-03-ak'' are shown. \figitem{A} The stimuli are thresholded, \eqnref{eq:threshold2}, by setting all negative values to zero. \figitem{B} Subsequent low-pass filtering, \eqnref{eq:dendrite}, attenuates the carrier and carves out the AM signal. \figitem{C} Intrinsic Gaussian white-noise is added to the signals shown in \panel{B}. Note the reduced internal noise amplitude in the noise split (iii) condition. \figitem{D} Spiking output of the LIF model, \eqnsref{eq:LIF}--\eqref{spikethresh}, in response to the sum of \panel{B} and \panel{C}. \figitem{E} Power spectra of the LIF neuron's spiking activity. Both, baseline activity (\panel[i]{E}) and noise split (\panel[iii]{E}), have the same peaks in the response spectrum at $r$, $f_{EOD} - r$, $f_{EOD}$, and $f_{EOD} + r$. With RAM stimulation (\panel[ii]{E}), the peak at the baseline firing rate, $r$, is washed out.} +\end{figure*} + +\subsection{Leaky integrate-and-fire models for P-units}\label{lifmethods} + +Modified leaky integrate-and-fire (LIF) models were constructed to reproduce the specific firing properties of P-units \citep{Chacron2001, Sinz2020, Barayeu2023}. Its basic components (static nonlinearity, low-pass filtering and spike generation) are equivalent to models of hair cells in auditory systems \citep{Eggermont1983}. The sole driving input into the P-unit model during baseline, i.e. when no external stimulus was given, is the fish's own EOD, modeled as a cosine wave +\begin{equation} + \label{eq:eod} + y(t) = y_{EOD}(t) = \cos(2\pi f_{EOD} t) + \end{equation} +with EOD frequency $f_{EOD}$ and an amplitude of one. + +To mimic the interaction with other fish, the EODs of a second or third fish with EOD frequencies $f_1$ and $f_2$, respectively, were added to the normalized EOD, \eqnref{eq:eod}, of the receiving fish according to their contrasts, $c_1$ and $c_2$ at the position of the receiving fish: +\begin{equation} + \label{eq:modelbeats} + y(t) = \cos(2\pi f_{EOD} t) + c_1 \cos(2\pi f_1 t) + c_2\cos(2\pi f_2 t) +\end{equation} +For two fish, $c_2 = 0$. + +Random amplitude modulations (RAMs) were simulated by first generating the AM as a band-limited white noise stimulus $s(t)$. For this, random real and imaginary numbers were drawn from Gaussian distributions for each frequency component in the range from 0 to 300\,Hz in the Fourier domain \citep{Billah1990,Skorjanc2023}. By means of the inverse Fourier transform, the time course of the RAM stimulus, $s(t)$, was generated. The input to the model was then +\begin{equation} + \label{eq:ram_equation} + y(t) = (1+ s(t)) \cos(2\pi f_{EOD} t) +\end{equation} +The contrast $c$ of the RAM is the standard deviation of the RAM relative to the amplitude of the receiving fish. + +First, the input $y(t)$ is thresholded by setting negative values to zero: +\begin{equation} + \label{eq:threshold2} + \lfloor y(t) \rfloor_0 = \left\{ \begin{array}{rcl} y(t) & ; & y(t) \ge 0 \\ 0 & ; & y(t) < 0 \end{array} \right. +\end{equation} +(\subfigrefb{flowchart}{A}). This thresholds models the transfer function of the synapses between the primary receptor cells and the afferent. Together with a low-pass filter +\begin{equation} + \label{eq:dendrite} + \tau_{d} \frac{d V_{d}}{d t} = -V_{d}+ \lfloor y(t) \rfloor_{0} +\end{equation} +the threshold operation is required for extracting the amplitude modulation from the input \citep{Barayeu2024}. The low-pass filter models passive signal conduction in the afferent's dendrite (\subfigrefb{flowchart}{B}) and $\tau_{d}$ is the membrane time constant of the dendrite. Dendritic low-pass filtering was also necessary to reproduce the loose coupling of P-unit spikes to the EOD while maintaining high sensitivity at small amplitude modulations. + +The dendritic voltage $V_d(t)$ is then fed into a stochastic leaky integrate-and-fire (LIF) model with adaptation, +\begin{equation} + \label{eq:LIF} + \tau_{m} \frac{d V_{m}}{d t} = - V_{m} + \mu + \beta V_{d} - A + \sqrt{2D}\;\xi(t) +\end{equation} +where $\tau_{m}$ is the membrane time-constant, $\mu$ is a fixed bias current, $\beta$ is a scaling factor for $V_{d}$, $A$ is an inhibiting adaptation current, and $\sqrt{2D}\;\xi(t)$ is a Gaussian white noise with strength $D$. Note, that all state variables, membrane voltages $V_d$ and $V_m$ as well as the adaptation current $A$, are dimensionless. + +The adaptation current $A$ follows +\begin{equation} + \label{eq:adaptation} + \tau_{A} \frac{d A}{d t} = - A +\end{equation} +with adaptation time constant $\tau_A$ \citep{Benda2003,Benda2005,Benda2010}. + +Whenever the membrane voltage $V_m(t)$ crosses the spiking threshold $\theta=1$, a spike was generated, $V_{m}(t)$ was reset to $0$, the adaptation current was incremented by $\Delta A/\tau_A$, and integration of $V_m(t)$ was paused for the duration of a refractory period $t_{ref}$ (\subfigrefb{flowchart}{D}): +\begin{equation} + \label{spikethresh} + V_m(t) \ge \theta \; : \left\{ \begin{array}{rcl} V_m & \mapsto & 0 \\ A & \mapsto & A + \Delta A/\tau_A \end{array} \right. +\end{equation} + +The P-unit models were integrated by the Euler forward method with a time-step of $\Delta t = 0.05$\,ms. For each trial of a simulation, $V_{m}$ was drawn from a uniform distribution between 0 and 1 and the initial value of $A$ was jittered by adding a random number drawn from a normal distribution with standard deviation of 2\,\% of its initial value. Then the first 500\,ms of any simulation were discarded to remove remaining transients. + +The eight free parameters of the P-unit model $\beta$, $\tau_m$, $\mu$, $D$, $\tau_A$, $\Delta_A$, $\tau_d$, and $t_{ref}$, were fitted to both the baseline activity (baseline firing rate, CV of ISIs, serial correlation of ISIs at lag one, and vector strength of spike coupling to EOD) and the responses to step increases and decreases in EOD amplitude (onset and steady-state responses, effective adaptation time constant, \citealp{Benda2005}) of recorded P-units. Model parameters of all 39 cells are summarized in file \texttt{models.csv} of our \texttt{punitmodel} repository at \url{https://github.com/bendalab/punitmodel/tree/v1}. + + +\subsection{Noise split} +\label{intrinsicsplit_methods} +Based on the Furutsu-Novikov theorem \citep{Furutsu1963,Novikov1965,Lindner2022,Egerland2020}, we split the total noise, $\sqrt{2D}\;\xi(t)$, of a LIF model, \eqnref{eq:LIF}, into two parts. The first part is the intrinsic noise term, $\sqrt{2D \, \alpha_{\rm noise}}\;\xi(t)$, whose strength is reduced by a factor $\alpha_{\rm noise}=0.1$ (\subfigrefb{flowchart}\,\panel[iii]{C}). The second part replaces the now missing intrinsic noise by a driving input signal $s_{\xi}(t)$, a RAM stimulus with frequencies up to 300\,Hz (\subfigrefb{flowchart}\,\panel[iii]{A}). The LIF model with noise split then reads +\begin{eqnarray} + \label{eq:ram_split} + y(t) & = & (1+ s_\xi(t)) \cos(2\pi f_{EOD} t) \\ + \label{eq:Noise_split_intrinsic_dendrite} + \tau_{d} \frac{d V_{d}}{d t} & = & -V_{d}+ \lfloor y(t) \rfloor_{0} \\ + \label{eq:Noise_split_intrinsic} + \tau_{m} \frac{d V_{m}}{d t} & = & - V_{m} + \mu + \beta V_{d} - A + \sqrt{2D \, \alpha_{\rm noise}}\;\xi(t) +\end{eqnarray} +Both, the reduced intrinsic noise and the RAM stimulus, need to replace the original intrinsic noise. Because the RAM stimulus is band-limited and undergoes some transformations before it is added to the reduced intrinsic noise, it is not \textit{a priori} clear, what the amplitude of the RAM should be. We bisected the amplitude of $s_\xi(t)$, until the CV of the resulting interspike intervals matched the one of the original model's baseline activity. The second-order cross-spectra, \eqnref{crossxss}, were computed between the RAM stimulus $s_{\xi}(t)$ and the spike train $x(t)$ it evoked. In this way, the effective signal-to-noise ratio can be increased while maintaining the total noise in the system. + + \section{Results} \begin{figure*}[t] @@ -495,188 +700,6 @@ The weakly nonlinear interactions in low-CV P-units could facilitate the detecta We have demonstrated pronounced nonlinear responses in primary electrosensory afferents at weak stimulus amplitudes and sufficiently low intrinsic noise levels. The observed nonlinearities match the expectations from previous theoretical studies \citep{Voronenko2017,Franzen2023}. The resulting nonlinear components introduce spectral components not present in the original stimulus, that may provide an edge in the context of signal detection problems at stimulus amplitudes close to threshold \citep{Schlungbaum2023}. Electrosensory afferents share an evolutionary history with hair cells \citep{Baker2019} and share many response properties with mammalian auditory nerve fibers \citep{Barayeu2023, Joris2004}. Thus, we expect weakly nonlinear responses for near-threshold stimulation in auditory nerve fibers as well. -\section{Methods} - -\subsection{Experimental subjects and procedures} - -Within this project, we re-evaluated datasets that were recorded between 2010 and 2023 at the Ludwig Maximilian University (LMU) M\"unchen and the Eberhard-Karls University T\"ubingen. All experimental protocols complied with national and European law and were approved by the respective Ethics Committees of the Ludwig-Maximilians Universität München (permit no. 55.2-1-54-2531-135-09) and the Eberhard-Karls Unversität Tübingen (permit no. ZP 1/13 and ZP 1/16). -The final sample consisted of 172 P-units and 30 ampullary electroreceptor afferents recorded in 80 weakly electric fish of both sexes of the species \lepto{}. Fish were obtained from a commercial supplier for tropical fish (Aquarium Glaser GmbH, Rodgau, Germany) and kept in tanks with a water temperature of $25\pm1\,^\circ$C and a conductivity of around $270\,\micro\siemens\per\centi\meter$ under a 12\,h:12\,h light-dark cycle. - -Before surgery, the animals were deeply anesthetized via bath application of a solution of MS222 (120\,mg/l, PharmaQ, Fordingbridge, UK) buffered with Sodium Bicarbonate (120\,mg/l). The posterior anterior lateral line nerve (pALLN) was exposed by making a small cut into the skin covering the nerve. The cut was placed dorsal of the operculum just before the nerve descends towards the anterior lateral line ganglion (ALLNG). Those parts of the skin that were to be cut were locally anesthetized by cutaneous application of liquid lidocaine hydrochloride (20\,mg/ml, bela-pharm GmbH). During the surgery, water supply was ensured by a mouthpiece to maintain anesthesia with a solution of MS222 (100\,mg/l) buffered with Sodium Bicarbonate (100\,mg/l). After surgery, fish were immobilized by intramuscular injection of from 25\,$\micro$l to 50\,$\micro$l of tubocurarine (5\,mg/ml dissolved in fish saline; Sigma-Aldrich). -Respiration was then switched to normal tank water and the fish was transferred to the experimental tank. - -\subsection{Electrophysiological recordings} -For the recordings fish were positioned centrally in the experimental tank, with the major parts of their body submerged into the water. Those body parts that were above the water surface were covered with paper tissue to avoid drying of the skin. Local analgesia was refreshed in intervals of two hours by cutaneous application of Lidocaine (2\,\%; bela-pharm, Vechta, Germany) around the surgical wounds. Electrodes (borosilicate; 1.5\,mm outer diameter; GB150F-8P; Science Products, Hofheim, Germany) were pulled to a resistance of 50--100\,\mega\ohm{} (model P-97; Sutter Instrument, Novato, CA) and filled with 1\,M KCl solution. Electrodes were fixed in a microdrive (Luigs-Neumann, Ratingen, Germany) and lowered into the nerve. Recordings of electroreceptor afferents were amplified and lowpass filtered at 10\,kHz (SEC-05, npi-electronics, Tamm, Germany, operated in bridge mode). All signals, neuronal recordings, recorded EOD, and the generated stimulus, were digitized with sampling rates of 20 or 40\,kHz (PCI-6229, National Instruments, Austin, TX). RELACS (\url{https://github.com/relacs/relacs}) running on a Linux computer was used for online spike and EOD detection, stimulus generation, and calibration. Recorded data was then stored on the hard drive for offline analysis. - -\subsection{Identification of P-units and ampullary cells} -Recordings were classified as P-units if baseline action potentials phase locked to the EOD with vectors strengths between 0.7 and 0.95, a baseline firing rate larger than 30\,Hz, a serial correlation of subsequent interspike intervals below zero, a coefficient of variation of baseline interspike intervals below 1.5 und during stimulation below 2. As ampullary cells we classified recordings with vector strengths below 0.15, baseline firing rate above 10\,Hz, baseline CV below 0.18, CV during stimulation below 1.0, and a response modulation during stimulation below 80\,Hz \citep{Grewe2017}. We here selected only those cells of which the neuron's baseline activity as well as the responses to band-limited white noise stimuli were recorded. - -\subsection{Electric field recordings} -For monitoring the EOD without the stimulus, two vertical carbon rods ($11\,\centi\meter$ long, 8\,mm diameter) in a head-tail configuration were placed isopotential to the stimulus. Their signal was differentially amplified with a gain factor between 100 and 500 (depending on the recorded animal) and band-pass filtered (3 to 1500\,Hz pass-band, DPA2-FX; npi electronics, Tamm, Germany). For an estimate of the transdermal potential that drives the electroreceptors, two silver wires spaced by 1\,cm were located next to the left gill of the fish and orthogonal to the fish's longitudinal body axis (amplification 100 to 500 times, band-pass filtered with 3 to 1\,500\,Hz pass-band, DPA2-FX; npi-electronics, Tamm, Germany). This local EOD measurement recorded the combination of the fish's own EOD and the applied stimulus. - - -\subsection{Stimulation}\label{rammethods} -Electric stimuli were attenuated (ATN-01M, npi-electronics, Tamm, Germany), isolated from ground (ISO-02V, npi-electronics, Tamm, Germany) and delivered via two horizontal carbon rods (30 cm length, 8 mm diameter) located $15\,\centi\meter$ laterally to the fish. -The fish were stimulated with band-limited white noise stimuli with a cut-off frequency of 150\,Hz (ampullary), 300 or 400\,Hz (P-units). The stimulus intensity is given as a contrast, i.e. the standard deviation of the noise stimulus relative to the fish's EOD amplitude. The contrast varied between 1 and 20\,\% (median 5\,\%) for P-units and 2.5 and 20\,\% (median 5\,\%) for ampullary cells. Only noise stimuli with a duration of at least 2\,s (maximum of 50\,s, median 10\,s) and enough repetitions to results in at least 100 FFT segments (see below, P-units: 100--1520, median 313, ampullary cells: 105 -- 3648, median 722) were included into the analysis. When ampullary cells were recorded, the white noise was directly applied as the stimulus. To create random amplitude modulations (RAM) for P-unit recordings, the noise stimulus was first multiplied with the EOD of the fish (MXS-01M; npi electronics). - -\subsection{Data analysis} Data analysis was done in Python 3 using the packages matplotlib \citep{Hunter2007}, numpy \citep{Walt2011}, scipy \citep{scipy2020}, pandas \citep{Mckinney2010}, nixio \citep{Stoewer2014}, and thunderlab (\url{https://github.com/bendalab/thunderlab}). The P-unit model parameters and spectral analysis algorithms are available at \url{https://github.com/bendalab/punitmodel/tree/v1}. - -\paragraph{Baseline analysis}\label{baselinemethods} -The baseline firing rate $r$ was calculated as the number of spikes divided by the duration of the baseline recording (median 32\,s). The coefficient of variation (CV) of the interspike intervals (ISI) is their standard deviation relative to their mean: $\rm{CV}_{\rm base} = \sqrt{\langle (ISI- \langle ISI \rangle) ^2 \rangle} / \langle ISI \rangle$. If the baseline was recorded several times in a recording, the measures from the longest recording were taken. - -\paragraph{White noise analysis} \label{response_modulation} -When stimulated with band-limited white noise stimuli, neuronal activity is modulated around the average firing rate that is similar to the baseline firing rate and in that way encodes the time-course of the stimulus. For an estimate of the time-dependent firing rate $r(t)$ we convolved each spike train with normalized Gaussian kernels with standard deviation of 2\,ms, if not mentioned otherwise, and averaged the resulting single-trail firing rates over trials. The response modulation quantifies the variation of $r(t)$ computed as the standard deviation in time $\sigma_{s} = \sqrt{\langle (r(t)-\langle r(t) \rangle_t )^2\rangle_t}$, where $\langle \cdot \rangle_t$ denotes averaging over time. - -\paragraph{Spectral analysis}\label{susceptibility_methods} -To characterize the relation between the spiking response evoked by white-noise stimuli, we estimated the first- and second-order susceptibilities in the frequency domain. For this we converted spike times into binary vectors $x_k$ with $\Delta t = 0.5$\,ms wide bins that are set to 2\,kHz where a spike occurred and zero otherwise. Fast Fourier transforms (FFT) $S(\omega)$ and $X(\omega)$ of the stimulus $s_k$ (also down-sampled to a sampling rate of 2\,kHz) and $x_k$, respectively, were computed numerically according to -\begin{equation} - \label{fourier} - X(\omega) = \sum_{k=0}^{N-1} \, x_k e^{- i \omega k} -\end{equation} -for $N=512$ long segments of $T=N \Delta t = 256$\,ms duration with no overlap, resulting in a spectral resolution of about 4\,Hz. Note, that for a real Fourier integral a factor $\Delta t$ is missing. For simplicity we use angular frequencies $\omega=2\pi f$ instead of frequencies $f$. - -In the experimental data the duration of the noise stimuli varied and they were presented once or repeatedly (frozen noise). For the analysis we discarded the responses within the initial 200\,ms of stimulation in each trial. To make the recordings comparable we always used the first 100 segments from as many trials as needed for the following analysis. - -In the simulations we generated for each trial a new realization of the noise stimulus. We discarded the first 500\,ms of the response and used the following 10 FFT segments for the analysis. This was repeated until $10^6$ or $10^7$ FFT segments were collected. - -The power spectrum of the stimulus $s(t)$ was estimated as -\begin{equation} - \label{powerss} - S_{ss}(\omega) = \frac{\Delta t}{N} \langle S(\omega) S^*(\omega) \rangle -\end{equation} -with $S^*$ being the complex conjugate of $S$ and $\langle \cdot \rangle$ denoting the average over the FFT segments. The factor in front is $\Delta t^2$ from the missing integration factors of the two Fourier transforms, \eqnref{fourier}, divided by $T=N \Delta t$, needed to make this a proper power spectral density. - -The power spectrum of the spike trains $S_{xx}(\omega)$ was estimated accordingly. Likewise, the cross-spectrum $S_{xs}(\omega)$ between stimulus and evoked spike trains was estimated according to -\begin{equation} - \label{crossxs} - S_{xs}(\omega) = \frac{\Delta t}{N} \langle X(\omega) S^*(\omega) \rangle -\end{equation} - -The first-order susceptibility (transfer function) -\begin{equation} - \label{linearencoding_methods} - \chi_{1}(\omega) = \frac{S_{xs}(\omega) }{S_{ss}(\omega) } -\end{equation} -was then computed from $S_{xs}(\omega)$ and $ S_{ss}(\omega)$. We report $\chi_{1}(\omega)$ in Hz/\%, i.e. firing rate per percent stimulus contrast. Multiplying $\chi_{1}(\omega)$ with the contrast of a sinusoidal stimulus in percent results in the amplitude of the evoked firing rate modulation in Hertz. - -The second-order cross-spectrum -\begin{equation} - \label{crossxss} - S_{xss}(\omega_{1},\omega_{2}) = \frac{\Delta t^2}{N} \langle X(\omega_{1} + \omega_{2}) S^*(\omega_{1}) S^*(\omega_{2}) -\end{equation} -quantifies nonlinear interactions that generate responses at the sum and difference (for negative $\omega_1$ or $\omega_2$) evoked by two stimulus frequencies $\omega_1$ and $\omega_2$. -The second-order susceptibility -\begin{equation} - \label{chi2} - \chi_{2}(\omega_{1}, \omega_{2}) = \frac{S_{xss} (\omega_{1},\omega_{2})}{2 S_{ss} (\omega_{1}) S_{ss} (\omega_{2})} -\end{equation} -normalizes the second-order cross-spectrum by the spectral power at the two stimulus frequencies. We report $\chi_{2}(\omega_{1}, \omega_{2})$ in Hz/\%$^2$, i.e. firing rate per percent stimulus contrast squared. -Throughout the manuscript we only show the absolute values of the complex-valued second-order susceptibility matrix and ignore the corresponding phases. - -\paragraph{Susceptibility index}\label{projected_method} -We expected to see a sharp ridge in the second-order susceptibility at $\omega_1 + \omega_2 = r$ \citep{Voronenko2017,Franzen2023}. To characterize this in a single number we computed a susceptibility index. First, we projected the absolute values of the second-order susceptibility matrix onto the diagonal by averaging over anti-diagonal elements. In this projection $D(f)$ we took the position of the maximum -\begin{equation} - \label{dmax} - f_{\rm peak} = {\rm argmax} D(r - 50\,\rm{Hz} \leq f \leq r + 50\,\rm{Hz}) -\end{equation} -within $\pm 50$\,Hz of the neuron's baseline firing rate $r$ as the position of the expected peak. For an estimate of the noise-floor surrounding this peak we averaged over 10\,Hz wide windows 10\,Hz to the left and right of the peak: -\begin{equation} - \label{dref} - D_{\rm ref} = \frac{1}{2}\left(\langle D(f_{\rm peak}-20\,\rm{Hz} \leq f \leq f_{\rm peak}-10\,\rm{Hz})\rangle_f + \langle D(f_{\rm peak}+10\,\rm{Hz} \leq f \leq f_{\rm peak}+20\,\rm{Hz}))\rangle_f \right) -\end{equation} -The size of the peak relative to this reference is then the susceptibility index -\begin{equation} - \label{siindex} - SI = \frac{D(f_{\rm peak})}{D_{\rm ref}} -\end{equation} -Values larger than one indicate a sharp ridge in the susceptibility matrix close to where the stimulus frequencies add up to the baseline firing rate. - - -\begin{figure*}[t] - \includegraphics[width=\columnwidth]{flowchart.pdf} - \caption{\label{flowchart} - Architecture of the P-unit model. Each row illustrates subsequent processing steps for three different stimulation regimes: (i) baseline activity without external stimulus, only the fish's self-generated EOD (the carrier, \eqnref{eq:eod}) is present. (ii) RAM stimulation, \eqnref{eq:ram_equation}. The amplitude of the EOD carrier is modulated with a weak (2\,\% contrast) band-limited white-noise stimulus. (iii) Noise split, \eqnsref{eq:ram_split}--\eqref{eq:Noise_split_intrinsic}, where 90\,\% of the intrinsic noise is replaced by a RAM stimulus, whose amplitude is scaled to maintain the mean firing rate and the CV of the ISIs of the model's baseline activity. As an example, simulations of the model for cell ``2012-07-03-ak'' are shown. \figitem{A} The stimuli are thresholded, \eqnref{eq:threshold2}, by setting all negative values to zero. \figitem{B} Subsequent low-pass filtering, \eqnref{eq:dendrite}, attenuates the carrier and carves out the AM signal. \figitem{C} Intrinsic Gaussian white-noise is added to the signals shown in \panel{B}. Note the reduced internal noise amplitude in the noise split (iii) condition. \figitem{D} Spiking output of the LIF model, \eqnsref{eq:LIF}--\eqref{spikethresh}, in response to the sum of \panel{B} and \panel{C}. \figitem{E} Power spectra of the LIF neuron's spiking activity. Both, baseline activity (\panel[i]{E}) and noise split (\panel[iii]{E}), have the same peaks in the response spectrum at $r$, $f_{EOD} - r$, $f_{EOD}$, and $f_{EOD} + r$. With RAM stimulation (\panel[ii]{E}), the peak at the baseline firing rate, $r$, is washed out.} -\end{figure*} - -\subsection{Leaky integrate-and-fire models for P-units}\label{lifmethods} - -Modified leaky integrate-and-fire (LIF) models were constructed to reproduce the specific firing properties of P-units \citep{Chacron2001, Sinz2020, Barayeu2023}. Its basic components (static nonlinearity, low-pass filtering and spike generation) are equivalent to models of hair cells in auditory systems \citep{Eggermont1983}. The sole driving input into the P-unit model during baseline, i.e. when no external stimulus was given, is the fish's own EOD, modeled as a cosine wave -\begin{equation} - \label{eq:eod} - y(t) = y_{EOD}(t) = \cos(2\pi f_{EOD} t) - \end{equation} -with EOD frequency $f_{EOD}$ and an amplitude of one. - -To mimic the interaction with other fish, the EODs of a second or third fish with EOD frequencies $f_1$ and $f_2$, respectively, were added to the normalized EOD, \eqnref{eq:eod}, of the receiving fish according to their contrasts, $c_1$ and $c_2$ at the position of the receiving fish: -\begin{equation} - \label{eq:modelbeats} - y(t) = \cos(2\pi f_{EOD} t) + c_1 \cos(2\pi f_1 t) + c_2\cos(2\pi f_2 t) -\end{equation} -For two fish, $c_2 = 0$. - -Random amplitude modulations (RAMs) were simulated by first generating the AM as a band-limited white noise stimulus $s(t)$. For this, random real and imaginary numbers were drawn from Gaussian distributions for each frequency component in the range from 0 to 300\,Hz in the Fourier domain \citep{Billah1990,Skorjanc2023}. By means of the inverse Fourier transform, the time course of the RAM stimulus, $s(t)$, was generated. The input to the model was then -\begin{equation} - \label{eq:ram_equation} - y(t) = (1+ s(t)) \cos(2\pi f_{EOD} t) -\end{equation} -The contrast $c$ of the RAM is the standard deviation of the RAM relative to the amplitude of the receiving fish. - -First, the input $y(t)$ is thresholded by setting negative values to zero: -\begin{equation} - \label{eq:threshold2} - \lfloor y(t) \rfloor_0 = \left\{ \begin{array}{rcl} y(t) & ; & y(t) \ge 0 \\ 0 & ; & y(t) < 0 \end{array} \right. -\end{equation} -(\subfigrefb{flowchart}{A}). This thresholds models the transfer function of the synapses between the primary receptor cells and the afferent. Together with a low-pass filter -\begin{equation} - \label{eq:dendrite} - \tau_{d} \frac{d V_{d}}{d t} = -V_{d}+ \lfloor y(t) \rfloor_{0} -\end{equation} -the threshold operation is required for extracting the amplitude modulation from the input \citep{Barayeu2024}. The low-pass filter models passive signal conduction in the afferent's dendrite (\subfigrefb{flowchart}{B}) and $\tau_{d}$ is the membrane time constant of the dendrite. Dendritic low-pass filtering was also necessary to reproduce the loose coupling of P-unit spikes to the EOD while maintaining high sensitivity at small amplitude modulations. - -The dendritic voltage $V_d(t)$ is then fed into a stochastic leaky integrate-and-fire (LIF) model with adaptation, -\begin{equation} - \label{eq:LIF} - \tau_{m} \frac{d V_{m}}{d t} = - V_{m} + \mu + \beta V_{d} - A + \sqrt{2D}\;\xi(t) -\end{equation} -where $\tau_{m}$ is the membrane time-constant, $\mu$ is a fixed bias current, $\beta$ is a scaling factor for $V_{d}$, $A$ is an inhibiting adaptation current, and $\sqrt{2D}\;\xi(t)$ is a Gaussian white noise with strength $D$. Note, that all state variables, membrane voltages $V_d$ and $V_m$ as well as the adaptation current $A$, are dimensionless. - -The adaptation current $A$ follows -\begin{equation} - \label{eq:adaptation} - \tau_{A} \frac{d A}{d t} = - A -\end{equation} -with adaptation time constant $\tau_A$ \citep{Benda2003,Benda2005,Benda2010}. - -Whenever the membrane voltage $V_m(t)$ crosses the spiking threshold $\theta=1$, a spike was generated, $V_{m}(t)$ was reset to $0$, the adaptation current was incremented by $\Delta A/\tau_A$, and integration of $V_m(t)$ was paused for the duration of a refractory period $t_{ref}$ (\subfigrefb{flowchart}{D}): -\begin{equation} - \label{spikethresh} - V_m(t) \ge \theta \; : \left\{ \begin{array}{rcl} V_m & \mapsto & 0 \\ A & \mapsto & A + \Delta A/\tau_A \end{array} \right. -\end{equation} - -The P-unit models were integrated by the Euler forward method with a time-step of $\Delta t = 0.05$\,ms. For each trial of a simulation, $V_{m}$ was drawn from a uniform distribution between 0 and 1 and the initial value of $A$ was jittered by adding a random number drawn from a normal distribution with standard deviation of 2\,\% of its initial value. Then the first 500\,ms of any simulation were discarded to remove remaining transients. - -The eight free parameters of the P-unit model $\beta$, $\tau_m$, $\mu$, $D$, $\tau_A$, $\Delta_A$, $\tau_d$, and $t_{ref}$, were fitted to both the baseline activity (baseline firing rate, CV of ISIs, serial correlation of ISIs at lag one, and vector strength of spike coupling to EOD) and the responses to step increases and decreases in EOD amplitude (onset and steady-state responses, effective adaptation time constant, \citealp{Benda2005}) of recorded P-units. Model parameters of all 39 cells are summarized in file \texttt{models.csv} of our \texttt{punitmodel} repository at \url{https://github.com/bendalab/punitmodel/tree/v1}. - - -\subsection{Noise split} -\label{intrinsicsplit_methods} -Based on the Furutsu-Novikov theorem \citep{Furutsu1963,Novikov1965,Lindner2022,Egerland2020}, we split the total noise, $\sqrt{2D}\;\xi(t)$, of a LIF model, \eqnref{eq:LIF}, into two parts. The first part is the intrinsic noise term, $\sqrt{2D \, \alpha_{\rm noise}}\;\xi(t)$, whose strength is reduced by a factor $\alpha_{\rm noise}=0.1$ (\subfigrefb{flowchart}\,\panel[iii]{C}). The second part replaces the now missing intrinsic noise by a driving input signal $s_{\xi}(t)$, a RAM stimulus with frequencies up to 300\,Hz (\subfigrefb{flowchart}\,\panel[iii]{A}). The LIF model with noise split then reads -\begin{eqnarray} - \label{eq:ram_split} - y(t) & = & (1+ s_\xi(t)) \cos(2\pi f_{EOD} t) \\ - \label{eq:Noise_split_intrinsic_dendrite} - \tau_{d} \frac{d V_{d}}{d t} & = & -V_{d}+ \lfloor y(t) \rfloor_{0} \\ - \label{eq:Noise_split_intrinsic} - \tau_{m} \frac{d V_{m}}{d t} & = & - V_{m} + \mu + \beta V_{d} - A + \sqrt{2D \, \alpha_{\rm noise}}\;\xi(t) -\end{eqnarray} -Both, the reduced intrinsic noise and the RAM stimulus, need to replace the original intrinsic noise. Because the RAM stimulus is band-limited and undergoes some transformations before it is added to the reduced intrinsic noise, it is not \textit{a priori} clear, what the amplitude of the RAM should be. We bisected the amplitude of $s_\xi(t)$, until the CV of the resulting interspike intervals matched the one of the original model's baseline activity. The second-order cross-spectra, \eqnref{crossxss}, were computed between the RAM stimulus $s_{\xi}(t)$ and the spike train $x(t)$ it evoked. In this way, the effective signal-to-noise ratio can be increased while maintaining the total noise in the system. - - -\section{Acknowledgements} -Supported by SPP 2205 ``Evolutionary optimisation of neuronal processing'' by the DFG, project number 430157666. -We thank Tim Hladnik, Henriette Walz, Franziska Kuempfbeck, Fabian Sinz, Laura Seidler, Eva Vennemann, and Ibrahim Tunc for the data they recorded over the years in our lab - - %\bibliographystyle{apalike}%alpha}%}%alpha}%apalike} \bibliography{journalsabbrv,references} % \bibliographystyle{apalike} %or any other style you like