27th Annual Computational Neuroscience Meeting (CNS*2018): Part Two

Reliable function of neuronal networks, especially of those that control vital motor functions, is crucial for the survival and health of organisms. Here

Reliable function of neuronal networks, especially of those that control vital motor functions, is crucial for the survival and health of organisms. Here, we investigate the role of Na+/K+ pump, which is the basic cellular engine maintaining the physiological gradients of Na+ and K+ ions across the membrane. Growing evidence indicates that the Na+/K+ pump current plays important roles in the electrical activity of neurons, contributing to functional and dysfunctional dynamics [1][2][3][4][5][6]. To understand how the pump dynamics contribute to neuronal activity, the kinetics of Na+/K+ pump function has to be quantitatively evaluated versus the dynamics of neurons and the stability of their activity regimes. The premise of this study is that the Na+/K+ pump contributes to the dynamics of neurons on the time scale of the period of their rhythmic bursting activity (6-10 s). In the leech heartbeat CPG, the basic building blocks are half center oscillators (HCO), which are pairs of mutually inhibitory HN interneurons producing alternating bursting activity. The neuropeptide myomodulin speeds up the period of the bursting pattern by increasing h-current and decreasing the Na+/K+ pump current in the heart interneurons (HNs), thus implicating the pump in the bursting dynamics of this CPG [6]. A definitive role of the Na+/K+ pump current in the bursting dynamics of HN neurons was revealed by using monensin, which stimulates the pump by diffusively increasing the intracellular Na+ concentration.
In the presence of h-current application of monensin decreases the period of a leech heartbeat HCO [4]. To make quantitative analysis possible, we developed a dynamic clamp implementation of a hybrid system with living neuron and mathematical model in real time in each of which we can manipulate pump parameters on the fly (Fig. 1). The implemented model had been described in the appendix of [4]. We explored parameter space of the full model and fine tune the real time model to support functional-like dynamics of HCO under variation of the key biophysical parameters. By focusing our analysis on a common network motif of mutually inhibitory units, we generated insights that generalize from our experimentally and computationally more tractable invertebrate system to less accessible vertebrate nervous systems. Synchronization in neurological systems is critical for the survival of many species. Vital functions such as locomotion and mastication, for example, depend upon mechanisms yielding robust and stable neuronal synchronization. Moreover, pathologies including Parkinson's disease and sleep disorders, for instance, are associated with neuronal synchrony deficiencies rendering patients incapable of leading a normal life. In this study we describe a transition recently observed in computer simulations of gap junction coupled neurons. The transition is mediated by a period-doubling cascade followed by chaos in synchronous neurons, with their firing regimes evolving from tonic (fast repetitive spiking) to bursting (periods of repetitive fast spiking followed by periods of quiescence) as a coupling parameter is increased. While tonic-to-bursting transitions play important roles, for instance, in thalamocortical neurons at sleeping transition states (Sherman, Trends Neurosci, 2001), and in sensory-motor nuclei that generate the typical tremors in Parkinson's disease (Llinas and Steriade, J Neurophysiol, 2006), little is known about the mechanisms regulating and controlling theses transitions at the level of the dynamic of the individual networked neurons. We use a Hodgkin-Huxley type model neuron (Rosa et al., Biosystems, 2015) to investigate the transition between tonic and bursting neuronal behaviors in small networks of electrically coupled neurons. Numerical simulations show that distinct neurons, one tonic and the other bursting, reciprocally coupled via gap-junctions, depending upon the individual characteristics of the two neurons, may synchronize either in the tonic or in the bursting regime, remaining in the state in which they first synchonized for extended increments in the strength of their coupling. However, we also found that in some cases, the two neurons synchronize initially in the tonic regime and with increased coupling strength, undergo a period doubling bifurcation cascade in route to chaos, go through chaos and then, still in synchrony, go into the bursting regime. Intriguing, we noticed that some peculiar common features of the independent single neurons are preserved when they are coupled and synchronized. For example, the characteristic firing rate at the border between tonic and bursting regimes for the individual neuron is passed on to the col- Cortical resting-state dynamics is organized on multiple spatiotemporal scales and involves cell-type-specific spike rates, slow and fast fluctuations, clustered inter-area correlations, and inter-area activity propagation. Simulations of large parts of cortex resolving the individual neurons and synapses enable studying how cortical network structure shapes this multi-scale activity, but have been limited by the available computational resources and simulation technology. Developments in the simulation technology of NEST and access to the JUQUEEN supercomputer have enabled us to overcome this barrier, and simulate a network of 32 vision-related areas of macaque cortex with each area represented by a 1 mm2 microcircuit with the full density of neurons and synapses [1], which avoids distortions due to downscaling [2]. The simulations rely on a recently derived connectivity map for the visual areas of macaque cortex that predicts the connection probability between any two neurons based on their types, areas, and layers [3]. This connectivity map integrates axonal tracing data with predictions from cortical architecture (neuron densities, layer thicknesses), inter-area distances, and neuronal morphologies. In line with models using simplified equations for the individual areas [4], our model predicts that cortex operates in a metastable state where slow activity fluctuations appear. In this regime, the power spectrum of simulated V1 spiking activity and the distribution of spike rates across V1 neurons agree well with those from parallel spike recordings in lightly anesthetized macaque [5]. Furthermore, the inter-area functional connectivity is similar to that from macaque resting-state fMRI [6]. The simulated neuronal activity propagates across areas mainly in the feedback direction, akin to LFP findings during sleep [7]. A mean-field-based analysis [8] shows that the order of activations of the areas is strongly associated with local stability properties, such that the most unstable areas are activated first. Our model reconciles microscopic and macroscopic accounts of cortical neural networks and provides a platform for further developments.
When we learn a new task, changes in our neural activity take place in order to accumulate and act upon relevant information. These changes can appear with different magnitudes in multiple brain areas. To understand the dynamics and ultimately the mechanisms of these changes, we follow mice as they learn to perform a visual change detection task and use wide-field GCaMP signaling to record their neural activity across the dorsal surface of the cortex. We also study random neural network models with cortical-resembling high-level area structures; by iteratively training these networks to perform the task we assess the similarities and differences in the mouse cortex and artificial recurrent networks. We find that initially, during the naïve behavioral stage, the visual cortex alone responds to the changing stimuli. As the learning progresses, frontal areas respond as well, and eventually, at the expert level, the whole mouse cortex responds to taskrelevant stimuli. Cortical activity becomes correlated across all areas, and responses in general become more stereotyped with precise temporal dynamics. Moreover, the dimension of this activity decreases as training progresses. Our artificial neural networks show similar learning-related phenomena. All together, we identify three cortexwide phenomena that emerge during learning of a basic sequential task: task-specific engagement of surprisingly widespread areas across cortex, an increase in the temporal precision and stereotypy of cortical activity, and a reduction of its dimensionality. These phenomena occur in both mouse cortex and in trained, minimally structured artificial neural networks, suggesting that they may recur across many learning systems and posing intriguing questions for further theoretical work.
Supercomputers are important tools for simulation and data processing. High-performance computing (HPC) workflows in neuroscience are applied in several areas: image processing, simulation, visualization, data storage, and parameter space exploration. Execution of workflows as coupled pipelines (Fig. 1) is desirable due to excessive intermediate data, the need for expert human live-interaction and interaction with multiscale components or live experiments. Launching applications on clusters requires batch scheduling, software setup and environment customization. To provide the scientific community with tools for reliable and efficient execution of complex, interactive workflows on HPC systems, we propose "Modular Science" (MS). MS is a software framework and social interaction contract for the deployment of such workflows. At its core lies the orchestration of scientific applications. The modular science orchestrator decomposes the execution of a workflow into process stages, enhancing robustness, debuggability and the probability of success with minimal effort or loss of resources for increasing number of components. MS will enable monitoring of basic workflow variables including data bandwidth, memory consumption and error tracking. Use case 1 includes interactive generation of neural network models using connectivity based on experimental data and executed with the NEST simulator [1]. Case 2 treats a simulation executed in NEST interacting with an Arbor simulation [2] and output processing with Elephant [3] and local field potential calculation [4]. Case 3 includes generation of a full brain simulation using TVB neural mass models (NMM) based on DTI data [5]. For these cases, simulation results must be compared against functional data to iteratively refine the model, interactively guided by an expert. In case 4, we plan to enable a multiscale full brain simulation combining systemic NMMs coupled to specific high-resolution regions in NEST. Firing rate output from the NMMs is coupled to spiking input for the neuron scale simulations, and output is then processed using Elephant for analysis. MS is in early stages of development. Here we present preliminary results from selected cases as a proof-of-concept for complex workflows deployed on the JSC's infrastructure. Our framework is open source and deployable on both local clusters and supercomputing centers and is compatible with most batch systems. This framework will support the reproducibility of large workflows and robust but efficient usage of HPC and data storage resources. Such a framework will be crucial to attack new, large scale neuroscientific problems requiring complex/multiscale workflows combining pluggable simulators and analytical tools. and impairment of a certain circuit across the brain can be translated into malfunction of physiology and behavior, which defines the clinical phenotypes of clinical states. To investigate the impact of PrPSc aggregates on the brain circuit, the amyloidogenic peptide PrP(106-126) derived from PrPC in either aggregated or non-aggregated state was challenged to organ-cultured brain section of wild type mice. The changes occurred in the brain section were monitored at the level of electrophysiology. For the functional connectivity analysis, mutual information were evaluated for each pair of 8 × 8 recording electrodes. For two CSD X = (x_1, x_2, …, x_N) and Y = (y_1, y_2, …, y_N) at each of the two channels under analysis, mutual information(MI) measures the statistical dependence between X and Y [1]. MI is similar in spirit to cross correlation, but is much more general, because MI is capable of capturing the nonlinear dependencies that the cross correlation might have missed. To estimate MI, we use the k nearest-neighbor approach [2]. In our study, Normalized Mutual information (NMI) NMI(X, Y) = (I(X, Y))/ (√(I(X, X)) √(I(X, Y))) is measured to scale the results between 0 and 1. The result indicates that aggregated, but not non-aggregated, PrP(106-126) affected CA2 and CA3 regions of hippocampal circuit with extraordinary activation of CA3 (Fig. 1). The region specific activation implies specific channel expression compared to other areas so pharmacologic study was applied.
Introduction: Voltage-gated sodium currents (VGSCs) are critical mediators of subthreshold oscillations and high frequency burst discharge in sensory neurons. The Nav1.6 type sodium channels mediate a complex form of sodium current that has three distinct components, namely, a transient (INaT), a resurgent (INaR) and, a persistent (INaP) current. Here we developed a unique three-component conductancebased model for this VGSC, and studied its role in the control of burst discharge in proprioceptive sensory neurons in the trigeminal Mesencephalic (Mes V) nucleus. Using model simulations, real-time closedloop dynamic clamp electrophysiology and theoretical analyses of bursting dynamics, we identified a unique push-pull modulation of burst timing by the INaP and INaR components. While INaP enhanced membrane resonance, burst frequency and duration, INaR prolonged the inter-burst intervals and decreased noise sensitivity. We further utilized the model Nav1.6 conductances to restore normal burst patterns using the dynamic-clamp approach, in the abnormally bursting Mes V neurons of a neurodegenerative mouse model. Our results highlight a key role of this VGSC in information processing in these sensory neurons.

