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

finished model in methods

Browse Source
This commit is contained in:
Jan Benda 2025-07-01 09:30:22 +02:00
parent 1e1f66b1ea
commit 28a6ef2bd8
1 changed files with 27 additions and 27 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

54
nonlinearbaseline.tex
View File

@ -728,7 +728,12 @@ When stimulated with band-limited white noise stimuli, neuronal activity is modu
% The average firing rate during stimulation, $r_s = \langle r(t) \rangle_t$, is given by the temporal average $\langle \cdot \rangle_t$ over the duration of the stimulus of the trial-averaged firing rate. To quantify how strongly a neuron was driven by the stimulus, we computed the response modulation as the standard deviation $\sigma_{s} = \sqrt{\langle (r(t)-r_s)^2\rangle_t}$ of the trial-averaged firing rate.
\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 according to $X(\omega) = \sum_{k=0}^{n} \, x_k e^{- i \omega k}$ for $n=512$ long segments of 256\,ms duration with no overlap, resulting in a spectral resolution of about 4\,Hz. For simplicity we use angular frequencies $\omega=2\pi f$ instead of frequencies $f$.
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} \, 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.
@ -739,7 +744,9 @@ The power spectrum of the stimulus $s(t)$ was estimated as
\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 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
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
@ -793,12 +800,12 @@ Values larger than one indicate a sharp ridge in the susceptibility matrix close
\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}. 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
Modified leaky integrate-and-fire (LIF) models were constructed to reproduce the specific firing properties of P-units \citep{Chacron2001, Sinz2020, Barayeu2023}. 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}
\carrierinput = y_{EOD}(t) = \cos(2\pi f_{EOD} t)
\end{equation}
with EOD frequency $f_{EOD}$ and an amplitude normalized to one.
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}
@ -829,9 +836,9 @@ the threshold operation is required for extracting the amplitude modulation from
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 + \alpha V_{d} - A + \sqrt{2D}\;\xi(t)
\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, $\alpha$ is a scaling factor for $V_{d}$, $A$ is an inhibiting adaptation current, and $\sqrt{2D}\;\xi(t)$ is a white noise with strength $D$. Note, that all state variables, membrane voltages $V_d$ and $V_m$ as well as the adaptation current, are dimensionless.
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}
@ -840,35 +847,28 @@ The adaptation current $A$ follows
\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$, and integration of $V_m(t)$ was paused for the duration of a refractory period $t_{ref}$ (\subfigrefb{flowchart}{D}):
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. The intrinsic noise $\xi(t)$ (\eqnref{eq:LIF}, \subfigrefb{flowchart}{C}) was added by drawing a random number from a normal distribution $\mathcal{N}(0,\,1)$ with zero mean and standard deviation of one in each time step $i$. This number was multiplied with $\sqrt{2D}$ and divided by $\sqrt{\Delta t}$. For each trial of a simulation, the variables $A$, $V_{d}$ and $V_{m}$ were drawn from a distribution of initial values, estimated from a 100\,s long simulation of baseline activity after an initial 100\,s long integration that was discarded as a transient.
%\begin{equation}
% \label{eq:LIFintegration}
% V_{m_{i+1}} = V_{m_i} + \left(-V_{m_i} + \mu + \alpha V_{d_i} - A_i + \sqrt{\frac{2D}{\Delta t}}\mathcal{N}(0,\,1)_i\right) \frac{\Delta t}{\tau_m}
%\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 remainig transients.
%\paragraph{Fitting the model to recorded P-units}
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 (table~\ref{modelparams}).
\notejb{add table with all 39 cells}
\begin{table*}[tp]
\caption{\label{modelparams} Model parameters of LIF models, fitted to 3 electrophysiologically recorded P-units \citep{Ott2020}.}
\begin{tabular}{lrrrrrrrr}
\hline
\bfseries cell & \bfseries $\mu$ & \bfseries $\beta$ & \bfseries $\tau_{m}$/ms & \bfseries $D$/$\mathbf{ms}$ & \bfseries $\tau_{A}$/ms & \bfseries $\Delta_A$ & \bfseries $\tau_{d}$/ms & \bfseries $t_{ref}$/ms \\\hline
2012-07-03-ak& $-1.32$& $10.6$& $1.38$& $0.001$& $96.05$& $0.01$& $1.18$& $0.12$ \\
2013-01-08-aa& $0.59$& $4.5$& $1.20$& $0.001$& $37.52$& $0.01$& $1.18$& $0.38$ \\
2018-05-08-ae& $-21.09$& $139.6$& $1.49$& $0.214$& $123.69$& $0.16$& $3.93$& $1.31$ \\
\hline
\end{tabular}
\end{table*}
In \figrefb{fig:noisesplit} and \figrefb{fig:model_full}, for each of the stimulus and response segments needed for the spectral analysis, \eqnsref{powerss}--\eqref{chi2}, a simulation was run. The first second was discarded and the analysis was based on the last second of the data. The resulting spectra thus have a spectral resolution of 1\,Hz. To speed up model simulations for \figrefb{fig:modelsusceptcontrasts} and \figrefb{fig:trialnr} we simulated up to $10^6$ trials of 3\,s duration. In each trial the first 0.5\,s were skipped, yielding 10 FFT segments of $T=0.25$\,s (4\,Hz resolution).
\notejb{tag git repo and insert reference to it}
% \begin{table*}[tp]
% \caption{\label{modelparams} Model parameters of LIF models, fitted to 3 electrophysiologically recorded P-units \citep{Ott2020}.}
% \begin{tabular}{lrrrrrrrr}
% \hline
% \bfseries cell & \bfseries $\mu$ & \bfseries $\beta$ & \bfseries $\tau_{m}$/ms & \bfseries $D$/$\mathbf{ms}$ & \bfseries $\tau_{A}$/ms & \bfseries $\Delta_A$ & \bfseries $\tau_{d}$/ms & \bfseries $t_{ref}$/ms \\\hline
% 2012-07-03-ak& $-1.32$& $10.6$& $1.38$& $0.001$& $96.05$& $0.01$& $1.18$& $0.12$ \\
% 2013-01-08-aa& $0.59$& $4.5$& $1.20$& $0.001$& $37.52$& $0.01$& $1.18$& $0.38$ \\
% 2018-05-08-ae& $-21.09$& $139.6$& $1.49$& $0.214$& $123.69$& $0.16$& $3.93$& $1.31$ \\
% \hline
% \end{tabular}
% \end{table*}
\subsection{Noise split}