Numerical uncertainty can critically affect simulations of mechanistic models in neuroscience

Understanding neural computation on the mechanistic level requires biophysically realistic neuron models. To analyze such models one typically has to solve systems of coupled ordinary differential equations (ODEs), which describe the dynamics of the underlying neural system. These ODEs are solved numerically with deterministic ODE solvers that yield single solutions with either no or only a global scalar bound on precision. To overcome this problem, we propose to use recently developed probabilistic solvers instead, which are able to reveal and quantify numerical uncertainties, for example as posterior sample paths. Importantly, these solvers neither require detailed insights into the kinetics of the models nor are they difficult to implement. Using these probabilistic solvers, we show that numerical uncertainty strongly affects the outcome of typical neuroscience simulations, in particular due to the non-linearity associated with the generation of action potentials. We quantify this uncertainty in individual single Izhikevich neurons with different dynamics, a large population of coupled Izhikevich neurons, single Hodgkin-Huxley neuron and a small network of Hodgkin-Huxley-like neurons. For commonly used ODE solvers, we find that the numerical uncertainty in these models can be substantial, possibly jittering spikes by milliseconds or even adding or removing individual spikes from the simulation altogether.Author summaryComputational neuroscience is built around computational models of neurons that allow the simulation and analysis of signal processing in the central nervous system. These models come typically in the form of ordinary differential equations (ODEs). The solution of these ODEs is computed using solvers with finite accuracy and, therefore, the computed solutions deviate from the true solution. If this deviation is too large but goes unnoticed, this can potentially lead to wrong scientific conclusions.A field in machine learning called probabilistic numerics has recently developed a set of probabilistic solvers for ODEs, which not only produce a single solution of unknown accuracy, but instead yield a distribution over simulations. Therefore, these tools allow one to address the problem state above and quantitatively analyze the numerical uncertainty inherent in the simulation process.In this study, we demonstrate how such solvers can be used to quantify numerical uncertainty in common neuroscience models. We study both Hodgkin-Huxley and Izhikevich neuron models and show that the numerical uncertainty in these models can be substantial, possibly jittering spikes by milliseconds or even adding or removing individual spikes from the simulation altogether. We discuss the implications of this finding and discuss how our methods can be used to select simulation parameters to trade off accuracy and speed.