Methods:
The Mes V neuronal model follows a conductance-based Hodgkin-Huxley formalism and incorporates a potassium leak current (Ileak), a slow 4-AP sensitive potassium current (IK), and the complex Nav1.6 VGSC with a fast/transient sodium current (INaT), a slowly inactivating persistent sodium current (INaP), and an unconventional resurgent sodium current (INaR). Model simulations were performed using MATLAB; bifurcation analysis was performed using XPPAUT. Real-time dynamic-clamp electrophysiology was performed in brainstem Mes V neurons in P8-P12 mice. Results: Our results using model simulations and dynamic-clamp approach demonstrate distinct effects of INaR versus INaPcomponents in modulating burst intervals in these neurons. While INaPdecreased inter-burst intervals (IBIs), INaRprolonged the IBIs, and, increased the burst refractoriness. The latter was found to be important for correcting burst irregularities simulated by random noise, as well as, in the SOD1G93A mouse model for Amyotrophic Lateral Sclerosis. Theoretical analyses of the steady-state model behavior revealed a push-pull negative feedback control of sodium resurgence by slow channel inactivation, as an underlying biophysical mechanism of burst duration and interval control. Finally, in a simplified monosynaptic reflex network model, we demonstrated that the burst interval control offered by these sodium currents in the Mes V sensory neurons, can modulate postsynaptic trigeminal motor neuron discharge patterns. Conclusions: Our results highlight a novel role for INaRin burst-timing control and noise modulation, with implications to other sensorimotor systems wherein, the Nav1.6 VGSCs underlie ectopic firing/bursting such as in neuropathic pain and spastic/epileptiform discharges. Processing of sequential information is crucial to encoding various inputs from the real world, thus working memory used to store sequential information in humans is of great interest. A number of studies have reported that subjects better memorize the first and last items in a sequence than the others, and refer to this as the primacy and recency effects [1][2]. However, the underlying mechanisms of these effects are still elusive. Here, we propose a novel model of these features of sequential memory performance, by introducing the concepts of sequential overwrite and non-uniform allocation of memory resources. First, for a precise investigation of sequential memory characteristics, we performed experiments where subjects were to memorize a series of visual patterns. As previously reported, we confirmed that the correct ratio of the first and last stimuli in the sequence was higher than those for the others. To explain this result, we modified the standard resource model with the assumption that memory resources of previous information are partially replaced by newly introduced information (overwrite , and that memory performance for each item is proportional to the amount of resources allocated. From this sequential overwrite model with a non-uniform resource allocation from data fit, we could readily explain the observed U-shaped memory performance in human psychophysical experiments. Next, based on our model that sequential overwrite is a key factor of memory performance, we predicted that memory performance can be affected by the modulation of memory overwrite. To test our idea, we designed an experiment where subjects perform sequential memory tasks under three different conditions: correct information, no information, and wrong information about the item numbers are presented before the task. We expected that these three conditions would vary the degree of memory overwrite of the sequential items and that this would affect the performance. As predicted, correct information improved memory performance while wrong information worsened it, and this effect was most significant in earlier sequences where overwrite was stronger. Model parameters fitted to the observed results suggested that the degree of overwrite was significantly different across conditions and well explained the performance. Our model suggests that sequential overwrite and non-uniform allocation of memory resources can explain the origin of the featured U-shape of sequential memory performance. Furthermore, our model suggests a possible mechanism of optimal memory allocation with prior information of the items to memorize. In higher mammals, the primary visual cortex (V1) is organized into various maps of visual functions such as ocular dominance, preferred orientation, and spatial frequency. It has recently been reported that the topography of these functional maps are geometrically correlated [4], such that the contours of the orientation and the spatial frequency maps intersect orthogonally [5]. This may imply an efficient tiling of processing units, but it is still unclear how this systematic organization can develop in the cortex. Here, we introduce our developmental model to suggest that the topography of the functional maps can be seeded altogether from the regularly structured retinal mosaic and that this shared origin results in topographical correlation among the maps. A previous model provides insight into our model, showing that a quasi-periodic orientation map can be seeded by themoiré interference between hexagonal lattices of ON and OFF retinal ganglion cells (RGCs) [6]. The key assumption was that the orientation tuning of a V1 neuron can be predicted bythe localalignmentof ON and OFF RGCs. This is supported by experimental observations that the structure of cortical functional maps is strongly correlated with the local organization of ON and OFF afferents [1]. Expanding this monocular model to binocular condition,we suggest that the local organization of ON and OFF RGCs in the retinal mosaic can also constrain the ocular dominance and spatial frequency preference. We found that the distance between ON and OFF RGCs could determine the separation of ON and OFF receptive field subregions of the connected V1 neuron, and could also change wiring strength to contralateral and ipsilateral feedforward circuits. With the notion that the ipsilateral connections are later developed to match the orientation preference through two pathways, our model showed that the phase difference between contra-and ipsilateral receptive fields of binocular V1 neurons could induce a preferencefor higher spatial frequency than that in the monocular region [4]. As a result, we successfully reconstructed the orthogonal relationships between orientation, ocular dominance, and spatial frequency maps, as observed in the experimental data [2,6]. Our results suggest a unified developmental model of various functional maps in visual cortex.
The basolateral amygdala (BLA) is known to be a core brain region for emotional function, such as fear memory. The thalamic reticular nucleus is the major source of inhibition to the thalamus. However, different thalamic nuclei in the rodent brain receive varying degree of inhibition from local interneurons, ranging from 15 to 20% of the neuronal population in the visual thalamus to < 4% in the somatosensory thalamus [1]. Despite the lower abundance of thalamic interneurons compared to excitatory thalamocortical (TC) cells, they have been shown to shape visual responses and to dynamically influence the extent of receptive fields [2]. As the morphological and electrophysiological properties of TC cells in firstorder thalamic nuclei show high degree of similarity across modalities (e.g., visual and somatosensory systems) [3], we hypothesized that local interneurons in different sensory circuits have similar cellular and synaptic properties, and explored them with the aid of data-driven computational models, in vitro patch-clamp recordings, and gene expression data. We characterized the properties of mice TC neurons of the ventrobasal (VB) nucleus and local interneurons by applying a standardized battery of electrical stimuli, biocytin staining and 3D morphological reconstruction. We qualitatively classified the passive responses and firing properties into different electrical types (e-types) and validated the classification by extracting electrical features from the voltage traces. We then used the 3D morphologies, electrical features, ionic current kinetics and distribution from experimental findings to constrain multi-compartmental models of the different e-types by using a multi-objective optimization strategy [4] and validated them with stimuli not used during model building. We complemented our data analysis and modelling pipeline by comparing the modelled e-types with single cell and synaptic properties systematically curated from the neuroscientific literature [5], along with gene expression data. The result of this analysis suggests that while some thalamic e-types are comparable to interneurons in cortical microcircuits, others are thalamus-specific and comparable to interneurons from the dorsal part of the lateral geniculate nucleus [6]. We investigate the role of synaptic connectivity in the cortical network that may lead to a better understanding of autism spectrum disorder (ASD). We have established a measure and studied properties of the integration of distinct stimuli in the cortical network model and investigated effects of network connectivity on this integration. Taken with ongoing experimental optogenetic studies [1], this model may assist to pave the way toward potential pharmacological targets for the treatment of ASD. Autism is a neurodevelopmental disorder for which there is no cure. It is characterized by impairments in social cognition and communication. Abnormalities in the ASD brain are not strictly localized, but involve multiple neural networks [2,3]. Numerous studies suggest a deficit in multisensory integration (MSI) in autism, both in human and animal models [4,5]. Specifically, it has been shown that ASD patients demonstrate a widened window of audio-visual temporal integration [6]. It has been proposed that deficits in the integration of multisensory cues leading to ASD are likely to be caused by dysfunctional connectivity in the brain [6,7,8].
Multisensory integration by neural populations in the cortex has been recently studied in many contexts [9,10,11]. Our model is an adaptation of the Traub model [12] of a single-column thalamocortical network, modified to be used in the GENESIS neuronal simulation environment [13]. The model comprises 14 types of cortical neurons each with its own compartmental morphology and electrophysiological properties connected in columnar structure. We apply pulse train stimuli to the cells at two different locations in the column and measure the local field potential (LFP) to characterize network activity. Our model demonstrates that the multiple distinct stimuli generate superadditive integrated LFP response, where the combined stimuli from two locations produces a larger response than the sum of the two individual ones. We then use this model to investigate the effect of network connectional parameters on temporal aspects of multisensory integration in the cortex. It is believed that the ASD condition may be associated with widened temporal windows of cortical integration. Existing ASD therapies concentrate on behavioral interventions that reduce symptoms and, to date, no drug therapy exists for ASD that would repair or strengthen brain circuits. With our model's measure of MSI and cholinergic impairment, it may be possible to gauge the efficacy of pharmacological agents whose action would ameliorate the ASD condition. Synaptic decoding of neural population activity at the single cell level presents a challenging question. One method to address this question in a rigorous way is to use detailed single neuron simulations with well-defined input patterns to study the input-output function of biophysically realistic neurons. In previous work we developed method to create artificial spike trains (ASTs) that can match spike train properties of cerebellar Purkinje cells in order to study the cerebellar corticalnuclear signal transformation [1]. Here we generalize this method to create well defined artificial spike trains (ASTs) made from templates of different types of recorded neurons and further test the method with surrogate data. The basic idea of our method is to use recorded neurons to construct rate templates of their activity using gaussians. Then we can draw gamma distributed spike trains from these rate templates to obtain ASTs with different regularity properties. We can scale templates to different firing rates, add a refractory period to the gamma distributions, and add well defined rate-correlations between multiple ASTs by manipulating the rate template. Here we first tested our method with constant rate templates, sinusoidal rate templates, and zap rate templates. We find that slow rate fluctuations (~ 1 Hz) can be well captured by individual ASTs, but that faster rate fluctuations require a population average of ASTs to recapture the rate template. The ability to capture faster rate fluctuations is a function of the regularity (kappa parameter of the gamma distribution) and the rate of the ASTs that are being generated. These properties parameterize fundamental limits of coding rate fluctuations with noisy spike trains. We then use pyramidal neuron and mossy fiber recordings from the cerebellar to test our algorithms for real data beyond the fast firing Purkinje cell populations previously used. Unlike cerebellar Purkinje cells which exhibit high firing rate and more regular spike trains, most pyramidal neurons and mossy fibers exhibit low firing rates and highly irregular and bursty firing pattern. In spite of these differences in firing rate and pattern, ASTs created using our method were able to match the statistical properties of spike trains in both these cell types. The ability to re-create original rate templates from such ASTs was limited by the same features as seen for our surrogate rate templates, and reveal limitations as to how faithfully low rate/high variability spike trains can communicate a rate code. Understanding the contributions of different layers to cortical processing within the visual and somatosensory systems has led to testable hypotheses about how multiple simple sensory features are combined within cortical columns of these areas. This understanding provides a concrete basis for bridging high-level models of information processing with proposed neurobiological components and recordings of neuronal activity. To complement these studies, insight into cortical oscillations provides a convenient framework for understanding population processing of sensory features. In particular, auditory cortex demonstrates an organizational hierarchy of rhythmic activity and these rhythms effect stimulus encoding [1]. Previous studies further implicate phase-resetting and reciprocal interlaminar interactions in feature selection [2,3]. Large-scale computational modeling studies of sensory systems have focused primarily on visual and somatosensory cortices. Despite recent advancements in characterizing anatomical and physiological properties of the auditory system, few models reconstruct fundamental differences between the auditory system and other modalities. Here we present a model that incorporates some features unique to auditory cortex and replicates various statistical response properties of a non-primary auditory area [4]. The model is a modification of a biologically realistic model of a thalamocortical network with multiple layers [5] converted for broader usage [6]. Simulations were performed within the NEURON simulation environment [7] using NetPyNe [8] for model handling and analysis.
The model provides predictions about how convergent inputs carrying distinct information are processed within a thalamocortical network and demonstrate relationships between laminar cortical oscillations and interlaminar processing of convergent inputs. We also investigate how phase-resetting reorganizes laminar processing and affects feedforward and feedback outputs. Our model simulations assume auditory processing relies on a hierarchical network structure. Feedforward inputs from earlier auditory populations (in area A1) are provided to the network based on statistical response patterns of multi-unit activity to a repertoire of auditory stimulation found in the same literature used to replicate non-primary auditory responses [4]. The model is modified and tuned to replicate the response properties of later areas under similar stimulation protocols. As the model uses a biophysically based multicompartmental formalism, we demonstrate how convergent extrinsic inputs shape local population activity reflected in both population firing and local field potentials. Biophysically realistic computational models are ideally suited for simulation of predictions that can be verified experimentally via chemical, surgical, and optogenetic manipulations. NeuroML is a modular, declarative, simulator-independent model language for describing and exchanging such models. NeuroML-DB.org is an online resource for rapidly locating NeuroML models by keyword search or ontological (anatomical or neurotransmitter) relations [1]. Rapidly locating models with desired dynamical and computational properties is an ongoing challenge. While NEURON simulator models can be examined via ModelDB, NEURON-based channel model dynamical properties can be viewed via ICGenealogy project [2], and NeuroML based models can be shared on the OpenSourceBrain.org website, a resource that allows viewing of systematic evaluations of NeuroML models has not been available. We have extended our previous work on NeuroML-DB and enabled interactive, visual inspection of temporal dynamics of channel models. Channel activation, deactivation, and inactivation voltage clamp protocols similar to those used in the ICGenealogy project were utilized to record the current and conductance responses of 133 NeuroML channel models available in NeuroML-DB. The resulting web application, accessible viaNeuroML-DB.org, allows visual inspection of, including zooming and panning, the temporal dynamics of channel currents and conductances during channel activation, inactivation, and deactivation voltage clamp protocols (Fig. 1, model from [3]). For channels with [Ca2+] dependent dynamics, the effect of [Ca2+] on the channel output response can be interactively examined in the interface. Additionally, files that reproduce the plots of each channel response can be downloaded via the interface. Using the jNeuroML library, the files can be converted to several other simulator formats. Finally, channel detail web pages display notes, if any, regarding numerical or simulation issues encountered while running the simulations. A channel model implemented in NeuroML can be converted to a wide range of formats and included in larger cell and network models, regardless of the original simulator used to implement the channel model. This work allows modelers to rapidly locate and visually inspect the dynamical properties of NeuroML channels and assess their suitability for inclusion in larger models. In ongoing work, we are utilizing the cell stimulation protocols [4] from the Allen Brain Atlas project to characterize the responses of NeuroML cell models. The cell protocols include ramp, step, threshold, and pink noise current injections and assess properties such as resting voltage, rheobase and threshold currents for all NeuroML cell models in NeuroML-DB, as well as quantify model run time complexity. Computational models of biological systems are rarely formally tested for agreement between model output and experimental data. SciUnit, a software framework for model validation, facilitates such rigorous testing. During model development, models can be continuously subjected to data-driven "unit tests" that quantitatively summarize model-data agreement, identifying modeling progress and highlighting output that fails to adequately reproduce observed data from the corresponding biological system. The OpenWorm Project is an international open-source collaboration to create a multiscale model of the organism C. elegans. At every scale, including subcellular, cellular, network, and behavior, this project employs one or more computational models that aim to recapitulate the corresponding biological system at that scale. This requires that the simulated behavior of each model be compared to experimental data both as the model is continuously refined and as new experimental data become available. We present data-driven OpenWorm model validation usingSciUnitat three model scales: 1) ion channels, 2) neurons, and 3) whole organism motor output. This workflow is publicly visible and accepts community contributions to ensure that modeling goals are transparent and well-informed. Model validation tests are executed continuously as the models are updated and refined, ensuring that development converges towards the ultimate design specification: agreement with the underlying biological system. Perceptual rivalry is the subjective experience of alternations between competing percepts when an individual is presented with an ambiguous stimulus. Mutual inhibition between pools of neurons encoding different interpretations of the stimulus is thought to underlie this disambiguation computation, where activity in one pool dominates and the corresponding percept is represented. A canonical cortical circuit model with mutual inhibition and fatigue can explain normalization, winner-take-all, rivalry, and various findings from psychophysics experiments [1]. However, this approach has yet to incorporate realistic spiking statistics. Meanwhile, balanced state theory has been used to explain why cortical neurons fire irregularly [2]. Researchers have modeled computing orientation selectivity and working memory in balanced networks, but competitive networks have only recently been investigated [3]. A recent study showed that alternations resembling rivalry result from random networks receiving stochastic but competitive inputs [4]. Here we explore a model of perceptual rivalry with realistic spiking. First, we show that normalization, winner-take-all, and rivalry behaviors can coexist with a realistic asynchronous-irregular state. Next, we compare the psychophysical properties of this model to those of a random network with stochastic input. Our model can explain Levelt's second and fourth propositions, a gamma distribution of dominance times, and can maintain a coefficient of variation of dominance times which is stable across changes in the input. By contrast, a random network cannot explain Levelt's fourth proposition, and does not reproduce a gamma distribution of dominance times.

Multiscale model validation with SciUnit
Thoracic sympathetic postganglionic neurons (tSPNs), innervated by preganglionic neurons in the spinal cord, are the last common motor output of the sympathetic nervous system, and directly control the vasculature and other internal organs. Dysfunction of tSPNs, such as hyperexcitability, has been observed after spinal cord injury, yet little is known about the cellular mechanisms that drive the excitability of tSPNs.
Combining electrophysiological data with computational modeling, we built the first physiologically-realistic single neuron model of tSPN in mice, and elucidated several cellular mechanisms that govern the tSPN dynamics. For example, we found that the postinhibitory rebound that has been observed in tSPNs ex vivo was induced by the sodium and potassium currents (INaand IKd) instead of the T-type calcium current. We also found that both the M-type potassium current (IM) and the calcium-dependent potassium current (IKCa) were necessary to replicate the spike rate adaptation. Together, we reproduced all the essential features of tSPNs ex vivo with eight types of ionic currents, which are INa, IKd, IM, IKCa, a fast transient potassium current (IA), a persistent calcium current (ICaL), a hyperpolarization-activated inward current (Ih), and a leak current (IL). Using this single neuron model, we employed an ensemble modeling approach to build a database of physiologically-realistic tSPN models [1][2][3], which enables a more comprehensive and rigorous examination of the range of tSPN responses to various synaptic inputs. Overall, this work lays the foundation to examine both the recruitment principles of synaptic inputs at tSPNs and the dysfunction of tSPNs after spinal cord injury in the future.
The mammalian visual cortex is composed of multiple areas that are organized in a hierarchical structure, with feedforward, feedback and horizontal connections. The bottom-up convergent connections generate larger spatial receptive fields and longer temporal integration windows at higher levels of the visual hierarchy. While the structure of these anatomical connections is relatively fixed on short timescales, the structure of functional interactions can rapidly change conformation due to changes in external stimuli and internal brain states. This flexibility is critical for selectively routing signals for perception, cognition, and behavior. Therefore, understanding how neurons form functional networks is fundamental for deciphering brain functions. Studies that investigate functional networks with resting-state fMRI can image the entire brain at once, but, due to the low temporal resolution of this method, they fail to uncover network dynamics at fast timescales that are important for many aspects of perception and decision-making. Studies that attempt to measure functional connectivity with electrophysiological recordings have typically been restricted to recordings from two brain areas at a time, with a limited number of simultaneously neurons in each dataset. Therefore, novel methods are needed for recording large populations of neurons with sufficient temporal resolution to study dynamic functional connectivity at a large scale. Here we make use of the newly developed Neuropixels probe, which contains 384 densely arranged recording sites along a linear shank. We built a platform to simultaneously record from 6 independent probes inserted in the mouse visual cortical areas including primary visual cortex (V1) and 5 higher-order visual cortical areas (LM, RL, AL, PM, and AM). The linear probes are inserted across the layers of the cortex in head-fixed awake mice. The high yields of the Neuropixels probe allow us to record simultaneously from more than 700 well-isolated neurons distributed across cortical layers and areas in the visual cortex of a single mouse. To maximize the probability of finding mono-synaptic functional connections, we mapped the retinotopy of each area with intrinsic signal imaging and specifically target regions with overlapping visual fields. This targeting is validated by receptive field mapping. To compare functional networks under the context of different sensory inputs, we studied activity during drifting gratings and natural movies, in addition to comparing with spontaneous nonstimulus driven activity. We used two methods to measure functional connectivity within the visual cortical network: fine timescale pairwise cross-correlogram (CCG) analysis and Granger causality analysis, both of which can reveal functional relationships of recorded neurons. We found that both the proportion of effective connections and the strength of the functional connection decay as a function of receptive field separation. The time delay of effective connections revealed layer-specific functional sub-networks, based on cortical layers estimated from current source density analysis. We also observed significant differences in functional connectivity between gratings, movies, and spontaneous activity. In sum, our platform provides a unique opportunity to directly study millisecond-timescale functional networks across 6 highly interconnected cortical areas.

P221
Building individualized dynamic brain models at high spatial resolution using In the current era of big data and ever improving methods for data acquisition, several initiatives have produced innovative whole-brain dynamical systems models of neural activity. Unfortunately, there remains a gap in (large-scale) systems neuroscience communities between generative models of brain dynamics vs. descriptive or statistical models of brain activity. The latter has proven highly amenable to data-driven characterizations of individual differences, while the former allows for the overt construction of mechanistic hypotheses regarding circuit function. In the present work we attempt to bridge this gap and present a procedure to fit independent large-scale dynamical-systems models to individual subjects using fMRI. Through a novel fitting procedure, we parameterize networks of hundreds to thousands of nonlinear neural masses (brain regions) each with their own intrinsic dynamics and connections for individual human subjects using resting state fMRI. We rigorously assess psychometric properties and sensitivity to preprocessing choices and a variety of nuissance factors. We perform multiple validity and reliability tests. Results demonstrate that our procedure is robust, reliable, and (for the cases tested) valid. We demonstrate that these models provide novel insight into both large-scale and local dynamics including a "hierarchy of intrinsic time scales" and spatial organization in a previously unexplored feature: transfer function curvature. By forward simulating individual subject's models we recreate individual subject's "dynamic functional connectivity" (dFC) patterns as well as higher-order phenomena such as dwell-time and transition probability between dFC patterns despite using time-invariant models ( Fig. 1).

P222
Revisiting efficient coding of natural sounds in the environment: unsupervised learning or task-based optimization? Efficient coding has been a leading computational principle for the sensory neuroscience. Following the visual system, Lewicki [1] argued that the auditory periphery can be explained by unsupervised learning of natural sounds. One of the study's claims is that the basis optimized to code human voice resembles the auditory nerve fibres, whose filter sharpness distribution is preserved across mammals. We were able to reproduce the matched distribution by applying the same algorithm to clean recordings of human voices. However, we also found that an efficient code for human voices recorded in the natural environment shows much sharper tuning than the auditory nerve fibres, even though the environment recording is closer to the sensory signal our ears receive and more natural than a studio recording. Our analysis showed that the waveforms are distorted on a short time scale comparable to the time window of the auditory nerve filters and that the mismatch can be reproduced by simulating environmental reverberations, suggesting that the primary factor in the mismatch is the environmental reverberations. How can we better model the auditory periphery including the environmental modulations? Inspired by a recent work on the visual hierarchy [2], we hypothesized that the auditory periphery is optimized to perform auditory tasks we face in the natural environment instead of unsupervised learning of the entire incoming signal. To test this, we built a deep convolutional neural network that receives reverberated waveform inputs. As a naturalistic task related to voices on a short time scale comparable to the time window of the auditory periphery, we chose phoneme classification. Waveforms and phone labels were taken from the TIMIT database. Each input had the length of 2000 data points, with the target phoneme at the centre. The input waveforms were convolved with an impulse response randomly chosen from a database [3]. After the training, the waveform filters learned in the first layer showed characteristics similar to the auditory nerve fibres, whereas a normal efficient code for the same input did not. This result does not depend on the speech dataset, since we could reproduce a qualitatively similar result by applying a different task, that is, classification of environmental sound recordings. Overall, the results suggest that the auditory periphery efficiently encodes taskrelated information in a reverberation-resistant manner rather than the entire incoming signal and that our understanding of sensory systems in the natural environment, not only the visual system, can be furthered by using a framework of task-based optimization.

P223
Emergence of auditory-system-like representation of amplitude modulation in a deep neural network trained for sound classification Takuya  Amplitude modulation (AM) of a sound contains essential information for auditory perception, and neural representation of AM has been targeted by a number of neurophysiological studies. These studies generally indicate that neurons in the auditory nervous system (ANS) exhibit tuning to AM frequency, usually characterized by the synchronization to the stimulus AM and the average spike rate during the stimulus . A model is designed to be optimal for reproducing a function of the sensory nervous system without assuming its anatomical or physiological properties. Thus characteristics emerging in the optimized model should reflect only the nature of the task and the input data. Especially, a deep neural network (DNN) is suitable for this purpose because it can take raw data as inputs and process them for behaviourally relevant tasks (Yamins, et al., 2014, PNAS). The present study attempts to understand the functional significance of the AM representation in the ANS by modelling the ANS with a DNN trained for sound classification. To directly compare the AM representation in the DNN with that in the ANS, we applied neurophysiological methods for analysing the DNN, and characterized the AM tuning of the units in each layer of the DNN. Similarly to the ANS, most units in the trained DNN exhibited low-pass or band-pass AM tuning. The activities of the units in the lower layers synchronized to high AM frequency, whereas those in the higher layers synchronized to low AM frequency. The average activities exhibited AM tuning only in the higher layers. We calculated the pairwise similarities of the layers in the DNN and the regions in the ANS, and found that the tunings in the lower/higher layers in the DNN were similar to the tunings in the peripheral/central regions in the ANS (Fig. 1). Such characteristics were not observed in the DNN before optimization, but gradually emerged during the optimization. By conducting the same analysis in the DNNs with different architectures, we showed that similarity to the ANS significantly correlated with the classification accuracy of the DNN. These results suggest that the AM representation in the ANS also emerged during optimization to natural sound recognition. This work was supported by JSPS KAKENHI Grant Number JP15H05915.

