role_of_noise/sigma.tex

181 lines
20 KiB
TeX

\subsection*{Determining noise in real world}
While in simulations we can control the noise parameter directly, we cannot do so in electrophysiological experiments.
Therefore, we need a way to quantify "noisiness".
One such way is by using the activation curve of the neuron, fitting a function and extracting the parameters from this function.
Stocks (2000) uses one such function to simulate groups of noisy spiking neurons:
\begin{equation}
\label{errorfct}\frac{1}{2}\erfc\left(\frac{\theta-x}{\sqrt{2\sigma^2}}\right)
\end{equation}
where $\sigma$ is the parameter quantifying the noise (figure \ref{idealizedactivation}). A neuron with a $\sigma$ of 0 would be a perfect thresholding mechanism. Firing probability for all inputs below the threshold is 0, and firing probability for all inputs above is 1. Large values mean a flat activation curve. Neurons with such an activation curve will likely fire even for some signals below the firing threshold, while it will sometimes not fire for inputs above the firing threshold. Its firing behaviour is influenced less by the signal, which indicates noisiness.
We also tried different other methods of quantifying noise commonly used (citations), but none of them worked as well as the errorfunction fit (fig. \ref{sigmasimulation} b)-d)).
\subsection*{Methodology}
The signal values were binned in 50 bins. The result is a discrete Gaussian distribution around 0mV, the mean of the signal, as is expected from the way the signal was created.
We have to account for the delay between the moment we paly the signal and when it gets processed in the cell, which can for example depend on the position of the cell on the skin. We calculate the cross correlation between the signal and the discrete output spikes. The position of the peak of the correlation is the time shift for which the signal influences the result of the output the most.
Then, for every spike, we note the value of the signal at the time of the spike minus the time shift. The result is a histogram, where each signal value bin has an associated number of spikes. This histogram is then normalized by the distribution of the signal. The result is another histogram, whose values are firing frequencies for each signal value. Because those frequencies are just firing probabilities multiplied with equal time steps, we can fit a Gaussian error function to those probabilities.
\subsection*{Simulation}
To confirm that the $\sigma$ parameter estimated from the fit is indeed a good measure for the noisiness, we validated it against D, the noise parameter from the simulations. We find that there is a strictly monotonous relationship between the two for different sets of simulation parameters. Other parameters often used to determine noisiness (citations) such as the variance of the spike PSTH, the coefficient of variation (CV) of the interspike interval are not as useful. In figure \ref{noiseparameters} we see why. The variance of the psth is not always monotonous in D and is very flat for low values of D.
%describe what happens to the others
%check Fano-factor maybe?
\begin{figure}
\includegraphics[width=0.45\linewidth]{img/ordnung/base_D_sigma}
\includegraphics[width=0.45\linewidth]{img/dataframe_scatter_D_normalized_psth_1ms_test_tau}
\includegraphics[width=0.45\linewidth]{img/dataframe_scatter_D_psth_5ms_test}
% \includegraphics[width=0.45\linewidth]{img/dataframe_scatter_D_cv_test}
\caption{a)The parameter \(\sigma\) as a function of the noise parameter D in LIF-simulations. There is a strictly monotonous relationship between the two, which allows us to use \(\sigma\) as a susbtitute for D in the analysis of electrophysiological experiments. b-d) different other parameters commonly used to quantify noise. None of these functions is stricly monotonous and therefore none is useful as a substitute for D. b) Peri-stimulus time histogram (PSTH) of the spikes with a bin width of 1ms, normalized by c) PSTH of the spikes with a bin width of 5ms. d) coefficient of variation (cv) of the interspike-intervals.}
\label{noiseparameters}
\end{figure}
We tried several different bin sizes (30 to 300 bins) and spike widths. There was little difference between the different parameters (see appendix).
\subsection*{Electrophysiology}
We can see from figure \ref{sigmafits} that the fits look very close to the data. For very weak and very strong inputs, the firing rates themselves become noisy, because the signal only assumes those values rarely. This is especially noticeable for strong inputs, as there are more spikes there, and therefore large fluctuations, while there is very little spiking anyway for weak inputs.
\begin{figure}
\includegraphics[width=0.4\linewidth]{img/ordnung/cropped_fitcurve_0_2010-08-31-ad-invivo-1_0.pdf}
\includegraphics[width=0.4\linewidth]{img/ordnung/cropped_fitcurve_0_2010-08-11-aq-invivo-1_0.pdf}
% cropped_fitcurve_0_2010-08-31-ad-invivo-1_0.pdf: 0x0 px, 300dpi, 0.00x0.00 cm, bb=
\caption{Histogram of spike count distribution (firing rate) and errorfunction fits. 50 bins represent different values of the Gaussian distributed input signal [maybe histogram in background again]. The value of each of those bins is the number of spikes during the times the signal was in that bin. Each of the values was normalized by the signal distribution. To account for delay, we first calculated the cross-correlation of signal and spike train and took its peak as the delay. The lines show fits according to equation \eqref{errorfct}. Left and right plots show two different cells, one with a relatively narrow distribution and one with a distribution that is more broad, as indicated by the parameter \(\sigma\). Different amounts of bins (30 and 100) showed no difference in resulting parameters.}
\label{sigmafits}
\end{figure}
%TODO insert plot with sigma x-axis and delta_cf on y-axis here; also, plot with sigma as function of firing rate, also absoulte cf for different population size as function of sigma.
When we group neurons by their noise and plot coding fraction as a function of population size for averages of the groups, we see results similar to what we see for simulations.
Noisier cells have a lower coding fraction for small populations. For increasing population size, coding fraction increases for all groups, but the increase is much larger for noisy cells. For large population sizes the noisy cells show a better linear encoding of the signal than the more regular cells.
\begin{figure}
\includegraphics[width=0.45\linewidth]{img/ordnung/sigma_popsize_curves_0to300}
\caption{Left: Coding fraction as a function of population size for all recorded neurons. Color are by \(\sigma\) from the fit of the function in equation \ref{errorfct}, so that there are roughly an equal number of neurons in each category. Red: \(\sigma = \) 0 to 0.5, pink: 0.5 to 1.0, purple: 1.0 to 1.5, blue: 1.5 and above. Thick colored lines are average of the neurons in each group. For a population size of 1, coding fraction descreases on average with increasing \(\sigma\). As population sizes increase, coding fraction for weak noise neurons quickly stops increasing. Strong noise neurons show better coding performance for larger popuation sizes (about 8 to 32 neurons). Right [missing]: Increase in coding as a function of sigma. y-axis shows the difference in coding fraction between N=1 and N=32,}
\label{ephys_sigma}
\end{figure}
\subsection*{Determining the strength of noise in a real world example}
While in simulations we can control the noise parameter directly, we cannot do so in electrophysiological experiments.
Therefore, we need a way to quantify the intrinsic noise of the cell from the output of the measured cells. Common measures to quantify noisiness of neuronal spike trains are not directly correlated with intrinsic noise strength (figure \ref{noiseparameters}). An example for such a measure is the coefficient of variation (cv) of the interspike interval (ISI)\citep{white2000channel, goldberg1984relation,nowak1997influence}. The ISI is the time between each consecutive pair of spikes. The coefficient of variation is then defined as the standard deviation of the ISI divided by the mean ISI. Even though it is frequently used, we find that the cv as a function the intrinsic noise in our LIF-simulations is not necessarily monotonously related. In addition, for different membrane constants, which determine how quickly a neuron reacts to inputs, the same intrinsic noise results in widely different cv-values. Refractory periods also have an influence on the cv.
\notedh{their/holt1996 cv2 looks interesting.}
Another measure which has been used before is the standard variation of the peri-stimulus spike histogram
\citep{mainen1995reliability} \notedh{can't find any paper which did something like we did here, even though Schreiber et al. 2003 A new correlation-based measure of spike timing reliability - claim it's frequently done with psth}. This approach also does not work well, as it also depends on the membrane constant and to a lesser extend the refractory period.
The approach used here uses the activation curve of the neuron, fitting a function to it and extracting the parameters from the fitted function. It is assumed that the neurons show Gaussian noise. The mean of the distribution is the activation threshold and the width of the Gaussian is a measure for noise.
The probability of spiking as a function of the input is then the integral over the Gaussian, i.e. an error function.
Stocks (2000) uses one such function to simulate groups of noisy spiking neurons:
\begin{equation}
\label{errorfct}\frac{1}{2}\erfc\left(\frac{\theta-x}{\sqrt{2\sigma^2}}\right)
\end{equation}
where $\sigma$ is the parameter quantifying the noise (figure ?) %\ref{idealizedactivation}). \notejb{$\sigma$ quantifies the noise in units of the stimulus!!! THis is why this approach might work!}
A neuron with a $\sigma$ of 0 would be a perfect thresholding mechanism. Firing probability for all inputs below the threshold is 0, and firing probability for all inputs above is 1. If $\sigma$ is greater than 0, a neuron with such an activation curve will fire even for some signals below the firing threshold, while it will sometimes not fire for inputs above the firing threshold. For large values of $\sigma$ the activation curve becomes flatter, meaning the probability for inputs below the theshold eliciting a spike is large and the probability that an input above the threshold does not lead to firing is also large. The firing behaviour of such a cell is influenced less by the signal, which indicates noisiness.
However, for strong noise $(>10^{-2} \frac{mV^2}{Hz})$, results are not monotonous anymore. This happens at a point where $\sigma$ becomes large. Therefore, we excluded all values of the unit-less \(\sigma\) larger than two from the following analyses.
\subsection*{Methodology}
The signal was binned according to its amplitude. The result is a discrete Gaussian distribution around 0mV, the mean of the signal, as is expected from the way the signal was created.
After accounting for time delays in signal processing, we make a histogram which contains the distribution of spikes according to signal amplitude.
This histogram is then normalized by the distribution of the signal.
The result is another histogram, where values are firing frequencies for each signal value. Because those frequencies are just firing probabilities multiplied with equal time steps, we can fit a Gaussian error function to those probabilities.
\subsection*{Simulation}
To confirm that the $\sigma$ parameter estimated from the fit is indeed a good measure for the noisiness, we validated it against D, the noise parameter from the simulations. We find that there is a strictly monotonous relationship between the two for different sets of simulation parameters.
%Other parameters often used to determine noisiness (citations) such as the variance of the spike PSTH, the coefficient of variation (CV) of the interspike interval are not as useful. In figure \ref{noiseparameters} we see why. The variance of the psth is not always monotonous in D and is very flat for low values of D.
%describe what happens to the others
%check Fano-factor maybe?
\begin{figure}
\centering
%\includegraphics[width=0.45\linewidth]{img/ordnung/base_D_sigma}\\
\includegraphics[width=0.23\linewidth]{img/simulation_sigma_examples/fitcurve_50hz_0.0002noi500s_0_capped.pdf}
\includegraphics[width=0.23\linewidth]{img/simulation_sigma_examples/fitcurve_50hz_0.001noi500s_0_capped.pdf}
\includegraphics[width=0.23\linewidth]{img/simulation_sigma_examples/fitcurve_50hz_0.1noi500s_0_capped.pdf}
\includegraphics[width=0.23\linewidth]{img/ISI_explanation.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_sigma_membrane_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_cv_membrane_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_psth_1ms_membrane_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_psth_5ms_membrane_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_sigma_refractory_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_cv_refractory_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_psth_1ms_refractory_50.pdf}
\includegraphics[width=0.23\linewidth]{img/cv_psth_sigma_compare/dataframe_scatter_labels_D_psth_5ms_refractory_50.pdf}
\caption{a)The parameter \(\sigma\) as a function of the noise parameter D in LIF-simulations. There is a strictly monotonous relationship between the two, which allows us to use \(\sigma\) as a susbtitute for D in the analysis of electrophysiological experiments.
b-e) Left to right: $\sigma$, CV and standard deviation of the psth with two diffrent kernel widths as a function of D for different membrane constants (4ms, 10ms and 16ms). The membrane constant $\tau$ determines how quickly the voltage of a LIF-neuron changes, with lower constants meaning faster changes. Only $\sigma$ does not change its values with different $\tau$. The CV (c)) is not even monotonous in the case of a timeconstant of 4ms, ruling out any potential usefulness.
f-i) Left to right: $\sigma$, CV and standard deviation of the psth with two diffrent kernel widths as a function of D for different refractory periods (0ms, 1ms and 5ms). Only $\sigma$ does not change with different refractory periods.
}
\label{noiseparameters}
\end{figure}
We tried several different bin sizes (30 to 300 bins) and spike widths. There was little difference between the different parameters (see appendix).
\section*{-----------------------}
%We can use $\sigma$ instead of D*firing_rate: $\sigma$ makes it ind. of fr!!
\subsection*{Electrophysiology}
We find that the fits match the experimental data very well (figure \ref{sigmafits}). For very weak and very strong inputs, the firing rates themselves become noisy, because the signal only assumes those values rarely. This is especially noticeable for strong inputs, as there are more spikes there, and therefore large fluctuations, while there is very little spiking anyway for weak inputs.
% fish_raster.py on oilbird for the eventplot
% instructions.txt enth\"alt python-Befehle um Verteilungen und scatter zu rekonstruieren
\begin{figure}
\centering
\includegraphics[width=0.4\linewidth]{img/sigma/example_spikes_sigma.pdf}
\includegraphics[width=0.28\linewidth]{img/ordnung/cropped_fitcurve_0_2010-08-11-aq-invivo-1_0.pdf}
\includegraphics[width=0.28\linewidth]{img/sigma/cropped_fitcurve_0_2010-08-31-aj-invivo-1_0.pdf}
\notedh{Daraus ergibt sich nicht direkt eine Intuition, wieso das noisy ist. H\"angt einfach sehr am Eingangssignal; wenn es im eventplot (un-)regelm\"a\ss{}ig w\"are, k\"onnten wir auch einfach cv nehmen...}\\
% \includegraphics[width=0.28\linewidth]{img/ordnung/cropped_fitcurve_0_2010-08-31-ad-invivo-1_0.pdf}
\includegraphics[width=0.4\linewidth]{img/fish/dataframe_scatter_sigma_cv.pdf}
\includegraphics[width=0.4\linewidth]{img/fish/dataframe_scatter_sigma_firing_rate.pdf}
\includegraphics[width=0.32\linewidth]{img/fish/sigma_distribution.pdf}
\includegraphics[width=0.32\linewidth]{img/fish/cv_distribution.pdf}
\includegraphics[width=0.32\linewidth]{img/fish/fr_distribution.pdf}
% cropped_fitcurve_0_2010-08-31-ad-invivo-1_0.pdf: 0x0 px, 300dpi, 0.00x0.00 cm, bb=
\caption{Histogram of spike count distribution (firing rate) and errorfunction fits. 50 bins represent different values of the Gaussian distributed input signal [maybe histogram in background again]. The value of each of those bins is the number of spikes during the times the signal was in that bin. Each of the values was normalized by the signal distribution. For very weak and very strong inputs, the firing rates themselves become noisy, because the signal only assumes those values rarely. To account for delay, we first calculated the cross-correlation of signal and spike train and took its peak as the delay. The lines show fits according to equation \eqref{errorfct}. Left and right plots show two different cells, one with a relatively narrow distribution (left) and one with a distribution that is more broad (right), as indicated by the parameter \(\sigma\). An increase of $\sigma$ is equivalent to an broader distribution. Cells with broader distributions are assumed to be noisier, as their thresholding is less sharp than those with narrow distributions. Different amounts of bins (30 and 100) showed no difference in resulting parameters.}
\label{sigmafits}
\end{figure}
% TODO insert plot with sigma x-axis and delta_cf on y-axis here; also, plot with sigma as function of firing rate, also absoulte cf for different population size as function of sigma.
When we group neurons by their noise and plot coding fraction as a function of population size for averages of the groups, we see results similar to what we see for simulations (figure \ref{ephys_sigma} a)):
Noisier cells (larger $\sigma$, purple) have a lower coding fraction for small populations. However, coding fraction mostly stops increasing with population sizes once a population size of about 16 is reached. The increase is much larger for noisy cells (orange). The averages of the coding fraction for the noisy cells does not increase above the coding fraction of the less noisy cells for the population sizes investigated here (N=128). In contrast to the more regular cells, coding fraction is still improving for the noisy cells, so it is plausible that at a certain population size the noisy cells can outperform the less noisy cells.
Indeed, if results are not averaged and single cells are considered, we find that for large population sizes the noisy cells show a better linear encoding of the signal than the more regular cells (figure \ref{ephys_sigma} b), red).
%figures created with box_script.py
\begin{figure}
%\includegraphics[width=0.45\linewidth]{img/ordnung/sigma_popsize_curves_0to300}
\centering
\includegraphics[width=0.65\linewidth]{img/sigma/cf_N_sigma.pdf}
%\includegraphics[width=0.45\linewidth]{img/sigma/cf_N_ex_lines}
\includegraphics[width=0.45\linewidth]{img/sigma/sigma_cf_quot.pdf}%
\includegraphics[width=0.45\linewidth]{img/sigma/check_fr_quot.pdf}%
\caption{Left: Coding fraction as a function of population size for all recorded neurons. Cells are grouped by \(\sigma\) from the fit of the function in equation \ref{errorfct}. Lines are averages over three cells each, with the shading showing the standard deviation. For stronger noise, coding fraction is far smaller for a single neuron. With increasing population size, coding fraction increases much faster for the noisy cells than for the less noisy cells.
Right: Examples for the two cells with lowest, intermediate and highest $\sigma$. For a population size of N=1, the cell with the largest $\sigma$ (brown) has the lowest coding fraction out of all the cells here. The coding fraction of that cell increases hugely with population size. At a population of N=128, coding fraction is second highest among the pictured cells.}
\label{ephys_sigma}
\end{figure}
%The value of $\sigma$ is not signal independent. The same cell can have different values for $\sigma$ for different input signals.