Acknowledgments
Grant-in-Aid for Scientific Research on Innovative Areas "Innovative SHITSUKSAN Science and Technology.

P224
Reproducing the cognitive function with the robustness against the brain structure and with the efficient learning algorithm Yoshihisa Fujita, Shin Ishii Kyoto University, Graduate School of Informatics, Kyoto, Japan Correspondence: Yoshihisa Fujita (fujita-y@sys.i.kyoto-u.ac.jp) BMC Neuroscience 2018, 19(suppl 2):P224 Computational modeling of biological neural networks which generate cognitive functions has difficulties compared to artificial neural networks. One major problem is the individual variability of the brain structure which includes extensively different connectivity patterns of neurons for the same function. Another problem is that the widely-used learning algorithms in artificial neural networks such as error backpropagation have high computational costs and are not biologically plausible. Therefore, the models of biological neural networks require both robustness against the network structure and an efficient, plausible learning algorithm. How to achieve highly cognitive functions under these constraints is largely unknown. To tackle this issue, we developed a neural network model based on Extreme Learning Machine (ELM), which includes random and fixed connections. We assumed that this randomness corresponds to the structural variability. ELM utilizes its random connections and thus its learning algorithm is quite simple. Using ELM, we tried to implement the function to recognize words from a string of letters. Since this function requires the process to integrate letters into words, we adopt Vector Symbolic Architecture (VSA), in which the patterns of neural activities for recognized objects are expressed as binary vectors and the integrating process is expressed as a vector operation. We developed a new learning model combining ELM and VSA. Unlike ordinary ELM, the balance between excitatory and inhibitory neurons was crucial in our model, which is consistent with biological findings. Within this balanced condition, it could successfully learn the vocabulary with low computational costs. We used this model to examine how neural representations of misspellings are different from those of correct spellings. Our model can provide a clue how the brain achieves cognitive functions efficiently under the structural variability. Common marmosets (Callithrix jacchus) are non-human primates with a small brain size (~ 8 g), mature quickly, and can be genetically manipulated. These features make them suitable for understanding the complex structure and connectome of the primate brain, and provide another animal for inter-species comparison. The complexities of the primate are also shown through individual structure and connectome differences. The RIKEN BSI team obtain diffusion magnetic resonance imaging (dMRI) data of a large number of marmoset brains (~ 50) to construct macro-scale connectomes and reveal the individuality between brains. Using dMRI, we can observe water diffusion phenomena within small voxels (200 μm isotropic) across the brain and use this to estimate the direction of fiber bundles through voxels. Observing the whole brain with dMRI and connecting these estimated fiber directions between adjoining voxels allows us to estimate long-range fiber bundles through the white matter of the brain. Here, we propose a pipeline for obtaining the macro-scale connectome of RIKEN's marmoset data and analyze their individual differences of the structure and connectome. The pipeline mainly consists of two components: (1) global reconstruction of fibers, and (2) improved parcellation of brain regions by a deep learning technique. (1) Because the fiber structure of the marmoset brain is complex, we often observe multiple and different directions of fiber bundles through a voxel of dMRI. We focus on the Bayesian method for global fiber reconstruction using dMRI (Reisert, et. al. 2011), which can successfully distinguish these different directions in Human dMRI analysis. (2) To make the connectivity matrix between brain regions, we parcellated individual brain regions and counted the number of fiber bundles between ROIs. We developed a structural registration method with the help of a recent deep learning technique for parcellating the white and gray matter (VoxResNet). We estimated the cortex, white matter, subcortical regions and cerebellum in the individual brains using this method. This rough parcellation of the brain can enhance the precise correspondence between the individual and a standard brain with atlas, in which a detailed anatomical structure was defined. We analyzed RIKEN's data with our pipeline and obtained the individual and average connectivity matrix of the marmoset brain. The graph structure of the connectivity matrices of marmosets show an exponential decay rule (EDR), which states that the strength of connectivity between two regions decays exponentially with respect to the distance between them. This rule is consistent with the graph structure of the connectivity matrices found in studies of other species (mouse and macaque monkey). We evaluated and visualized the individuality of marmoset connectivity matrices using t-SNE; a low-dimensional manifold mapping method. We found that we could discriminate between all individuals through their connectivity matrices, despite variations in the experimental settings of the dMRI and b-values. We can recognize rapidly and effortlessly complex visual scenes. Such amazing ability in visual recognition needs the effective processing of visual information along the multiple stages of visual pathways. Neurophysiological experiments have provided evidence for a "simple-to-complex" processing model based on a hierarchy of increasing complex image features, performed along the feedforward pathway of the ventral visual system. On the other hand, visual system has abundant feedback connections, whose number is even larger than the feedforward ones. Li et al. [1] showed that top-down signals allowed neurons of the primary visual cortex (V1) to engage stimulus components that are relevant to a perceptional task and to discard influences from components that are irrelevant to the task. They showed that V1 neurons exhibited tuning curves modulated depending on the task context. We demonstrated what kinds of top-down signals generate the tuning curves [2]. However, it remains unclear how top-down signals reflecting task behaviors emerge and how they modulate the tuning curves of V1 neurons. To address this issue, we develop a model of visual system that consists of networks of V1 a higher visual area, and a recognition area. We consider one of the perceptual tasks used by Li et al., or bisection task. Neurons of the higher visual area receive the top-down signal reflecting a decision of the task, as well as the feedforward inputs from V1 neuron, and feed the outputs back to V1 neurons. We also use a reinforcement learning to acquire an adaptive behavior to the task. The synaptic weights from neurons of the higher visual area to those in the recognition area were determined by Mirror Decent (MD) method [3] on the basis of error rate of behavior. The synaptic weights involved in top-down signals in three areas were shaped by Hebbian learning, concurrently with the learning of the feedforward connections. We show here how the feedforward and feedback connections involved in the three areas are formed by the learning of adaptive behavior. Top-down signals are generated concurrently with the acquisition of adaptive behavior. We also show that the tuning modulations of V1 neurons are caused by the change in activity through long-range connections of V1 neurons, elicited by top-down signals from the recognition area to V1 via the higher area. Our model provides the results on the tuning properties of V1 neurons that are compatible with the experimental results by Li et al. These results provide insights into understanding how behavior affects information processing in early sensory areas. Computational models in neuroscience typically contain a number of parameters that are uncertain, either because they vary between cells or dynamically within a cell, or because they are difficult to measure accurately. Uncertainty quantification is a means to quantify the uncertainty in the model output that arise from uncertainty in the model parameters, while sensitivity analysis is the process of quantifying how much of the output uncertainty each parameter is responsible for. Unfortunately, uncertainty quantification and sensitivity analysis are not standard practices in the field of neuroscience, and models are commonly presented without any form of uncertainty quantification. To help alleviate this problem we have created Uncertainpy (https ://githu b.com/simet enn/uncer tainp y), an open-source Python toolbox, tailored to perform uncertainty quantification and sensitivity analysis of neuroscience models. Uncertainpy aims to make it easy for users to perform uncertainty quantification and sensitivity analysis without requiring detailed prior knowledge. The toolbox allows uncertainty quantification and sensitivity analysis to be performed on already existing models, and does not require changes to be made to the model implementation. Uncertainpy primarily bases its analysis on polynomial chaos expansions [1], which are faster than the more standard Monte-Carlo based approaches. Polynomial Chaos expansions are obtained from the previously developed package Chaospy [2]. Uncertainpy does not merely perform an uncertainty analysis of the "raw" model output (e.g. membrane voltage traces), but is tailored for neuroscience applications by an built-in capability of identifying characteristic features in the model output. Uncertainpy then performs an uncertainty analysis of these features. For example, the toolbox can quantify the uncertainty and sensitivity of salient model response features such as spike timing, action potential width, mean interspike interval, and other features relevant for various neural and neural network models. Uncertainpy comes with several common neuroscience models and features built in, and including custom models and new features is easy. We here present Uncertainpy, and demonstrate its broad applicability by performing an uncertainty quantification and sensitivity analysis of three case studies relevant for neuroscience: the original Hodgkin-Huxley point-neuron model [3], a multi-compartmental model of a thalamic interneuron [4] implemented in the NEURON simulator, and a sparsely connected recurrent network model [5] implemented in the NEST simulator. A preprint of this work is available at bioRxiv [6]. The feed forward propagation of visual information in rate-coded neural networks discards information about which low-level features are driving high-level transformation invariant features. In particular, when multiple visual stimuli are present, the network has no way of assigning which low-level features belong to high-level features. Information about the configuration of high-level features, i.e. the spatial composition of the high-level feature in terms of its low-level features, is therefore lost. In visual psychology, this is known as the feature binding problem. Capsule networks [1] demonstrate a new type of artificial neuron called capsules. The activity of a capsule is represented by a vector rather than a "firing rate" value. The magnitude of the vector represents the probability that its preferred feature is present, whilst the direction of the vector represents the configuration of the feature transform. Capsules therefore provide a simultaneous representation of the presence of a feature and its configuration in terms of lower level features, somewhat analogously to feature binding. But these are not a plausible model of brain function. We show how a spiking neural network can give rise to emergent spatiotemporal spike patterns and feature binding representations. An earlier modelling study [2] showed how synchronized activity can emerge over a series of layers. Building on this work, we show how incorporating randomized axonal delays leads to the emergence of spatiotemporal patterns of spikes (polychronization). This is an inductive process over a series of layers. Such spike patterns emerge even when the input neurons have randomized spike times. These spatiotemporal spike patterns carry information relating to the hierarchical binding relations between lower and higher features. Our simulations demonstrate that neurons can learn to respond invariantly over a range of transformations of a highlevel feature, whilst simultaneously representing the configuration of the feature transform using the precise timings of their spikes. The relative timings between the spikes of neurons representing a high-level feature vary continuously and monotonically as the high-level feature undergoes transformation. Such a representation may be how the brain forms a similar representation to that of capsules and could be part of the brain's solution to the feature-binding problem. Our spiking network models can also represent the hierarchical binding relations between lower and higher level features through the emergence of binding neurons [3], which fire if and only if a neuron encoding a lower-level feature is participating in firing a neuron representing a higher-level feature. This implies that the low level feature is part of the high level feature. Such binding neurons develop through visually guided learning with STDP. Excitatory plasticity has long been the focus of learning in spiking neural networks. From the earliest pairwise Spike-Timing Dependent Plasticity (STDP) rules to triplet STDP rules and beyond, excitatory learning rules have been explored both experimentally and theoretically. Inhibitory plasticity has only more recently been appreciated and shows promise in network stabilisation, homeostasis, and predictive coding [1,2]. Very recently, investigations into the properties of inhibitory plasticity in decorrelating excitatory responses and shaping excitatory synaptic weights have emerged in both experimental and theoretical studies [3,4]. This study aims to highlight these effects and make the claim that under a given inhibitory plasticity rule [1], inhibitory spatio-temporal receptive fields play a crucial role in the development of excitatory receptive field structures. Decorrelation of excitatory cells has been identified as a function of inhibitory neurons and plasticity [3]. This description can be misleading. In fact, excitatory neurons in a network with a correlative inhibitory plasticity rule have activity which is decorrelated with their inhibitory input activity. The result is two-fold. First, excitatory receptive fields can only cover those stimuli which are not predicted/explained by their inhibitory inputs. In the case of balance, inhibitory inputs can "explain away" all incoming excitatory stimulation. However, if this balance is incomplete, it leads to excitatory cells responding to only those stimuli which its inhibitory inputs do not cover. The connectivity (and receptive field tuning widths) of inhibitory cells therefore determine the corresponding excitatory cell response characteristics. Curious effects also occur if excitatory and inhibitory cells are active on different input time scales. This is most clear when we consider dynamic stimuli with inhibitory neurons and inhibitory synapses acting on a timescale significantly faster or slower than excitatory cells/synapses. Under these conditions, excitatory cells compete to form a receptive field on the timescale of the incoming inhibition. The excitatory receptive field forms a peak at this timescale (e.g. close to or far from stimulus onset) and learning at all other timescales reflect features that predict activation at that peak. These effects are studied in a spiking neural network model. In particular, this study shows the emergence of structure in the inhibitory weights (related to the correlation of pre and post synaptic cell responses) and how this affects the emergence of excitatory spatiotemporal receptive fields.

Inhibitory plasticity moulding excitatory spatio-temporal receptive fields in a spiking neural network model
Neurons in the brain are connected in intricate arrangements, communicating with each other mostly through chemical synapses. The collective dynamics of neuronal networks-responsible for the cognitive and other functions of the brain-is strongly influenced by their connection topology. Motivated by reports [1,2] that at least some parts of the brain may be wired in a modular fashion, we first show that the occurrence of modularity promotes the dynamical balance between excitation and inhibition in neuronal networks. In other words, the range of values for the fraction of inhibitory neurons in the network that will result in a high probability of persistent activity is considerably amplified by modular organization of the connection topology. We then explore the possibility that this modularity arises in a self-organized manner when a set of neurons interacting with each other spontaneously organize themselves into densely inter-connected communities with sparser connections to other communities. For this purpose, in parallel with action potential generating models of individual neurons, we investigate spike-time dependent plasticity (STDP) synaptic dynamics [3] for evolving the connectivity organization between neurons. Note that STDP is considered to be the principal mechanism by which neurons strengthen or weaken connections, and thereby "learn" according to the Hebbian paradigm ("neurons that fire together, wire together"). We show that, beginning from a homogeneous globally connected network, the interplay of neuronal dynamics and STDP results in the spontaneous emergence of modules in the network. Furthermore, along with synapses, neurons are also connected through gap junctions in real systems. The precise role played by the gap junctions in a network is not yet fully understood, although it is believed to be important for synchronization and rhythmic activity. We consider a system of neurons with two distinct types of interactions (synapse and gap junctions) and such systems can be better described using a multiplex network [4] (Fig. 1). Using this network paradigm, we aim to uncover new emergent dynamical behaviors that arise due to coexistence of two types of connections and get further insights into how gap-junctions enhance the robustness and efficiency of real neuronal systems in general.

P231
Uncovering the mesoscopic organisation of the macaque brain Anand Pathak, Shakti N. Menon, Sitabhra Sinha The Institute of Mathematical Sciences, Theoretical Physics, Chennai, India Correspondence: Anand Pathak (anandpathak31@gmail.com) BMC Neuroscience 2018, 19(suppl 2):P231 In order to unravel how the brains of higher organisms carry out cognitive and motor functions, it is crucial to understand the structural organization of the neurons at different levels of hierarchy. At an anatomical level, the mammalian brain is compartmentalized into different regions (lobes, gyri, nuclei etc.) Brain regions, each having millions of neurons, are connected to each other through axonal bundles projecting out from their respective neurons. One approach to study the structural connectome of the brain is to consider the brain network at the scale of brain regions. This approach has become possible to implement in brains of higher mammals, e.g., the macaque monkey, for which a large amount of data on different brain areas and their connectivity have been collated in the online repository CoCoMac. Using as our starting point a previous study [1] that organized the data available from the CoCoMac database of macaque brain connectivity, we have reconstructed an unambiguous and comprehensive brain network with regions covering the entire brain cortex as well as subcortical regions (Fig. 1). Mesoscopic analysis of a network pertains to those substructures that occur at a level much higher than few nodes but at a scale lower than whole network. In a brain network, knowing the mesoscopic organization can characterize its basic structural and functional make up. In particular we consider the modular and hierarchical organization of the macaque brain. Modular networks consist of sub-networks called modules that are densely connected within themselves but have sparse inter-modular connections. A high degree of structural modularity in a brain network reveals functional compartmentalisation Fig. 1 Schematic representation of the neuronal system as studied in this work and co-ordination among the brain regions. A hierarchical network on the other hand is organized into layers that have dense connections between consecutive layers and relatively sparse connections between non-consecutive layers. Uncovering the hierarchical organization of a brain network not only illuminates the structural and functional compartmentalization but also the directionality of information flow. Analysing these two distinct types of mesoscopic organization of the macaque brain network reveals the larger plan of the macaque brain architecture. Our analysis reveals that the macaque brain network exhibits a highly modular structure, which in spite of being in accordance with the known functional organization, provides interesting new insights and implications about the functioning of many less studied brain regions. The arrangement of modules is spatially contiguous except one saliently fragmented module that is particularly intriguing and suggestive in its functionality. This study strongly suggests that even though macaque brain connectivity is clearly governed by the spatial configuration of the brain regions, the modular structure is essentially independent of the geometry and shape of the brain, and hence the emergence of modules seems to be a more fundamental attribute. An even more surprising observation is that the macaque brain network also has a highly hierarchical structure. Our new original approach for determining the hidden hierarchical structures in a network opens up a whole range of possible analysis to further understand the structure and function of neuronal networks in general. Recordings of extracellular electrical, and later also magnetic, brain signals have been the dominant technique for measuring brain activity for decades. The interpretation of such signals is however nontrivial [1], as the measured signals result from both local and distant neuronal activity. In volume-conductor theory the extracellular potentials can be calculated from a distance-weighted sum of contributions from transmembrane currents of neurons. Further, given the same transmembrane currents, the contributions to the magnetic field recorded both inside and outside the brain can also be computed [2]. This allows for the development of computational tools implementing forward models grounded in the biophysics underlying the different measurement modalities [1]. LFPy ( [3], LFPy.github.io) incorporated a now well-established scheme for predicting extracellular potentials of individual neurons with arbitrary levels of biological detail. It relies on NEURON ( [4], neuron.yale.edu) to compute transmembrane currents of multicompartment neurons which is then used in conjunction with an electrostatic forward model [5]. We have now extended its functionality to populations and networks of multicompartment neurons with concurrent calculations of extracellular potentials and current dipole moments. The current-dipole moments are used to compute non-invasive measures of neuronal activity, like magnetoencephalographic (MEG) signals [2,6] and, when combined with an appropriate head-model, electroencephalogram (EEG) scalp potentials. One such built-in head-model is the 4-sphere model including the different electric conductivities of brain, cerebral spinal fluid, skull and scalp [6,7]. Modeling studies of cortical network dynamics frequently aim to include realistic assumptions on structural and effective connectivity  [4,6] to achieve a qualitative reproduction of experimentally observed neuronal activity. Here, we develop a quantitative validation approach where mean-field theory [2] guides the adaptation of a generic pointneuron network model to macaque motor cortex. We describe the characteristics of the experimental data extracted and used for comparison and present preliminary results for the generic network model. The underlying network model is an upscaled version of the Potjans &Diesmann [4] layered spiking network model extended to a size of 4x4mm2and a total of ~ 1.2 million leaky integrate-and-fire neurons [3]. In contrast to the original model this mesocircuit model uses lateral distance-dependent connection probabilities derived from cortical neuroanatomical data. To compare the output with observations we subsample single unit activities from the corresponding layer in the simulated network with the same number of neurons and with the same spatial arrangement of the recording array as in the experimental data. The model describes a system in ground, idle or resting state with uncorrelated input. In order to perform a quantitative comparison with experimental data we therefore conducted a resting state experiment with macaque monkeys not given any specific task or stimulus. We recorded neuronal activity from premotor and motor cortex using a chronically implanted 4x4mm2Utah array with 100 electrodes [1,5]. A video of the monkey was used to differentiate between periods of rest and spontaneous movements. The experimental single unit activities (~ 140 neurons) are subdivided into putative excitatory and inhibitory neurons based on their spike widths. We find that a) putative inhibitory and excitatory activity is in a balanced state, b) spike counts increase during movement, c) inhibitory units contribute more strongly to firing rate modulations than excitatory units, d) they also tend to be more strongly correlated among each other and e) the dimensionality of cortical activity is decreased during movement. Our results are to a large degree in accordance with mean-field theoretic predictions and may thus allow us to infer constraints on the parameter space of the mesocircuit model. Environment stimuli are continuously processed by the central nervous system (CNS) to better adjust, adapt, and learn new responses that optimize our benefits. At neural level, the external stimuli are coded as spikes of electric activity, called action potentials (APs). Neurons respond to changes in the environment by altering their firing speed, or phase, which means that instead of firing at a regular pace the neuron starts firing faster or they slow down. The amount of change, or resetting, in their firing period is determined by the timing, duration, and strength of the external stimulus. However, neurons connect with each other and create large networks capable of elaborated firing patterns that drive the response of the organism. We modeled the neural network as hierarchically-organized layers of neurons and in each layer the neurons' response is dictated by its own phase resetting behavior. We successfully generalized mathematically and then checked numerically that knowledge of how one isolated neuron responds to a stimulus can help predicting the response of a larger network to complex stimuli.

Generalized phase resetting and phase-locked mode prediction in biologically-relevant neural networks
Prefrontal cortex plays a crucial role in advanced cognitive functions.
Previous experimental observation has shown that glutamatergic inputs to the basal dendrites of cortical pyramidal neurons activate AMPA and NMDA receptors which can bring the dendrites into a long-lasting depolarized state: a dendritic plateau potential. These sustained depolarizations push the cell body towards spike threshold and reduce the membrane time constant. In such a "Prepared" state, the pyramidal cells can respond to other sparse synaptic inputs more quickly and easily, facilitating synchronization of firing. During the plateau depolarization, a neuron can tune into ongoing network activity and synchronize spiking with other neurons to provide a coordinated "Active" state (robust firing of somatic action potentials), which would permit "binding" of signals through coordination of neural activity across a population. Under this scenario, Active cells are recruited from cells in the Prepared state, and therefore the transient Active ensemble is embedded in the longer-lasting Prepared ensemble of neurons. We hypothesize that "embedded ensemble encoding" may be an important organizing principle in networks of neurons, explaining how electrical signaling endows central nervous system with capacity to form large number of neural ensembles. Also, embedded ensemble encoding pulls together two concepts (rate coding vs. temporal coding) that are typically seen to be in opposition. We have developed a Fig. 1 Modeling of the glutamate evoked dendritic plateau potentials morphologically-detailed model reconstructed from a cortical Layer 5 prefrontal pyramidal neuron in the NEURON simulator. Both synaptic AMPA/NMDA and extrasynaptic NMDA inputs are placed on basal dendrites to model the induction of plateau potentials (Fig. 1a-d).
The active properties of the cell are tuned to match the amplitude and duration of experimentally observed plateau potentials utilizing voltage-sensitive dyes in dendrites and whole-cell patch recording in soma (Fig. 1e). In addition, the effects of input location, receptor conductance, calcium-activated potassium channels and voltageactivated calcium channels are explored in detail. These findings help us to better understand the implications of dendritic plateaus at the cellular and network level. In the future, this detailed individual cell model can be used to develop cortical meso-scale network models for exploring the hypotheses pertaining to the recruitment of neurons into neural ensembles.

P236
Integrating large brain networks and network analysis to understand the epileptogenic zone Over 20 million people in the world suffer from medically refractory epilepsy (MRE). Approximately 50% of MRE patients have focal MRE, meaning that a small focal region in the brain, the epileptogenic zone (EZ), is the source of the seizures. For patients with focal MRE, treatment by surgical resection of the EZ can be effective, provided the EZ is reliably identified and entirely removed. Identification of the EZ often requires a surgical implantation of subdural grid or stereotactic depth electrodes (SEEG) electrodes, followed by a visual inspection of hundreds of EEG signals during seizure events that occur over several days to weeks. Clearly, surgical outcome relies heavily on precise localization of the EZ. We aim to integrate efforts from computational modeling and data analysis of SEEG recordings to better localize the EZ in an epilepsy patient. The Virtual Brain (TVB) is a computational platform that can integrate patient-specific information, such as brain connectivity derived from MRI and clinician's EZ hypotheses, to form personalized brain models capable of simulating realistic functional signals (i.e. SEEG). From a data analysis perspective, we used a novel network-based algorithm, coined the "fragility algorithm", that has demonstrated capabilities of localizing the EZ by analyzing network stability with respect to nodal perturbations. The fragility algorithm, unlike single channel frequency analysis, looks at the SEEG as a network and can efficiently analyze a SEEG network to create a heatmap of predictions on the EZ. The fragility algorithm determines which nodes within the epileptic network (i.e. SEEG channel) are the most fragile, i.e., nodes that if connections from it are perturbed slightly will destabilize the network. Fragility weights for each node are then used to predict the EZ. We built personalized brain models for two temporal focal MRE patients to determine the algorithm's predictions in two different situations. We simulated for each patient in silico: (1) inside: placement of the EZ at the clinical hypothesis and (2) outside: placement of the EZ outside the resection region. We then applied the fragility algorithm on the simulated and actual SEEG data from each patient to see if the fragility maps of the simulated EZ scenarios resemble the fragility maps derived from actual recordings. With TVB integrated with fragility analysis, we can hypothesize where the true EZ might be for a given MRE patient and whether or not it was correctly identified by clinicians using standard visualization methods. In one patient who had a successful surgery, we assume the EZ lies within the resected region and we found that the predicted EZ in the real SEEG and simulated data as identified by the algorithm, matched the clinically annotated EZ. In contrast, the patient with failed surgery, we assume the EZ lies outside the resected region and we found that the predicted EZ in the real SEEG data and the simulated data does not match the clinical EZ. These results suggest that the failed epilepsy surgery was due to the fact that the EZ was not within the resected region, while in the success case it was. These results outline how personalized brain models can help determine sensitivity of EZ localization algorithms to locations of the EZ and it can integrate with data analysis to validate whether the EZ is properly localized in a surgical resection. Knowledge of mesoscopic brain connectivity is important in understanding inter-region communication and information processing. Models of structural connectivity have been used to investigate the relationship with functional connectivity, to compare brain structures across species, and more [5,2,4,6]. Early models of connectivity were constructed with the assumption that the regions are homogeneous [3]. This assumption is useful, but it depends on predefined regional parcellations and describe connectivity at a region-limited level of resolution. Here, we go beyond the regional approach and construct a model of the whole brain connectivity at the scale of 100 micron voxels, extending the previous work of [1]. We use the Allen Mouse Brain Connectivity Atlas, a large scale dataset measuring with two photon tomography the projections of sets of neurons infected with a viral tracer to approximately 5x10^5 target voxels [3]. While this dataset is large, the 429 sources we have in wild type C57BL/6 mice pale in comparison to the 2.5 × 10^5 source voxels (all experiments were conducted in one hemisphere). To meet this challenge, we propose a model which relaxes the assumption of homogeneity of connections within a region and instead assumes smoothness across major brain divisions. It is a simplification of the structured regression model proposed in [1], which outperformed a regional model in the visual cortex. We model the connectivity at each source voxel as the kernel-weighted average of the projection patterns of nearby injections (Fig. 1). The kernel is a radial basis function with bandwidth hyperparameter fit in the inner loop of a nested K-fold cross validation procedure. The voxel-scale model strongly outpredicts the previous regional modeling both in relative mean squared error in cross-validation, and when compared to a human-curated dataset. We plan to release the code to perform construct this network. As our model performs an interpolation, it inherently produces a dense connectivity. While the true sparsity of the mouse connectome is debated [7], a manual analysis of true positive and true negative signals in the data [3], poses the question of the relevance of very weak connections. Thus we plan to use this model to examine the whole brain distribution of connection weights, and the importance of very weak connections on common graph theoretical models.
The new voxel-scale model has been used in several applications including analysis of the modularity of the mouse cortical network. This new voxel-scale model of the mouse connectome permits researchers to extend their previous analyses of structural connectivity to unprecedented levels of resolution, and allows for natural comparison with functional imaging and other datasets. Fig. 1 An illustration of the voxel-scale model projected into a cortical surface view. We predict the connectivity at each voxel in the brain using a kernel weighted average of the injections in each major brain division Convolutional neuronal networks have had great success in many applications [1] as well as in describing neural responses in the ventral stream of primates [2]. However such networks require labeled data to train and are quite brittle: a single pixel change can with high confidence change the category prediction of the classifier [3]. In contrast to what is known from biology, these networks focus on feedforward visual processing, largely ignoring the influence of lateral connections. They also focus on supervised rather than unsupervised learning. Can we use knowledge from biology to improve the robustness of these systems by including recurrent lateral connections learned in an unsupervised manner? Visual processing in the brain makes use of recurrent lateral connections. Surround suppression, contour integration, and figure-ground segmentation are all examples of forms of contextual modulation that cannot be simply explained by feedforward mechanisms, and are often included as part of the extra-classical receptive field. Experimental and modeling studies suggest different excitatory and inhibitory cell types are important for mediating these lateral connections. The brain is also able to learn structured representations in a largely unsupervised manner with little to no labeled data.
Here, we explore how representations learned in an unsupervised way can be used in conjunction with features learned in a supervised manner to improve the robustness of convolutional neuronal networks. We leverage recent work [4] demonstrating how optimal lateral connections can be learned through a modified Hebbian learning rule. We combine optimal context integration with deep learning to not only learn different features within a convolutional neuronal network, but also the optimal lateral connections between the neurons encoding these features. These connections can be implemented using a collection of cell types and lateral connections which match well biological observations in cortical circuits. We demonstrate the influence of lateral connections by testing our model on two standard image datasets: MNIST and CIFAR-10. In additional to the original images, we also generate noisy versions of the images by adding additive white gaussian noise and salt-and-pepper noise in increasing levels. These noisy images are used to test the robustness of the trained convolutional neural networks. In contrast to traditional deep learning approaches, we combine supervised and unsupervised training in our model. We then test models with and without optimal lateral connections on the original images in each dataset, as well as on the noise-corrupted images. As expected, increasing levels of noise degrades model performance for both the MNIST and CIFAR-10 datasets; however, models with optimal lateral connections are more robust to this noise and can achieve higher classification accuracies on both datasets, especially at higher noise levels. This result demonstrates the potential usefulness of combining supervised and unsupervised learning techniques in real-world vision tasks. Although the natural scene statistics of MNIST and CIFAR-10 differ greatly, the same unsupervised learning rule worked for both systems. Our results suggest that the integration of lateral connections into convolutional neural networks is an exciting avenue of future research (Fig. 1). Neuronal responses in early visual cortex are primarily driven by inputs to the classical receptive field and are influenced by stimuli in the receptive field surround. This type of spatial contextual effect is thought to arise due to the statistical structure present in natural scenes, with the surround providing context for the information in the classical receptive field. Lateral connections between neurons in the same cortical area are generally thought to be responsible for transmitting information in the near surround. In a number of recent experimental studies, excitatory neurons have been demonstrated to have like-to-like connectivity with neurons coding for the same feature (e.g. orientation) preferentially connecting to each other with higher probability and/or strength and specific rules for connectivity of inhibitory neuron types have been described. On the other hand, normative models of lateral interactions, relying on sparsity and saliency in the optimal representation of natural Starting from the assumption that each excitatory neuron represents the probability of a feature being present in the sensory stimulus, we hypothesize that lateral connections serve to optimally (in a Bayesian sense) integrate evidence from the surround. We show that such optimal integration of contextual information can be implemented by a neuronal network. Using natural scene statistics obtained from the Berkeley Segmentation DataSet and in vivo electrophysiological data from awake mouse V1 neurons, we compute the synaptic weights between neurons in our network resulting from the optimal integration of contextual information. We show that this network has like-to-like connectivity between excitatory neurons, in agreement with experimental observations. However, in our model, even neurons with non-overlapping classical receptive fields can have strong connections if they code for features which often co-occur in natural scenes. The distance dependence of connections is similar to those observed experimentally and is found to be heavier-tailed than an exponential decay. These results generalize to other classes of receptive fields including gabors and sharp on-off band-like receptive fields. The network also needs multiple types of inhibition -local normalization, surround inhibition and gating of inhibition from the surround -which we map to the parvalbumin (PV), somatostatin (SST and vasoactive intestinal peptide (VIP) expressing interneuron cell classes respectively. We compared our local circuit model with a correlation-based model, in which the lateral connectivity between cortical neurons is determined by the correlation between their classical receptive fields. We find that the correlation-based model results in a much steeper decay of both like-to-like connectivity and distance dependence, compared to our model. We also show that, compared to a feedforward network, the presence of this local network structure increases the capacity to reconstruct a natural scene from the activities of neurons in the network under noisy conditions. We hypothesize that optimal integration of context is a general computation of cortical circuits, and the local network rules constructed for mouse V1 generalize to other areas and species. Retinal ganglion cells, which transmit the output of the retina to the brain, typically have a concentric center-surround receptive field and consist of ON and OFF types. Although the retinal responses to natural [3] and artificial [1] scenes have been well characterized, it remains unclear how efficient the retinal code is at transmitting visual information. Here we seek to identify the biological constraints of the visual system and the redundancies in natural scene statistics that have shaped the retinal code, by varying respectively the architectural constraints of and statistics of inputs to artificial neural networks trained to classify objects in natural images. We find that when we allow an overcomplete representation in the early layers (i.e. many more neurons in each layer than pixels in the image), the trained network exhibits oriented receptive fields in these layers. However, when we severely restrict the number of neurons in these early layers (guided by the intuition that the optic nerve can only contain a limited number of fibers), concentric center-surround receptive fields emerge. We also find that the response patterns of these neurons naturally cluster in two functional groups analogous to ON and OFF cells in the retina. Moreover, the receptive fields of the first overcomplete layers downstream of our artificial retina are oriented, like receptive fields of V1. Examining the connection weights between the two layers, we find that these oriented receptive fields are generated by drawing from a few center-surround neurons along an axis, mirroring Hubel and Wiesel's hypothesis about simple cells in visual cortex. These results suggest that the retinal code, with its two main cell populations and its concentric receptive fields, could be optimized to transmit relevant visual information to the brain with a limited number of neurons. This interpretation of the constraints underlying the retinal code is an alternative to that of Karklin and Simoncelli [2], which suggests that the retinal code is optimized for metabolic and noise constraints. Finally we provide evidence that the utility of center-surround encoding for compression arises from covariances among orientation-selective neuron activations in response to natural scenes. Bandlimited power in cerebral cortex occurs in the form of transient events, a fact often obscured by spectral averaging in classical empirical studies of brain rhythms. Events in the beta band (15-29 Hz) may be locally rhythmic but are brief and sparse, while occurring in scalp-level and mesoscale electrophysiology across diverse neocortical regions in humans, nonhuman primates and rodents. Beta is implicated in nearly every aspect of cortical function, varying during attention, perception and decision-making as well as motor control. However, we presently lack amechanisticmodel that might afford beta events a functional role susceptible to pharmacological, electromagnetic or genetic controls. We are motivated by studies in primary somatosensory cortex (SI), where beta is consistently reported to have suppressive effects on the perception of subsequent tactile stimulation-even a full second after an event occurs. To investigate the underlying mechanisms, we use a geometrically reduced biophysical circuit model consisting of motifs universal to neocortex, which fits both perceptual and mean magnetoencephalographic (MEG) tactile evoked responses (ERs) after event occurrence in SI, while offering unprecedented predictivity of ER features. Previous modeling indicates that beta events are generated by calcium-mediated bursts incident on L1 with one primary source in nonlemniscal thalamus. Recent studies show that such inputs act simultaneously on superficial pyramidal and interneuron subnetworks. In our model, incoming bursts prime pyramidal dendrites while recruiting neurogliaform (NGF) cells either directly or through electrical coupling after interneuron synchronization. Modeling shows how both mechanisms act in sequence, with NGF suppression occurring more often due to the slow timescale of GABABinhibition (> 250 ms), explaining the net behavioral effect. Before NGF is recruited, however, incoming burstsamplifystimuli through subthreshold facilitation -similar to previous oscillation models but with significant improvements in precision. Our model predicts the correct beta phase at stimulus onset in "hit" trials, and also makes two precise predictions of poststimulus bandlimited phase coherences verified at 160-200 Hz at 40 ms (from model L2/3) and 90-110 Hz at 60 ms (from model L5). After this priming phase, L2/3 NGF cells act through GABAB channels on pyramidal somata (in L2/3) and on L5 middle apical dendrites. This mutes L5 pyramidal bursts and clamps L2/3 pyramidal somatic voltages to the potassium reversal. The latter effect is essential for a conspicuous ER feature generated by a completely novel mechanism, BMC Neurosci 2018, 19(Suppl 2):65

P241
where high-voltage downward-propagating dendritic spikes collide with a soma fixed at its minimum voltage. This causes a biophysically maximal local current ideal for driving signaling cations into the soma, raising the possibility that GABAB suppression gates learning after early-prediction error. Our modeling circumscribes nearly the entire known phenomenology on beta, generating predictions verified at the scalp level with invasively testable analogues. In the process, we unite previously unlinked findings among theoretical, electrophysiological, and anatomical domains. Our results shed particular light on neurogliaform GABABaction in neural computation, with concrete manipulations accessible through a suggestive pharmacopeia including alcohol, opiates, serotonin and neuropeptide Y. In our companion poster, Multiplexed coding using differentially synchronized spikes: Part 1, theory and simulations, we demonstrated in silico that the rate of asynchronous spikes can encode the intensity of a slow signal while the timing of synchronous spikes simultaneously encodes abrupt changes in stimulus intensity. This suggests that a single set of neurons can represent distinct features of an external stimulus using differentially synchronized spikes. Here we tested in vivo whether this occurs in real neurons. To this end, we recorded extracellularly from single units in the primary somatosensory cortex of sedated rats. Computer-controlled mechanical stimuli were applied to the whisker pad as discrete steps of increasing force. We asked whether we could decode (1) the force during the sustained phase of the step based on the rate of asynchronous spiking and (2) the timing of the onset and offset of the step based on the timing of synchronous spikes. Using responses from 17 neurons deemed to be responsive to our mechanosensory input, we constructed a firing rate histogram (FRH) in two ways (Fig. 1). By constructing the FRH with a broad Gaussian kernel (500 ms), we found that the magnitude of the firing rate tracks the intensity of the force. By constructing the FRH with a narrow Gaussian kernel (5 ms), we found that abrupt changes in the force were reflected in abrupt increases in the firing rate caused by transient synchronization. Applying a synchrony threshold to the latter FRH yielded 86% sensitive and 100% specific detection of the stimulus transients. These findings definitively demonstrate that synchronydivision multiplexing occurs in somatosensory cortex. Specifically, we have shown that stimulus intensity (a first-order stimulus feature) is encoded by the rate of asynchronous spiking while abrupt variations in stimulus intensity (a second-order stimulus feature) are simultaneously encoded by the timing of synchronous spikes. The results help reconcile apparently contradictory evidence for rate or temporal coding by showing how these coding strategies can operate together.

P243
The Neurons use action potentials, or spikes, to encode information.
Proper neural coding thus relies on the proper control of spike initiation. The spike initiation process reflects the highly nonlinear interaction between different ion channels, which means that subtle variations in ion channel expression or function can dramatically impact excitability. Yet excitability is normally very stable, which raises the question of how excitability is regulated so robustly. Emerging data argue that the biophysical basis for excitability is highly degenerate, meaning that many different combinations of ion channel conductances can yield equivalent excitability. This degeneracy is thought to facilitate the robust regulation of excitability by allowing changes in any one ion channel to be compensated for by changes in many other ion channels. We hypothesized that parameters that are able to compensate for one another will be correlated. We further hypothesized that the strength of pairwise correlations will be weakened as the degree of degeneracy (i.e. the number of parameters that can compensate for one another) increases. To test these hypotheses, large sets of conductance-based Morris-Lecar models were generated with randomly chosen parameter values describing different conductance densities or activation properties. A different number of parameters was allowed to vary for each set of models; all other parameters were held constant at their baseline values. From these sets, we identified model neurons with comparable excitability and, using only those models, determined the pairwise correlation between parameters that had been randomly determined. Correlations were observed between some parameters and, as further predicted, the strength of correlation decreased as the number of randomly varied parameters increased. Based on these results, we expect that highly degenerate systems will exhibit only weak pairwise correlations in their parameter values. Conversely, strong correlations may suggest that a system has only modest degeneracy and that it may, therefore, be less able to compensate in the face of strong perturbations. Multiplexing refers to the simultaneous transmission of multiple signals through a single communication channel. In engineered systems, multiplexing is often implemented by partitioning different signals to different frequency bands (frequency-division multiplexing) or to different temporal epochs (time-division multiplexing). Mounting evidence suggests that the brain also multiplexes but it remains unclear how this might occur. We hypothesized that the brain can form multiplexed representations of first-and second-order stimulus features (i.e. stimulus intensity and abrupt variations therein, such as occur at edges) using spikes that are differentially synchronized across a set of neurons receiving common input. To test our hypothesis, we built a feed-forward neural network comprising Morris-Lecar (ML) model neurons. All neurons received a common mixed input constructed from two distinct signals, slow and fast, plus uncorrelated fast noise. The slow and fast signals represent input from upstream sensory neurons tuned to first-or second-order stimulus features based on their low-or high-pass filter properties, respectively. The two sensory streams converge on the ML model neurons. According to our hypothesis, slow and fast signals are independently encoded by different types of spikes. Specifically, the rate of asynchronous (Async) spikes encode the slow signal whereas the timing of synchronous (Sync) spikes encode the fast signal. To assess the feasibility of the multiplexed coding scheme, we fit linear-nonlinear (LN) rate models to PSTHs from our conductance-based spiking models. In a conventional LN model, input passes through a linear filter and then through a static nonlinearity whose output is firing rate. We constructed a multiplexing LN model with two parallel streams; the same mixed signal is presented to both filters but the output of each filter passes through a different nonlinearity. Unlike the two input streams, which represent input from two differently specialized sets of sensory neurons, the two streams within the LN model represent two operating modesintegration (low-pass filtering) and coincidence detection (high-pass filtering)-used by a set of neurons operating in a hybrid mode. The two-stream LN model more accurately predicted true firing rate (using the PSTH of spiking models as reference) than the one-stream LN model, especially for synchronous spikes. Our results demonstrate that a set of cortical pyramidal cells can implement multiplexing by simultaneously encoding slow and fast features of a mixed signal through a multi-modal filter. These results are further validated experimentally, as presented in our companion poster Multiplexed coding using differentially synchronized spikes: Part 2, experiments. Mathematical modelling significantly contributes to our understanding of the mechanisms underlying neuron and network behaviour. Hodgkin-Huxley (HH) equations are frequently used to model neuron conductances. The majority of the existing HH conductance modelling workflows are based on data acquired under current clamp or voltage clamp (VC) recording conditions. While current clamp provides phenotypically rich information on a neuron's firing properties, it is often difficult to disentangle the influence of the multiple conductances at play in a real neuron. Isolation of a single conductance by pharmacological or heterologous expression simplifies analysis; however, exhaustive exploration of the kinetic and voltage properties using VC is needed for HH parameterisation that requires weeks of recording time. A more experimental time efficient approach would not only benefit current efforts to develop specific HH models for the entire set of voltage gated channels expressed in the brain, but also free resources to explore these ion channels under pathological or pharmacological conditions. We present an improved HH conductance modelling workflow that uses data derived from dynamic action potential clamp (DAPC) recordings to extract conductance model parameters. A major difference in DAPC over VC recording is that DAPC is an action potential (AP) weighted systems identification approach and is more aligned to the common final use of HH models, which is the building of single neuron and network AP firing models. In this study, first, using fully simulated conditions, we show that with as little as one second of DC recording time, we can produce parameters with an average error of less than 4%. When deployed into simple neuron models, these parameters produced firing rates that approached 100% accuracy in fully simulated experiments. Second, we undertake a real-world test using NaV1.2 channels and show that training our model with five or fewer APs could produce a HH conductance model that predicted the subsequent AP firing with 97% firing rate accuracy. Further, the AP traces overlapped with 94% accuracy. We conclude that DAPC based workflows can be as, or even more, accurate than VC based workflows for extracting HH conductance model parameters. Importantly, this accuracy can be obtained with considerably less recording time and effort providing additional opportunity for exploration of HH conductance equations in different experimental conditions and positioning this approach to have the efficiency to exploit current advances in ion channel genetics and precision drug discovery. The ability to rapidly learn new objects and environments is a critical task for the brain and the locations of landmarks is crucial to learning new environments. Grid cells encode local position information through regular, tiled firing fields that are anchored to each learned environment [1]. Grid cells are organized in modules of cells that share the same scale and orientation. Individual modules are ambiguous over larger spaces but sets of grid cell modules with varied scale and orientation accurately and unambiguously encode locations over many large environments [3]. Experimental results show that landmarks and sensory cues are critical for anchoring grid cells and learning new environments, but the exact mechanisms are unclear. It has been proposed that each grid cell module operates like a residue number system [2] and can encode relationships between layers [4]. In this paper we exploit a core property of residue number systems that allows a small, fixed set of connections to perform arithmetic operations. We propose a biologically-plausible model with two grid cell layers that encodes the relative locations of landmarks in environments and show that it reliably encodes environments in the presence of noise. Grid cells in each layer are organized into modules, where each module represents different relative scales and orientations. The "location layer" encodes the location of each sensed landmark. Grid cell modules in the "relative location layer" are paired with a corresponding module in the location layer. Relative location operations are performed only between these paired modules without access to other modules in the layers. The difference (subtraction) between the phase of currently active neurons and the neurons corresponding to landmarks seen in a recent period is encoded in the relative location layer. Locations for recently active landmarks are updated at each step in the location layer with a motor command. However, the relative location layer always represents the relative positions of each landmark pair and is invariant to absolute position within the environment. The model thus learns a set of representations of relative positions of landmarks that is stable for a given environment, and enables efficient disambiguation of previously seen environments (Fig. 1). Our model is compared with two other models: (1) a bag-of-features model that only compares landmarks without locations and, (2) an ideal model that exhaustively examines all environments to find the BMC Neurosci 2018, 19(Suppl 2):65 best match. Using the relative positions of landmarks, our model is able to achieve perfect accuracy when there is little noise and lags the ideal model slightly for very noisy test cases. The bag-of-features model is no better than chance when a small pool of five landmarks is used. Further research will explore the generalization ability of the model and the addition of an unsupervised temporal clustering layer that can reinstate learned relative location representations in order to predict input sensations that have not recently been seen.

Unsupervised learning of relative landmark locations using grid cells
Recognizing and tracking the direction of visual objects from the environment is crucial to the control of many animal behaviors. During development, the synaptic circuitry of V1 requires anatomical and functional refinement essential for the emergence of direction selectivity (DS). Under such refinement, synapses strengthen or weaken depending on the relative timing of spikes, a phenomenon known as spike timing-dependent plasticity (STDP). In addition to STDP, a ubiquitous feature in neural circuits is short-term plasticity (STP), where the strength of synapses varies from milliseconds to seconds as a result of recent activity. Early in development, Layer2/3 (L2/3) neurons recorded in V1 of ferrets exhibit spatial clustering, where neighboring neurons exhibit a weak preference for direction [1]. When the network achieves functional maturity, weak direction biases (WDBs) acquired early in development become accurate predictors of strong DS. Despite the crucial role of visual experience needed for proper establishment of DS maps, circuit mechanisms underlying the development of DS remain elusive. In studies of developing neural circuits, STDP has featured prominently in explaining DS through the formation of asymmetric intracortical connections, with weak connections outnumbering strong ones [2]. Amongst infrequent but strong connections, neurons exhibit correlated responses to visual stimuli moving in a particular direction of motion. Given the implication of STDP in explaining DS, what potential mechanism may set the initial condition in driving the onset of WDBs? Using a bio-inspired LIF model, we examine the role of STP as a complementary mechanism to STDP in explaining DS. We simulate a network of excitatory connections with initial weights reflecting the distribution of early L2/3 connectivity. The model captures responses of DS neurons during circuit refinement observed in V1 data (Fig. 1a). While the network receives inputs, the STP mechanism introduces direction-dependent modifications in synaptic strength, which in turn breaks the symmetry of the excitatory connections, introducing an asymmetric connectivity scheme that drives the onset of WDBs. In turn, STDP refines these connections by progressively increasing the asymmetry, making weak connections more prominent, and strong connections a rare commodity, as observed in L2/3 neurons (Fig. 1b). This connectivity scheme reflects preference for direction that matches the WDB introduced by STP. In combining both STP-STDP, and by treating both mechanisms in isolation, we compare versions of the model to analyses of V1 data [3], where increases in DS during subsequent days are primarily due to decreases in spiking responses to the null direction and orthogonal orientation. Learning-related changes in responses are quantified by direction circular variance (Fig. 1c). Overall, results highlight the role of STP as a complementary mechanism to STDP in explaining DS, and propose that STDP is necessary but insufficient to account for the predictive basis of DS.  Sense of agency, i.e., the feeling that oneself caused something to happen, is fundamental to the experience of volition, self-consciousness and responsibility for one's own actions, and the degradation of this experience characterizes certain psychiatric disorders. Despite its irrefutable significance, the literature still lacks a mathematical exposition of the computational principles that underlie it. We theorize sense of agency as theconfidencein one's perception of action-outcome effect to be consistent with the hypothesis that the self is thecommon sourcebehind this effect. We adapted the Bayesian inference model of Sato, Toyoizumi and Aihara [1] that was originally used to explain the ventriloquism effect as a Bayesian estimate of a common source behind the audio-visual stimuli. Formalizing sense of agency by this Bayesian principle distinguishes our theory from existing works. Intentional binding, i.e., the perceived shortening of the time interval between voluntary action and its outcome, has been reported as an implicit measure of sense of agency. Yet, the exact nature of this link is far from comprehension [2]. Our Bayesian model gives a simple coherent account of this link: the shorter perceived interval between the actionoutcome timings is more consistent with the causal role of one's action in producing the immediate outcome, and thus increases the confidence of the Bayesian estimate, modeled as sense of agency. We compared the predictions of our model to the results of two pertinent intentional binding studies. The first follows the seminal experiment reported by Haggard, Clark & Kalogeras [3] that showed voluntary actions produced intentional binding effects but involuntary actions produced the prolonged opposite perception of the action-outcome interval. The second case follows the study of Wolpe, Haggard, Siebner & Rowe [4] that investigated the contribution of sensory uncertainty to intentional binding by manipulating the intensity of outcome tones. They showed that when the outcome reliability was reduced, action binding was diminished and tone binding was increased. Our Bayesian psychophysics model reproduces these empirical results based on a computational principle. This study investigates dynamical network models designed to reproduce the electrophysiological responses of the human auditory cortex to mismatch negativity eliciting stimulations. To this end, we first focused on modeling On and Off responses of tonal stimuli via recurrent circuits as they exist in the auditory cortex [1][2][3]. In the simulation, the recurrent circuits are represented by a two-layer network, where the input stimuli from the thalamus reach the 1stlayer and indirectly affect the 2nd-layer activities through recurrent inter-layer connections. With a stream of stimuli fed to the 1stlayer (input), various types of On/Off responses can be reproduced in the 2ndlayer (observation) given proper inter-layer connections. The simulation results account for relevant properties of cortical On/Off responses and provide thereby clues about the underlying physiological mechanisms. (1) A subtle change in interlayer connections switches the response type between On, On and Off, and Off. Furthermore, it can also switch between 'sustained' and 'suppressed' activity during the stimulus presentation. Hence, the diverse On/Off responses observed at different locations in auditory cortex [2] may reflect diverse inter-layer connections between the input and the observation layer. Interestingly, symmetric interlayer connections do not give rise to On/Off responses, underlining the importance of asymmetric forward-backward interactions for the change detection function at the cortical level. (2) The distinct onset-and offset-frequency-receptive-fields (FRF) observed in A1 neurons in [3] can be accounted for by the two-layer scheme. We conclude that the tonotopically organized input layer has distinct equivalent inter-layer connections with the observation layer. (3) Furthermore, the simulation demonstrates that generation of Off response in the 2ndlayer relies on tonic inhibition in the 1stlayer during the stimulus. This nicely matches physiology, as the reduced Off response by NMDA receptor antagonists observed in [1] is due to the reduced inhibition during the stimulus, because excitatory synapses on inhibitory neurons are more sensitive to the NMDA receptor antagonists [4]. To summarize, the recurrent circuits in our model provide a parsimonious solution for the change detection function at the cortical level. This model for On/Off responses will be further used to reproduce mismatch responses to omitted and deviant stimulus in the oddball paradigms. Various techniques to scale down models of biological neural networks have been developed. However scaling can never preserve all properties of a network, especially not the high number of connections between neurons observed in the brain. Therefore, simulating large-scale biological neural network models remains important and doing so in a reasonable time is one of the major technical challenges in computational neuroscience. Conventionally large-scale simulations are executed on High Performance Computing (HPC) clusters and the tools to distribute neural network simulations across such systems are now relatively mature [1]. However, HPC systems are expensive and not well suited to real-time simulation. Bespoke 'neuromorphic' hardware has been developed to address these problems, but they come with their own challenges and limitations. Graphics Processing Units (GPUs) were initially designed to accelerate the rendering of 3D graphics, but have evolved into massively parallel accelerators, programmable at a relatively high-level of abstraction. This has led to their use in accelerating a wide range of HPC applications and proved vital in the practical application of deep learning. The computational properties of biological neural networks simulations are not necessarily a natural fit for GPU acceleration but modern GPUs can still achieve significant speed-up across a wide range of realistic models. In this work we present the latest developments in GeNN [2], our open source library for generating optimised neural network simulation code for NVIDIA GPUs. GeNN does not require that users understand the details of GPU programming; instead it lets users describe neuron and synapse models using C-like code and then generates the code for running the resultant models on GPUs. To demonstrate the power of this approach we re-implemented a point neuron model of a cortical column [3], which was previously simulated at 0.2× realtime using a 23-node HPC cluster [4]. Using GeNN and a single GPU, we are able to simulate the same model at only 0.25× real-time. This model also demonstrates how GeNN allows the initialisation of a network to be offloaded onto the GPU, significantly reducing overall runtime. GeNN can be used directly from C++, allowing users full control of their simulation and making it easy to integrate GeNN with other libraries. This is particularly useful, now that GPUs such as the Jetson TX1 are available in form factors small enough to use on board autonomous robots. We demonstrate this using a recent model of the bee central complex [5] which we re-implemented in GeNN and used as a controller for a wheeled robot equipped with a Jetson TX1. Finally, for those requiring simulator-agnostic model creation, we have created interfaces to use GeNN from SpineML and the Brian 2 simulator-making it trivial for existing users of these tools to take advantage of GPU acceleration. We demonstrate this with a SpineML model of optical flow calculation inspired by the honeybee [6] (Fig. 1). In many applications, it is important to understand how quickly and reliable a neuron receiving noisy background input would respond to an external input signal. Here, we consider a leaky integrate-and-fire neuron receiving Gaussian, exponentially correlated noise, and obtain analytically the probability that the neuron fires at least one action potential when an external input signal, in the form of post-synaptic potentials (PSPs), is applied. We derive expressions for the probability of firing in three regimes, when the dynamical time scale of the external input is much quicker than the membrane time constant, when it is much slower, and when it is of comparable speed. In case of fast inputs, the firing probability can be obtained by assuming that the probability distribution of the free membrane potential is translated according to the size of the signal, In the case of slow inputs, we can use the adiabatic approximation, assuming that the membrane potential distribution is quasi-stationary, and obtain the firing probability using previously developed methodologies [1,2]. In the most challenging case of comparable timescales of membrane potential and input fluctuation, we derive a novel correction to the adiabatic case by applying a Taylor expansion on the noise-induced component of the membrane potential. To determine the spiking probability for arbitrary input waveforms, we heuristically combine the three approximations depending on the time derivative of the input. Our analytic approximations match well with the results obtained from numerical simulations for input signals of a wide range of sizes and shapes, and under both strong and weak background noise (See Fig. 1). Our method can be used to estimate the response latency of sensory neurons and to study the reliability of neural responses to temporally filtered and jittered synchronous input.    [2,3]. Here, we present a MB model that captures these data to compute reward prediction errors (RPEs) for learning, thus implementing the Rescorla-Wagner model. Current models posit that the alpha− (A) and beta− (B) lobes of the MB encode the signed valence of reward information and actions [1]: DANs in the A-lobe (D− in Fig. 1a) are excited by negative (−ve) rewards, and depress active KC synapses onto MBONs that bias actions toward approach (M+ in Fig. 1a); DANs in the B-lobe (D+ in Fig. 1a) are excited by positive (+ve) rewards, depressing active KC synapses onto MBONs that bias actions toward retreat (M− in Fig. 1a). If MBONs provide excitatory feedback to their respective DANs, the learned reduction in MBON firing can offset the excitatory reward signal arriving at that DAN. Thus, D+ and D− may both encode RPEs in the signed (+ve or -ve) reward valence. Moreover, the difference in MBON firing rates, m+-m−, signals the learned net valence associated with a sensory cue. We first show two problems with this model: (1) it cannot learn reward magnitudes above an upper bound; (2) it learns only when KC-DAN excitation is minimal or absent, in contrast with experiments [2]. We propose a solution, in which D+/D− neurons are instead inhibited by -ve/+ ve reward signals, and in which KC-DAN excitation is required (Fig. 1a). We also derive a plasticity rule for KC-MBON synapses that performs gradient descent on the RPE, and that resembles experimentally observed rules [4]. We call this model the Signed Valence Circuit (SVC). As before, DANs encode RPEs in the signed reward valence (Fig. 1d), and the difference in DAN firing rates, d+-d−, yields the net RPE (Fig. 1e). The SVC can learn rapid changes to reward contingencies in just 5-10 trials (Fig. 1c). In the SVC, D+/D− respectively signal RPEs for -ve/+ ve rewards, so do not actually contribute to learning +ve/−ve valences, counter to experimental evidence [1]. However, in a dual version of this circuit-in which D+/D− are driven by +ve/ve rewards-D+ no longer signals decrements in -ve rewards, again in contrast with experiments [5]. We therefore combine the SVC and its dual to produce the Signed RPE Circuit (SRC; Fig. 1b), in which the lobes encode the signed RPE of both +ve and -ve reward signals (Fig. 1G). Lastly, the SRC performs well in a traplining task (Fig. 1H-I)-repeating learned routes and minimizing the distance traveled between feeding areas-a behavior exhibited by bees [6] and other species, and a foraging analogue of the travelling salesman problem.

Firing probability for a noisy leaky integrate-and-fire neuron receiving an arbitrary external input signal
Synapses do not merely relay spikes from presynaptic to postsynaptic neurons. Instead, they edit outputs of presynaptic neurons depending on presynaptic spike history. During stimulation period, synaptic transmission can facilitate or depress [1][2][3]. However, the exact functions of this 'short-term synaptic plasticity' remains contentious partly due to its diversity. Thus, it is necessary to characterize short-term synaptic plasticity to comprehend its functions.
To this end, we analyzed the data from Allen Institute's large scale multi-patch pipeline project, which focuses on studying synaptic connections in mouse primary visual cortex (V1). As short-term synaptic plasticity depends on both pre-and post-synaptic neuron classes [4][5][6], we defined the synapse classes using pre-and postsynaptic neuron classes and studied their properties. To characterize short-term synaptic plasticity in synapse classes, we developed descriptive models. Specifically, our models considered five temporal dynamics of synaptic transmission: depression (Eq. 1 in Fig. 1a), facilitation (Eq. 2 in Fig. 1a), use-dependent modulation (Eq. 3 in Fig. 1a), desensitization of postsynaptic receptors (Eq. 4 in Fig. 1a) and slow modulation of release probability (Eq. 5 in Fig. 1a). The amplitude of postsynaptic potentials (PSPs) is proportional to resource (n), utility (P) and sensitivity of receptor (S); see Eq. 6 in Fig. 1a. Depression and facilitation were modeled by the dynamical systems proposed by the two groups, and the rest were modeled by those presented in an earlier work [7].
In the experiments, PSPs were measured while presynaptic neurons were stimulated with spike trains at 10, 20, 50, 100 and 200 Hz, and we fitted the time courses of PSPs to our synapse models. Specifically, we first used x-means clustering to identify homogenous synapses in the class. Then, we fitted the average PSPs from homogeneous synapses to the model using 'LMFIT' , an open-source package developed for flexible non-linear least-square minimization [8]. So far, we constructed ~ 10 synapse models in V1 that can capture short-term synaptic plasticity observed in V1 (See Fig. 1 for an example). These models suggest 1) that most synapse classes depress in V1 and 2) that short-term synaptic plasticity depend mainly on presynaptic neurons. We believe that these synapse models would allow us to better understand the neural basis of visual perception. The following points should be underscored. First, our synapse models will be further refined, as more data become available. Second, we are currently building network models of V1 with our synapse models to study functions of short-term synaptic plasticity in V1 in visual perception. For instance, we are studying short-term synaptic plasticity's contribution to the stimulus-specific adaptation in V1. The inferior olive is thought to contribute to the generation of timing and error signals for motor control. However, the specific role of its distinctive spatiotemporal activity patterns is still controversial. Synchronous states have been observed both in olivary subthreshold oscillations in vitro and in complex spikes of cerebellar Purkinje cells in vivo, which serve essentially as a one-to-one readout of spiking in olivary climbing fiber output. These observations show two types of spatial structure: zero-lag correlations that decay with distance (localized synchrony), and synchrony at nonzero phase lags that increase with distance (phase waves). Although evidence suggests that electrical synapses between olivary cells underlie this synchronization, the details of spatial pattern formation remain largely unexplained. We present a model that explains how both types of synchronized states arise from the interplay of subthreshold membrane potential oscillations and spatially constrained electrical coupling. Although the olive is known to have the highest density of electrical coupling (via dendro-dendritic gap junctions) of any brain region, olivary cells also exhibit a high level of heterogeneity in the intrinsic oscillation frequency, which opposes the synchronizing effect of the coupling. We capture these features in a network model of phase oscillators, derived from conductance-based models of olivary neurons following the theory of weakly coupled oscillators (Fig. 1). We show that a network of phase oscillators with the appropriate balance of heterogeneity and local coupling will produce spatial phase gradients which are stable and robust to perturbations of the network. This may be a generic mechanism for brain to harness random heterogeneity to create stable patterns of timing. We next consider the patterns of spiking triggered by external input to a network in the subthreshold phase gradient state. Depending on the spatial structure of input and the degree to which spike output depends on phase of oscillation, external input can evoke different types of patterns, reproducing both zero-lag and spatial-lag observations of complex spike synchrony. We present examples of both states in a spiking network model of olivary neurons, and discuss the potential functional role of these structured spiking states in terms of generating motor timing and filtering error signals. Mammals selecting actions in noisy contexts quickly adapt to unexpected outcomes to better resolve uncertainty in future decisions. Such feedback-based changes in behavior rely on plasticity within cortico-basal-ganglia-thalamic (CBGT) networks, driven by dopaminergic (DA) modulation of cortical inputs to the direct (d) and indirect (I) pathways of the striatum. DA error signals favor the D pathway over the I pathway for rewarding actions with the opposite tendency for aversive ones, effectively encoding the values of alternative actions. It remains unclear how changes in action value influence the mechanisms of the action selection process itself. Here we use a biologically plausible spiking model of CBGT networks to illustrate (1) how feedback-driven DA signals modify the strength of D and I pathways in accordance with a simple reinforcement learning model and (2) how asymmetries in D/I efficacy, resulting from the learning process, impact the accumulation of evidence for alternative actions. Simulations of corticostriatal synapses showed that DA feedback leads to asymmetrical weights in the D and I pathways within a given action channel and the ratio of these weights (w_D/w_I) effectively encodes the action's expected value (Q). We then simulated the full CBGT network in the context of a simple 2-choice value-based decision task under different weighting schemes for cortical inputs to the D and I pathways (high, medium, and low w_D/w_I) for one of the action channels. The simulated response times from these simulations were fit with two variants of a drift-diffusion model (DDM), leaving either the drift-rate or the boundary height free to vary with the w_D/w_I ratio. As w_D/w_I increases, the speed of information accumulation in the decision process also increases, providing a direct mapping between network level properties of CBGT systems and cognitive decision processes. Finally, we have incorporated the corticostriatal plasticity module into the CBGT network model to form an integrated learning and decisionmaking network. Fits of the DDM to integrated network outputs will provide novel predictions about the mapping between CBGT and DDM parameters-drift-rate, boundary height, accumulation onset time, bias, and others-that best captures RTs associated with variable reward schedules in human experiments performed in our lab. This framework also allows us to explore how particular basal ganglia network features, such as tonic dopamine levels and changes in synaptic connection strengths, relate to changes in decision-making strategies, including those driven by behavioral parameters such as expectation and motivation.

P256
The mean-field theory of dynamically balanced neuronal networks Takashi  Coherence among activities of individual neurons and local field potential oscillations has been suggested as a clue to the mechanisms underlying information integration in the brain [1][2][3][4][5][6]. Experiments have also revealed that balanced excitatory and inhibitory synaptic inputs to neurons underlie the local-field potential oscillations [7]. However, despite recent pioneering studies of oscillations in neuronal networks [8][9][10][11][12][13][14][15], how the local field potential oscillations emerge as a result of balanced excitatory and inhibitory inputs and how individual neuronal activities become coherent with those oscillations remain to be understood theoretically. In the present study, we investigate a simple neuronal-network model on a dynamical balance between excitatory and inhibitory recurrent inputs, developing an analytical method that extends a previous theory [16] and describes this type of networks theoretically for the first time see [17] for a preprint,. In this network, microscopic dynamics of a small number of neurons are amplified by the strong excitation and inhibition and reflected in the macroscopic dynamics of the mean synaptic input over the network which have been considered as the origin of local field potentials. Conversely, the macroscopic dynamics of the mean synaptic input constrain the microscopic fluctuations in the activities of individual neurons. As a result of these bidirectional interscale interactions, oscillatory patterns of the mean synaptic input similar to local field potential oscillations spontaneously emerge. As the magnitude of balanced excitation and inhibition is increased, the mean synaptic input and the neuronal activities become coherent. This type of coherent states can also be induced by applying external stimuli to a small number of neurons in the network. The above behaviour of the network model is predicted by our theory with good quantitative agreement between the theory and direct simulations. Numerical results further suggest that the coherent states allow selective and reproducible read-out of information from the network. In conclusion, our results suggest a novel form of neuronal information processing that accounts for the emergence of local field potential oscillations, their coherence with neuronal activities, and the role of coherent dynamics in information processing in the brain. We also expect our results to provide a foundation for designing artificial neuronal networks for reservoir computing and beyond.
is routinely collected using a variety of techniques. One such comparably easy-to-perform technique is measurements of extracellular potentials by insertion of electric probes into neural tissue. However, the interpretation of the low frequency part of the signal, the local field potential (LFP), is hard because the measured signals result from both local and remote neural activity. Applications of CNNs for LFP analysis are not yet widespread, in particular with regards to detecting and classifying neural events or states that may not readily be detected using conventional methods. Here, we ask the question: Can CNNs be trained to estimate the underlying model parameters of spiking neuron networks from the LFPs they generate? We apply a recently developed hybrid scheme for computing extracellular potentials from spiking point-neuron network models [2] to a cortex-like, sparsely connected network model consisting of one excitatory and one inhibitory population of leaky integrate-and-fire (LIF) neurons. The network is simple enough to allow for detailed analysis of its state space [3]. We systematically vary different network parameters (for example connection strengths and amount of external input), run each simulation and compute the corresponding 'virtual' LFP signals as if measured at different depths through the neuronal populations. We then train CNNs set up using Tensorflow on subsets of LFP data, and explore to what extent model parameters can be estimated by the CNNs for the remaining LFP data. We indeed find that these CNNs can, based on the generated LFP, accurately identify the model parameters underlying the simulations by this relatively simple spiking network. This work contributes to a better understanding of what information is available in the LFP signal. It is also a first step in the direction of new analysis methods applicable to experimental LFP data that can be used to obtain more detailed information about the underlying neurons and neural networks. As multimodal sensory information proceeds to the cortex, it is intercepted and processed by the nuclei of the thalamus. The main source of inhibition within thalamus is the reticular nucleus (TRN), which collects signals both from thalamocortical relay neurons and from thalamocortical feedback. Within the reticular nucleus, neurons are densely interconnected by connexin36-based gap junctions, known as electrical synapses. Electrical synapses have been shown to coordinate neuronal rhythms, including thalamocortical spindle rhythms, but their role in shaping or modulating transient activity as is propagates through the brain is far less understood. Our recent findings on the plasticity of electrical synapses led us to investigate the impact of changes in electrical synapse strength on the circuit that embeds them. We constructed a four-cell single-compartment Hodgkin-Huxley model comprising thalamic relay and TRN neurons, and used it to investigate the impact of electrical synapses on closely timed inputs delivered to thalamic relay cells (Fig. 1). We showed that the electrical synapses of the TRN assist in cortical discrimination of these inputs through effects of truncation, delay or inhibition of thalamic spike trains. Electrical synapses lead to increased thalamocortical spiking separation and independence when inputs to the thalamus are dissimilar in strength and arrival timing. Conversely, electrical synapses can result in fusion of thalamocortical spiking by masking smaller strength and/or temporal differences between the incoming synaptic inputs. Thus, electrical synapses within the thalamocortical circuit strongly influence what the cortex receives from thalamus, and whether the information relayed to cortex is amenable to discrimination. We expect that these are principles whereby electrical synapses play similar roles in regulating the processing of transient activity in excitatory neurons across the brain. Our simulations provide specific predictions regarding the impact of electrical synapses and plasticity in thalamocortical processing, which we are currently testing in experiments in vitro.

Electrical synapses between inhibitory neurons shape the responses of principal neurons to transient inputs in the thalamus: a modeling study
Introduction: Traumatic brain injury (TBI) occurs when an external force results in structural damage to the brain, typically in white matter regions. In cases of mild to moderate TBI, this damage often goes undetected with conventional imaging techniques. However, since the structural damage invokes axon shearing, magnetic resonance imaging (MRI) with the application of network science, may improve detection and ultimately uncover the underlying dysfunction in TBI. Furthermore, research using Diffusion Tensor Imaging (DTI) has shown that diffusion properties, as well as connectivity patterns, can depict properties of TBI networks, namely, decreased fractional anisotropy (FA), increased mean diffusivity (MD), higher small-worldness, higher modularity, and lower global efficiency. While these metrics provide insights into the properties of the structural network, the specific attributes of the network that are disrupted after TBI are still unknown.
Here, we aim to further advance the knowledge about the spatial attributes of TBI-related network dysfunction by applying novel network science methods.

Methods:
The study enrolled 29 controls and 43 TBI patients who underwent an MRI scan (sagittal T1-weighted volumes and axial diffusion-weighted volumes acquired in 32 directions) within 4 ± 2 days from injury. The Human Connectome Project Multimodal Parcellation, providing 180 regions per hemisphere, was utilized for node definitions in the structural graph with each region corresponding one node. Further resolution of these graphs was achieved at 2-fold, 5-fold, and 10-fold splitting, via k-means with biological constraints, of each region. Edges in the structural graph were represented by streamlines seeded from each white matter voxel in the cerebrum, thresholded by anisotropy and curvature, and calculated using probabilistic Bayesian tractography and deterministic FACT algorithm tractography. Streamlines were retained if two different nodes were connected, the connection included the seed, and was at least 10 mm. The resulting edge definitions for weighting the structural graphs include: streamline counts, mean FA, along with corrections for node volume and streamline length. Adjacency matrices were constructed using the aforementioned node and edge definitions. These matrices were analyzed for graph measures of segregation, integration, and influence with subsequent group analysis via participation coefficient. The final analysis applied spatial machine learning algorithms for assessing network dysfunction. Results: Previous research has shown sensitivity in results from structural network construction. Here, we comprehensively construct a variety of graphs for each subject and utilize each graph as a valid representation of the structural network. First, the diffusion properties in TBI subjects showed similar patterns for alterations in FA and MD, and specific track related decreases in FA will be presented. Second, graph metrics for segregation, integration, and influence show interesting changes in the acute TBI case, most notably were changes in network efficiency. Next, feature extraction was implemented to find indications of disconnections in the TBI structural network followed by spatial machine learning algorithms to show spatial attributes of these network difference. The results present a novel step towards understanding the structural network dysfunction in acute mild to moderate TBI. A distinguishing feature of synaptic and endocrine secretory vesicle exocytosis is the high degree of variability in response latency. Stochastic Ca2+ channel gating is a major source of this stochasticity [2][3][4]. Another reason for high variability is that only a small number of Ca2+ ions enter the cell through a single channel during an action potential, as well as the stochasticity in the Ca2+ binding to Ca2+ buffers and sensors. This leads to a widely-held assumption that solving mass-action reaction-diffusion equations for buffered Ca2+ diffusion does not provide sufficient accuracy for modeling Ca2+-dependent cell processes. However, several recent comparative studies show a surprising close agreement between deterministic and trial-averaged stochastic simulations of Ca2+ diffusion, buffering and binding [1,2], as long as Ca2+ channel gating is not strongly Ca2+ dependent [2][3][4]. We present further analysis and comparison of stochastic and mass-action simulations, focusing on Ca2+ dynamics downstream of Ca2+ channel gating. Smoldyn (www.smold yn.org) is used for stochastic simulations, while CalC (Calcium Calculator) is used for deterministic simulation (www.calci umcal culat or.org), with Ca2+-binding scheme of Felmy et al. [5]. We show that the discrepancy between deterministic and stochastic approaches can be surprisingly small even when only as few as 30 ions enter per single channel-vesicle complex, assuming that Ca2+ binding depletion is accurately implemented, as described previously [1,2]. We argue that the reason for the close agreement between stochastic and mass-action simulations is that the discrepancy between the two approaches is determined by the size of correlations between reactant molecule numbers rather than their variance. Further, contrary to naïve intuition, mass-action reactiondiffusion description provides an estimate of the vesicle release latency (i.e. first passage time) probability density function, since it is identical to the system of equations for the first moments of the Ca2+-sensor system derived from the underlying master equations under the moment-closing assumption of negligible correlations (Fig. 1). This may explain the close match between stochastic and mass-action simulations of Ca2+ dependent exocytosis, despite the high variability in Ca2+ diffusion, buffering and binding. Potential explanations for the small size of molecule number covariances will be discussed. Naturalistic odour stimuli have a rich temporal structure. It has been hypothesised that this structure contains information about the olfactory scene, for example the distance to an odour source [1,2]. Furthermore, it has been suggested that animals might exploit this structure and extract this information in order to find odour sources [3]. As some of this information may lie in the frequency content of the stimuli [2], we studied input frequency dependent responses of mitral cells (MCs) in the olfactory bulb (OB), the first processing stage in the mammalian olfactory system. Specifically, we investigated whether MCs show frequency tuning and, if they do, how different components of the glomerular layer circuitry shape and determine the tuning. We used a model of the OB (modified from [4]) containing periglomerular cells (PGCs) and MCs, thus focusing on the recurrent and feed-forward inhibition in the glomerular layer. Simple sinusoidal currents of varying strengths and frequencies were used as input to the model. We constructed frequency tuning curves, extracted the peak resonance frequencies and looked at how these changed for different parameter combinations. We also considered the strength of the tuning, measured as (max firing rate − mean firing rate)/mean firing rate. We found that the resonance frequency decreased as the excitation of PGCs (both from the input and from the MCs) increased, whereas the strength of the PGC inhibition onto MCs did not seem to have a strong effect. Furthermore, the resonance strength increased with the strength of the excitatory connection between MCs and PGCs when the PGCs received sufficient external input from olfactory stimuli. These results suggest that the MCs can indeed show frequency tuning and that this depends on the strength of the excitatory synaptic input to PGCs, which provide inhibitory input to the MC. However, the observed frequency tuning occurred in a narrow range (19.5-33.0 Hz). Future work should investigate how the OB could use this frequency tuning to obtain information about the surrounding olfactory scene.

P264
The combined effect of homeostatic structural and inhibitory synaptic plasticity during the repair of balanced networks following deafferentation Although a number of previous experimental and theoretical studies have investigated network reorganisation following deafferentation down to the level of synaptic elements [4], the mechanisms that are involved in this process are still not completely understood. We examined the dynamics of the repair mechanism by incorporating activity dependent homeostatic structural plasticity [1] into a spiking neural network model balanced by inhibitory synaptic plasticity [6]. Results from our simulations suggest that the process of reconfiguration of lateral connectivity following sensory deprivation is extremely sensitive to the balance of excitation and inhibition (E-I) in the network. We find that while fast homeostatic inhibitory synaptic plasticity is able to re-establish the E-I balance in neurons outside the lesion projection zone (LPZ), it prevents them from transferring excitatory activity to the deprived neurons in the LPZ. On the other hand, uncontrolled disinhibition by suppression of homeostatic inhibitory synaptic plasticity initially allows deprived neurons to regain activity but fails to stabilise the network back to a functional balanced state. These observations are in accordance with findings that indicate that inhibition plays a critical role in network rewiring [2] seemingly by stimulating structural plasticity mechanisms seen during development [5]. The sprouting of inhibitory axons outwards from the LPZ, opposite to excitatory axons has also been observed, possibly to re-inhibit neurons outside the LPZ [3]. Therefore, we hypothesise that the ratio of excitation and inhibition must follow a specific trajectory in the different regions of the network to enable successful repair as has been observed in various studies. The model of structural plasticity implements the dynamics of synaptic elements as dependent on intrinsic properties of individual neurons only [1]. The configuration of the network, by the formation and removal of synapses therefore depends solely on the numbers of various synaptic elements. Our current work extends this model by considering other factors that affect network rewiring, such as the activity dependent stability of synapses, and inhibition gradient guided axonal sprouting [4], to build a more faithful simulation of the underlying dynamics. This will enable us to study the effects of network reorganisation after deprivation on its computational functions, such as associative memory.
Since synchronized neuronal activity might underlie efficient communication in the brain, alterations thereof, as found in EEG/MEG studies of schizophrenic patients, might contribute to the symptoms characterizing schizophrenia [1]. A robust finding is a deficit in the gamma band auditory steady-state response (ASSR) [2]. Fast-spiking PV+ interneurons seem to be a major contributor to gamma oscillations. However, this class of inhibitory interneurons can be divided into at least two subgroups: basket cells (BCs) and chandelier cells (ChCs) [1]. Interestingly, for both subtypes cellular/molecular alterations have been identified in schizophrenia [1]. However, the role these two subgroups play during the generation of gamma oscillations, and during abnormal oscillations in schizophrenia, remains unresolved. We use a simple model, consisting of three populations of theta neurons: (1) pyramidal cells (PCs), (2) BCs and (3) ChCs (based on [3,4]). We assume that the prolonged GABAergic decay time at ChC synapses is a major contributor to gamma and beta band ASSR deficits in schizophrenia [3] and model this by increasing the decay time constant for ChCs. We then explore the model behaviour in response to oscillatory inputs in the beta and gamma range, for different ratios of BCs vs. ChCs (BCs are known to be more numerous than ChCs [1]), different strengths of inhibition of the ChCs onto PCs (ChCs might exert powerful inhibition because of their synapses directly targeting the axon hillock of PCs [1]) and reductions in the strength of inhibition of BCs (a possible result of genetic alterations in schizophrenia [2]). At realistic BCs/ChCs ratios, increased ChC inhibition, due to increased decay times is not sufficient to strongly reduce gamma power as it has been described for schizophrenia patients. Under the assumption that they exert much more powerful control over PC firing stronger reductions were observed. However, the model did not reproduce other deficits that have been described in schizophrenia, such as an increase in beta power for 20 and 40 Hz stimulation [3]. Simultaneously reducing BC inhibition did not change this overall behaviour. Interestingly, prolonged decay times at BC-PC synapses led to both a strong decrease of gamma power and an increase in beta power, matching experiments more closely. We conclude that changes in the dynamics at ChC-PC synapses might not be a major contributor to gamma and beta band ASSR deficits in schizophrenia. Our results suggest that the more numerous BCs are likely to dominate the influence inhibitory interneurons exert on the PC population during oscillatory entrainment.

P266
Modeling rod-cone parallel processing in the retina Parallel processing underlies computation in many neural circuits. Several common circuit motifs control how parallel processing contributes to circuit function: (1) divergence of common inputs to parallel circuits; (2) distinct linear shaping of signals in different parallel circuits; and (3) location of key circuit nonlinearities relative to the convergence points of signals from different parallel circuits. Interactions between rod and cone mediated signals in the retina provide an excellent opportunity to investigate these computational elements. Vision relies on inputs from both rod and cone photoreceptors across light conditions ranging from moonlight to dawn, and visual perception is strongly influenced by interactions between the resulting signals. To understand how retinal mechanisms contribute to these perceptual interactions, we aim to develop a model that predicts retinal output in response to temporally and spatially modulated images in dim and intermediate light. We will use direct retinal recordings from cells across the primate retina to constrain the architecture of our model and test its ability to capture key features of rod-cone interactions. The model will then be used to predict neural responses to novel stimuli, specifically focusing on identifying stimuli that highlight the importance of specific circuit features in shaping retinal outputs; significant discrepancies between predictions and empirical measurements will be utilized in finetuning the model. Ultimately, this model will improve both our understanding of how perceptually-relevant computation operates in parallel circuits and our ability to incorporate relevant computational features into devices (e.g. retinal prosthetics) that aim to replicate retinal function. Avalanches have been suggested to reflect a scale-free organization of cortex. It is hypothesized that such an organization may relate to a particularly effective form of activity propagation which is balanced between failure (activity fails to reach a target area) and overactivation (activity reaches a target area via many routes leading to wasted activity or epileptiform activity patterns). We electrically stimulated a computer model of mouse primary motor cortex (M1) and analyzed signal flow over space and time. Initially we stimulated a 300 μm × 600 μm slice of M1 using a 10 μm × 10 μm 0.5 nA stimulus across all 6 layers of cortex (1350 μm) for 100 ms. Waves of activity swept across the cortex for a half a second after the end of the electrical stimulus. We extracted avalanches from the data by counting events, spikes, occurring within 1 ms frames. An avalanche of length N was defined as N consecutively active frames, preceded by a blank frame and followed by a blank frame. A graph of the cortical slice above, with the 0.5 nA stimulus, displayed a bimodal distribution. We observed 18 avalanches in total with 4 single neuron avalanches and all the other avalanches containing more than 1000 neurons each. The largest avalanche contained 7000 neurons. Studies have generally shown avalanche activity to show a linear log-log graph starting highest from small avalanches and decreasing as the avalanches get larger. We looked at responses of M1 to lower amplitude stimuli between 0.05 and 0.5 nA to see if they may fit a classic inverse power-law curve. We graphed M1 response to a 500 ms electric stimulus at various amplitudes and found particularly clear inverse power-law responses to stimuli between 0.16 and 0.18 nA. In the 300 μm × 300 μm slice of M1 for 500 ms using 0.16nA we observed 90 avalanches from as small as a single neuron action potential in isolation to 13 neurons spiking. A large proportion were SOM neurons participating in the avalanches but they also included IT neurons at this level of stimulation. The simulation of avalanches in cortex offers advantages for analysis that are not readily done experimentally in in vivo or in vitro. We have been able to record from every neuron in our M1 slice and follow activity from cell to cell. In the future we will analyze how avalanches take place within and between layers.

Simulation of avalanches in mouse primary motor cortex (M1)
Experimental data is accumulating at an unprecedented-and accelerating-rate. However, as the BRAIN Initiative 2025 report points out: (1) even excellent quality data will not yield solid conclusions unless it is adequately integrated and interpreted, and (2) turning experimental knowledge into understanding inevitably requires rigorous theory and modeling. Biophysically realistic modeling provides a tool to integrate, organize and bridge data at multiple scales and develop hypothesis about the biological mechanisms underlying physiological and pathological brain function. NetPyNE (www.netpy ne.org) is a high-level Python interface to the widely used NEURON simulator [1]. It provides high-level declarative language designed to facilitate the definition of data-driven multiscale models, e.g., a concise set of connectivity rules vs. millions of explicit cell-to-cell connections. The user can then easily generate NEURON network instances from these specifications, run efficient parallel simulations (with predefined setup for supercomputers), and exploit the wide array of built-in analysis functions (e.g. connectivity matrix, voltage traces, raster plot, information transfer measures). A recent feature provides the ability to place extracellular LFP recording electrodes at arbitrary 3D locations and plot the LFP signal, power spectra or spectrogram. All this functionality is also accessible via a graphical user interface (GUI) based on the state-of-the-art Geppetto technology. The GUI provides an intuitive way to define the model, including an interactive Python console and full synchronization with the underlying Python-based model (Fig. 1). The user can visualize the 3D network, run simulations and choose from the available analysis plots. NetPyNE's standardized format clearly separates model parameters from implementation and can be exported/imported to NeuroML, thus making it easier to understand, reproduce, reuse and share models. This has motivated the conversion of several published models to NetPyNE specifications, including the Potjans&Diesmann cortical network, the Traub thalamocortical network, the Cutsuridis CA1 microcircuit and the Tejada dentate gyrus network. The tool is also being used to develop a variety of new models exploring mouse M1 microcircuits [2], the claustrum network, cerebellum circuits, transcranial magnetic stimulation (TMS) in cortex, or the underlying biophysics of EEG recordings. We expect the NetPyNE tool to make data-driven biophysically-detailed network modeling accessible to a wider range of researchers and students, including those with limited programming experience, and encourage further collaboration between experimentalists and modelers. Fig. 1 The NetPyNE GUI provides an intuitive way to generate, simulate and analyze data-driven biophysically detailed network models. NetPyNE GUI showing 3D representation of (simplified) M1 model and simulation results (raster plot, statistics, voltage traces and power spectra) and 3D intracellular reaction-diffusion models [1]. This has been used to probe intracellular calcium dynamics in both physiological and pathological conditions. Originally rxd provided only limited extracellular support with an isolated space around each segment. We have extendedrxdto include coarse-grained macroscopic models of the extracellular space [2]. NEURON thus allows detailed neuron models to be embedded in a 3D macroscopic model of tissue. Extracellular diffusion is implemented using the Douglas-Gunn alternating direction implicit method, an efficient scheme which supports parallelization. Reactions are now implemented using Just-In-Time compilation, allowing numerical integration to use faster compiled code rather than slow interpreted code. The macroscopic tissue model is based on a volume averaging approach, allowing the user to specify both the free volume fraction (the proportion of space in which species are able to diffuse) and the tortuosity (the average multiplicative increase in path length due to obstacles). These tissue characteristics can be spatially dependent enabling the modeler to account for differences in brain region or pathological effects of injury. We applied the rxd simulation framework to develop a model of ischemic stroke, which required multiscale coupling of electrophysiology with intracellular molecular alterations, and consideration of network properties in the context of bulk tissue alterations mediated by extracellular diffusion [3]. We initially modeled spreading depression triggered by elevated potassium in a cube of tissue (Fig. 1). Occlusion of a blood vessel in the brain triggers a cascade of changes, including: (1) synaptic glutamate release, related to excitotoxicity; (2) elevated extracellular potassium, leading to spreading depression; (3) cell swelling, reducing the extracellular volume and increasing the tortuosity; (4) production of reactive oxygen species, which give rise to inflammation. These cascades occur over multiple time-scales, with the initial rapid changes in cell metabolism and ionic concentrations triggering several damaging agents that may ultimately lead to cell death.
The NEURON simulator (neuron.yale.edu) provides a computational framework for studying not only networks of neurons but also the interplay between electrophysiology and chemical dynamics, (both intracellular and extracellular reaction-diffusion models). The models underlying these studies can be specified, simulated, and analyzed using both Python and graphical tools. NEURON's graphical tools previously focused on supporting pure electrophysiology models. We describe a new integrated graphical toolset, powered by wxPython 4.x, for specifying and visualizing NEURON models incorporating both reaction-diffusion dynamics and traditional electrophysiology simulation. In comparison to electrophysiology models, these models feature new types of regions (1D and 3D, intracellular organelles, extracellular space, etc.), new types of kinetics, etc. Our toolset includes an expanded RxDBuilder supporting recent enhancements to NEURON's reaction-diffusion capabilities, including extracellular and 3D intracellular simulations. The intracellular 3D graphical tools provide a detailed view of the cells morphology, enabling the modeler to select a region of interest over which to plot relevant intracellular concentrations.
With the extracellular space, the GUI allows the modeler to choose to view the concentration dynamics for: a single voxel in the extracellular space, an average around the cell of section of interest, or over the whole extracellular space. To allow model changes from both the console and the GUI, the graphical tools are run in a separate thread that periodically polls the internal state; a function is provided to allow arbitrary wxPython windows to be run in the same thread, allowing user customization. For performance reasons, state variables are recorded in C++ during simulations; visualization occurs via Python at a user-specifiable interval. A session consisting of; the models, their current state and the graphical tools may be saved and loaded for future reuse. We demonstrate the utility of these model construction and visualization tools with a 3D intracellular calcium wave model and an extracellular model of spreading depression.
In recent years, nonlinear Hawkes processes implemented as point process Generalized Linear Models (GLMs) have been proven to be a useful tool for analyzing microelectode array recordings. On a more theoretical side they are related to well known models of neuronal dynamics such as the spike response model and can capture the spiking temporal patterns of Izhikevich canonical models. Unlike ODEbased neuron models, point process GLMs can be fitted directly to the spike times data. For most relevant models, the models can be easily fitted using standard optimization tools, as the likelihood function is strictly convex. Despite the acknowledged utility of nonlinear Hawkes process GLMs, the dynamics of fitted models has attracted attention only recently. In particular, simulation of fitted models can often produce unphysiologically high firing rates, despite passing many goodness-of-fit tests. Here, "unphysiologically high rate" means "rate, close to 1/absolute refractory period, " reflecting "runaway excitation". To make nonlinear Hawkes process GLMs useful for long-term prediction of neuronal activity and simulation studies, it is important to understand which model features can lead to runway excitation. The mathematical theory of nonlinear Hawkes processes is not fully developed. Prior studies (e.g. Bremaud and Massoulie, Ann Prob, 1996) have focused either on the mere existence of finite stationary firing rates (and do not allow to estimate their actual values in general), or on the theoretical examination of infinite neuronal networks with some degree of homogeneity in the connectivity, whereas actual recordings typically contain not more than several hundreds of neurons and lack homogeneity. The question we consider here is how to predict runway excitation for an arbitrary finite network of Hawkes processes. Several recent theoretical approaches, based on statistical physics, allow to approximate stationary firing rates for nonlinear Hawkes processes. Those include mean field approximations, 1-loop approximation (Ocker et al., PLoS CB, 2017), quasi-renewal (QR) approximation (Gerhard et al., PLoS CB, 2017) and the regular firing rate test. These approaches are quite different conceptually, were introduced in different settings and have limitations in different directions. E.g. mean field approximation does not work for neurons with absolute refractory period without additional adjustments; 1-loop approximation inherits the same issue but also shows poor accuracy for strong nonlinear functions (at least for some networks), whereas QR approximation is primarily designed to work for exponential nonlinearity only. Moreover, The QR approximation can lead to prediction of multiple "fixed points" (stationary firing rates) that may relate to the actual dynamics in a nontrivial manner. To summarize, so far the strengths and limitations of these different approaches have not been compared systematically. Also their application to real data has been limited. We present a study that compares how the above approaches work for simple single and multiple-neuron nonlinear Hawkes process GLMs and compare their predictions with simulations. We identify model features that make some approaches work much better than others. We show that, in some cases, the different approaches can complement each other. Finally we demonstrate how the different approaches work when being applied to multivariate nonlinear Hawkes process GLMs fitted to actual spiking data. Methods to reliably predict seizures in patients with epilepsy have been sought for several decades. Reliable seizure prediction would have a major impact in the quality of life of people with pharmacologically intractable seizures, allowing for new seizure prevention therapies based on warning and closed-loop electrical stimulation systems [1]. While most seizure prediction systems have relied upon EEG and/ or ECoG, the predictive value of intracortical neural signals remain little explored [2]. Here, we demonstrate that seizures can be predicted early in advance from the neural activity of small neocortical patches distal from the identified seizure onset areas. We used multiunit activity and local field potentials recorded via microelectrode arrays (Blackrock Microsystems, Salt Lake City, Utah) plus machine learning algorithms to show that interictal and preictal activity in people with focal seizures can be discriminated. Intracortical signals were recorded in 5 patients undergoing neuromonitoring for resective surgery from a neocortical area distal to identified seizure onset areas [3]. Preictal periods were defined as the one-hour period leading to a seizure with a 5-minute interval between the preictal period and the seizure onset time. Interictal periods excluded the four hours preceding any seizures [2]. This setting attenuates potential errors and uncertainty in the determination of actual seizure onset times and the separation or interictal and preictal periods. Long short-term memory (LSTM) recurrent neural networks were used to assess the predictive power of the different features extracted from the recorded neural activity signals. Substantial predicted power, as assessed by the area under the receiver operating characteristic curves, was achieved with a 90% score for at least one type of feature in each patient. Importantly, we show that successful prediction can be achieved based exclusively on the multiunit activity of recorded neurons detected by thresholding high-pass filtering the electric potentials. This result indicates that neural activity in the recorded local neocortical patch exhibited preictal changes not only in subthreshold postsynaptic potentials that could be driven by the distal epileptogenic areas, but also changes in the local neuronal spiking activity in the recurrent neocortical networks.

Multiunit activity patterns in neocortex predict upcoming seizures in human focal epilepsy
Our findings indicate that large-scale neuronal networks are engaged beyond the identified epileptogenic seizure onset areas towards the onset of a seizure, and open new perspectives for seizure prediction and control by emphasizing the contribution of multiscale neural signals in these networks. Long-range communication in the nervous system is carried out with the propagation of action potentials along the axons of nerve cells. While typically thought of as being unidirectional, it is not uncommon for axonal propagation of action potentials to happen in both directions. This is the case because action potentials can be initiated at multiple 'ectopic' positions along the axon. Axons are endowed with ionotropic and metabotropic receptors for transmitters and neuromodulators that can alter membrane excitability, and initiate ectopic action potentials [1]. Action potentials generated at distinct sites, and traveling toward each other, will collide. Recently, it has been suggested that some biological axons may show properties of crossing action potentials [2], and that Hodgkin-Huxley type models may be inadequate for representing some axons. However, this view has been challenged in a subsequent study using a reduced Hodgkin-Huxley model [3]. As neuronal information is encoded in the frequency of action potentials, the rate of action potential collision and annihilation may affect the way in which neuronal information is received, processed and transmitted. Additionally, action potential collision and annihilation can be of relevance in the treatment of spinal cord injury along with chronic pain of peripheral origin [4]. Here we present numerical simulations and experimental results aimed at helping to elucidate the subject of colliding action potentials [5]. We introduce an axonal multicompartmental model with the compartments represented by Hodgkin-Huxley equations reciprocally connected to each other by diffusive coupling. The numerical simulations are capable of mimicking low frequency ectopic spiking with orthodromic and antidromic action potential propagation. They predict that colliding action potentials traveling in opposite directions annihilate and do not cross. We further discuss this matter in the context of axonal excitability and supernormality in the wake of action potential generation for neurons of type I and type II. We also present results of experimental work performed on the earthworm ventral cord and on the crustacean stomatogastric nervous system. Both numerical simulations and experimental outputs clearly and unambiguously demonstrate that annihilation is inevitable. Consolidation of synaptic changes in response to neural activity is thought to be fundamental for memory maintenance over a time scale of hours. In experiments, synaptic consolidation can be induced by repeatedly stimulating presynaptic neurons. However, the effectiveness of such protocols depends crucially on the repetition frequency of such stimulations and the mechanisms that cause this complex dependence are unknown.

Optimal stimulation protocol in a bistable synaptic consolidation model
Here we propose a simple mathematical model that allows us to systematically study the interaction between the stimulation protocol and synaptic consolidation. We show the existance of optimal stimulation protocols and we explain the results using phase-plane techniques. Our model explains why the temporal structure of LTP induction protocols, in particular the repetition frequency, is important for the outcome of plasticity experiments (Fig. 1). Biological neurons exhibit an extended richness of biophysical mechanisms besides passive input integration. Among these, spike frequency adaptation has received great interest due to its richness in time scales [1] that allows sensory neurons to optimally transmit information [2]. Understanding the effect of such history-dependent processes on the dynamics of a recurrent neural networks has proven to be a hard task. Recent developments in mean-field approaches in the presence of such history-dependent processes [3,4] allow to compute the mean and fluctuations of the network activity. To obtain the temporal structure of the recurrently generated fluctuations however, one has to solve the system self-consistently for the fluctuations. This can be done, in the large-N limit of non-adaptive, randomly connected networks of rate units, using Dynamical Mean-Field Theory (DMFT), thanks to which it was first shown the existence of a quiescent phase and a chaotic phase in the network dynamics [5]. However, this technique was until now restricted to non-adaptive networks.

Oscillations and chaos in adaptive neural networks
Here we apply DMFT to a randomly connected network of adaptive rate neurons. The technical challenge emerging from this setting is that the resulting mean-field system is two-dimensional, causing standard DMFT techniques to be not applicable. We propose an iterative method that allows fast computation of the mean power spectral density of the network activity. We show that in a large portion of the adaptation parameter space, the dynamics of the adaptive neural network is qualitatively different from the non-adaptive one. Besides the purely chaotic and purely quiescent phases, the adaptive network features two new phases. For strong recurrent connectivity and strong adaptation, the chaotic dynamics exhibit a finite-width peaked power spectral density, which means that in the DMFT limit the system can be described as a stochastic oscillation. For lower connection strength, a bistable phase also emerges, in which a stable fixed point coexists with limit cycles, even in the large-N limit. Finally, we extend the wellknown result for the eigenvalue spectrum of Gaussian random matrices [6], to the adaptive case. This allows us to compute the stability of the zero fixed-point, and we show that it approximately predicts the separation between pure chaos and stochastic oscillations and the oscillation frequency at the criticality. The additional dynamical richness of adaptive neural networks could explain the better performance in learning tasks that require integration of information over long time scales [7]. This could highlight a novel important role of adaptation that arise through network-level interactions.  . The model includes spatial integration from V1 to MT, V1 iso-orientation surround suppression (IOSS) and normalization stages in V1 and MT. First, we implemented single-stage gain adaptation on the raw motion energy signals in twelve direction channels at each spatial location. The adapted signals are then used to compute the surround suppression signal and a spatially local classical untuned normalization signal. The normalized signals pass through a V1 opponency stage before being normalized and integrated (across space&direction) at the MT stage to form CDS and PDS units. Second, we implemented a recently proposed form of adaptation in which normalization weights between units are updated by a learning rule that aims to achieve pairwise response-product homeostasis (RPH; Westrick et al., 2016, J Neurosci), extending this mechanism from orientation to direction channels. Results: For single-stage gain adaptation, we found that the effect on V1 tuning of prolonged adaptation to a single direction depended on the size of the adapting and test stimuli (drifting grating patches optimized for the unit under study) in a manner qualitatively consistent with electrophysiological results in terms of response amplitude (Patterson et al., 2013, J Neurosci). However, tuning curves showed attractive shifts in the presence of untuned normalization when flank adaptation was limited to the classical receptive field (CRF), unlike repulsive shifts reported in the literature. When untuned normalization was omitted, there were no attractive shifts and no shifts in tuning for MT CDS cells, but there were repulsive shifts for PDS cells, contrary to the literature. For RPH normalization, we were able to achieve the desired repulsive shifts in V1 direction tuning for flank adaptation, but direction channels not driven by the adapter showed implausibly large increases in gain.

Conclusions:
Our results so far suggest that simple combinations of mechanisms believed to be fundamental to processing along the V1-to-MT pathway are not sufficient to account simultaneously for the influences of adaptation across a diverse set of stimulus configurations in the CRF and surround. To remedy this, we are implementing twostage gain adaptation and exploring alterations to RPH normalization that can better account for physiological data.
How features encoded by neurons are reorganized during development is not well understood. We examined this question in the context of visual processing in primary visual cortex (V1) of mice. Recent experimental studies on mice suggest the distribution of preferred orientations and directions of neurons in primary visual cortex (V1) are balanced after eye-opening, through an activity-dependent mechanism [1,2]. At eye-opening, the orientation preferences of V1 neurons are highly biased towards the horizon axis; as time goes on, neuronal preferences to other orientations emerge via visual experiences. This developmental balancing is considered to be important for optimal processing of the visual environment [3][4][5]. To achieve balancing of preferred orientations, two general conditions must be met. (1) The encoded feature should be plastic to achieve balancing during development while stable enough to ensure normal visual function.
(2) There exists intrinsic signaling mechanisms that specifically guide neurons to develop under-represented orientations. We propose that inhibitory neurons in the network are capable to provide such mechanism. The first condition guarantees the onset of reorganization during visual development, while the second condition guarantees the reorganization process to be bias balancing. We proposed a computational model to unveil the underlying mechanism for the balancing of orientation preference in V1. We built a 2-layer network consists excitatory and inhibitory neurons. The first layer is composed of Poisson excitatory neuron, which provides tuned input from external visual stimuli. The second layer consists inhibitory neurons as well as two types of excitatory neurons: senior excitatory neurons, which have the initial biased preference to the horizon axis, freshman excitatory neurons, which have no preference before learning. All excitatory neurons receive feedforward currents conveying orientations of external inputs as well as recurrent connections between themselves. These connections are subject to spiking-time-dependent-plasticity (STDP), which shapes the tuning property of excitatory neurons. In the critical period, inhibitory neurons become mature, and through reciprocal connections, they mediate mutual inhibition between excitatory neurons, leading the freshman neurons to encode orientations unrepresented in the network. For instance, when the horizon direction is presented, because of the strong responses of senior neurons, inhibitory neurons become active, which suppress freshman neurons, so that freshman neurons are less likely respond to the horizon direction; on the other hand, when an underrepresented orientation is presented, inhibitory neurons are less active, and freshman neurons learn to encode this orientation via STDP. After onset of visual experience, the network undergoes a balancing process that equalizing the initial bias in distribution of preferred orientations. Our model successfully reproduced the balancing of orientation preferences in V1, and unveiled the important role of inhibitory neurons in the visual development (Fig. 1). The striatum is a structure of the basal ganglia that is critical for reinforcement learning. In the striatum, cholinergic tonically active neurons (TANs) are thought to gate the dopaminergic input to medium spiny neurons during their involvement in action selection and reinforcement. TANs exhibit a context-dependent pause in their activity, during which the dopamine (DA) concentration in the striatum varies to encode reward prediction error (RPE), i.e. the difference between the expected and obtained reward. Although this mechanism has been the subject of many experimental studies, the role of TANs in motor learning is not well understood. However, it is known that TANs generate a short burst in response to a stimulus, which is followed by a pause in TAN activity for several hundred milliseconds. During the pause, dopaminergic neurons modulate DA release into the striatum to encode the RPE and thus induce learning. After this pause, TANs return to normal tonic firing and striatal dopamine concentration stabilizes at its baseline levels.

Cholinergic modulation of reinforcement learning in the striatum
The duration of the TAN pause depends on dopaminergic inputs to TANs through activation of D2 receptors. During baseline tonic firing, TANs-being cholinergic-control the output of dopaminergic neurons by releasing acetylcholine (ACh) that binds to the nicotinic receptors of the latter. When a reward is presented, TANs receive a short stimulus from the thalamus. This short increase in TAN activity is then followed by inhibition via a slow after-hyperpolarization (sAHP) current, which lasts several seconds, inducing a pause in the tonic firing of TANs. Another current, the hyperpolarization-activated cation current (h-current) allows quick recovery from sAHP. The h-current in these neurons is down-regulated by DA via D2 receptors. Therefore, the TAN pause is produced by the slow AHP current, and the length of the pause is modulated by the faster h-current. Thus, encoding the RPE depends on the dynamic interactions between DA and ACh release mechanisms in the striatum. In this study, we constructed a mathematical model of ACh-DA interactions to clarify TANs' role in reinforcement learning. We fit our model to the data obtained in electrophysiological experiments. Furthermore, we integrated the model ACh-DA interactions into our previously published model of the reward-based motor learning during center-out reaching movements. We simulated the effects of the striatal dopamine deficiency as observed in Parkinson's disease patients. Additionally, we simulated and mechanistically explained the effects of administration of L-DOPAa common treatment for early phases of Parkinson's Disease-clarifying the mechanism by which L-DOPA recovers learning in these patients.
In simulations, our model shows that both the baseline DA concentration and phasic DA release positively correlate with the duration of the TAN pause. Therefore, in the case of striatal DA deficiency, the loss of learning is associated not only with lower DA concentration but also with a shorter TAN pause, which means there is a shorter time period for learning to occur. We simulated L-DOPA administration by increasing the baseline concentration of DA in the striatum, which did allow partial recovery of motor learning functionality even though the magnitude of phasic DA release was not affected. Our model explains this recovery by L-DOPA-mediated prolongation of the TAN pause, which increases learning efficiency. The respiratory and cardiovascular systems work together to oxygenate tissues and remove carbon dioxide and are physiologically integrated. Central neural circuits that control the respiratory and cardiovascular functions are located in brainstem and receive sensory feedback to maintain the gas homeostasis. Respiratory and cardiovascular physiologic outputs are partially synchronized/modulated by each other, and the respective brainstem neuronal networks have reciprocal synaptic connections. However, no quantitative mechanistic description was suggested to explain specific aspects of the cardio-respiratory interactions and their alterations in certain pathophysiological conditions. Two major markers of cardio-respiratory interactions were previously identified: cardio-ventillatory coupling (CVC) and respiratory sinus arrhythmia (RSA). CVC is usually interpreted as a form of partial synchronization between cardiac and respiratory rhythm which is characterized by varying probability of a heartbeat to occur at different phases of the respiratory cycle. RSA is a phenomenon concerned with changes in heartrate at different respiratory phases, which is usually represented by the dependence of the inter-heartbeat interval (R-R interval) on the respiratory phase with R-R interval shortened during inspiration and prolonged during expiration. Due to similar representation, CVC and RSA are often confused. However, there is substantial experimental evidence that independent mechanisms mediate the two phenomena. Here, we introduce a closed loop model of the integrated respiratory and cardiovascular control system to describe mechanisms for both CVC and RSA. The model combines and extends our previous data-driven models that incorporated mechanisms of cardiovascular input to the respiratory system or respiratory input to the cardiovascular system. In this model, CVC is mediated by the pulsatile inputs from arterial baroreceptors to neurons of the respiratory central pattern generator (rCPG) with pulses corresponding to the increases and subsequent relaxations in arterial pressure caused by heart contractions. We implement baroreceptor input to the rCPG Fig. 1 The emergence of balanced orientation selectivity in the network. The upper panel displays the distributions of neuronal orientation preferences at different developmental stages, which are: early stage (S1, just after eye-opening), middle stage (S2, during the development), and late stage (S3, mature). The lower panel displays the orientation tuning of the inhibitory neuron pool at different developmental stages. The inhibitory neuron pool is connected evenly to all excitatory neurons. At S1, its tuning property is determined by senior neurons. At S3, it no longer has preference, since the tunings of excitatory neurons cover the whole orientation space

Brainstem mechanisms of cardio-respiratory coupling
model. This success is demonstrated through (i) inferred spectral densities that are consistent with nonparametric spectral estimates using Welch's method-see the Fig. 1-and (ii) changes in inferred parameter values with anaesthesia level that are consistent with mechanistic effects of isoflurane [3]. These results open the way for future work wherein one may be able to estimate the (sub)cellular effects of psychoactive agents merely from recordings of electrocortical activity. Many neuroscientists have used computer simulation. The accuracy of computer simulation is affected by not only a computational model and a numerical method but also floating point precision used in the simulation. The floating point types in C language are float (single precision), double (double precision), and long double (extended double precision). Generally, without consideration, we will select the double precision floating point for computer simulation. In recent years, a large-scale simulation and a real-time simulation of the neural system has been extensively attempted. For a large-scale simulation and a real-time simulation, not only a supercomputer but also a desktop workstation are used. In a simulation using a desktop workstation, a graphics board or an accelerator board with a GPU achieves acceleration of the simulation. When we perform a simulation using a graphics board, computational time with the single precision floating point is shorter than that with the double precision floating point. Furthermore, using single precision floating point can also reduce data transfer time. If the floating point precision has little effect on the accuracy of a simulation result, we can use single precision without worry and perform an efficiently accelerated simulation. In this study, we investigate the effect of the single precision, the double precision, and the extended double precision floating points on the dynamics of the neuronal activity in the computer simulation.  Fig. 1a which contains 3 parameters. Namely, a gain of the Direct pathway G1, a gain of the Indirect pathway G0, and a time constant of the velocity storage mechanism (VSM) in the Indirect pathway H. This simple model can reproduce OKR eye velocity well. Several studies have reported that the cerebellum is involved in the adaptation of OKR gain defined as eye velocity/visual stimulus velocity [2]. However, its contribution to generation and adaptation of the Direct and Indirect components has been controversial. Presently, we performed cerebellectomy before or after OKR gain adaptation in goldfish and evaluated the model parameters.

Decomposing adaptable elements of optokinetic response into cerebellar and non-cerebellar contributions by modeling and cerebellectomy approach
Methods: Goldfish were gently fixated at the center of a white cylindrical water tank with eye coils binocularly sutured on the cornea around the pupils for eye position measurement. Visual stimulus was projected on the wall of the water tank, and rotated in the clockwise direction at 20 deg/s for 8 s and stopped for 8 s repeatedly to generate horizontal OKR. This stimulation was continued for 3 h to induce OKR gain adaptation. Two kinds of experiments were conducted: Acute and chronic cerebellectomy. In the former experiment, normal goldfish (n = 8) underwent the 3-hour OKR training, then the cerebellum was acutely removed. In the latter experiment, cerebellectomy was conducted at least one week before the experiment, and the same 3-hour training was applied (n = 8). From eye velocity data recorded during these experiments, the parameters G1, G0and H of the OKR model implemented on MATLAB Simulink were estimated by fitting the unit step response of the model to the experimental data with a nonlinear optimization method (lsqnonlin function). Results: In the acute cerebellectomy experiment (Fig. 1b, Left), G1and G0increased during the 3-hour training. H increased for the initial 20 min, but decreased thereafter and reached back to its pre-training value. After acute cerebellectomy, the increased G1and G0went back to their pre-training values. By contrast, the parameter H, which increased and then decreased during the training, increased after acute cerebellectomy. In the chronic cerebellectomy experiment (Fig. 1b, Right), G1and G0did not change during the 3-hour training while H increased gradually and did not show significant decrease unlike in the acute cerebellectomy experiment, and reached to a value comparable to that after acute cerebellectomy. Conclusion: Changes in both Direct and Indirect components of OKR eye velocity represented by G1and G0in the model are totally cerebellum dependent. By contrast, change in VSM time constant represented by H consists of cerebellar and non-cerebellar contributions. These results suggest that OKR adaptation, specifically the changes in the VSM contains cerebellar and non-cerebellar adaptable elements.