In vitro large-scale experimental and theoretical studies for the realization of bi-directional brain-prostheses

Brain-machine interfaces (BMI) were born to control “actions from thoughts” in order to recover motor capability of patients with impaired functional connectivity between the central and peripheral nervous system. The final goal of our studies is the development of a new proof-of-concept BMI—a neuromorphic chip for brain repair—to reproduce the functional organization of a damaged part of the central nervous system. To reach this ambitious goal, we implemented a multidisciplinary “bottom-up” approach in which in vitro networks are the paradigm for the development of an in silico model to be incorporated into a neuromorphic device. In this paper we present the overall strategy and focus on the different building blocks of our studies: (i) the experimental characterization and modeling of “finite size networks” which represent the smallest and most general self-organized circuits capable of generating spontaneous collective dynamics; (ii) the induction of lesions in neuronal networks and the whole brain preparation with special attention on the impact on the functional organization of the circuits; (iii) the first production of a neuromorphic chip able to implement a real-time model of neuronal networks. A dynamical characterization of the finite size circuits with single cell resolution is provided. A neural network model based on Izhikevich neurons was able to replicate the experimental observations. Changes in the dynamics of the neuronal circuits induced by optical and ischemic lesions are presented respectively for in vitro neuronal networks and for a whole brain preparation. Finally the implementation of a neuromorphic chip reproducing the network dynamics in quasi-real time (10 ns precision) is presented.


INTRODUCTION
Millions of people worldwide are affected by neurological disorders that disrupt connections between brain and body, causing paralysis, or impair cognitive capabilities. This number is likely to increase in coming years, yet current assistive technology is still limited. Over the last decade Brain-Machine Interfaces (BMIs) and neuro-prostheses (Nicolelis, 2003;Hochberg et al., 2006Hochberg et al., , 2012Nicolelis and Lebedev, 2009) have been the object of extensive research and offer the promise of treatment for such disabilities. These devices could profoundly improve the quality of life for affected individuals, and could have a more widespread impact on society.
Neural interfaces have mainly been devoted to restoring motor function that is lost due to injuries at the level of the spinal cord (Collinger et al., 2013), or to recover sensorial capacities, e.g., artificial retinal or cochlear implants (Chader et al., 2009). However, recent interest has also focused on neural prostheses for restoring cognitive functions. For example, a hippocampal prosthesis for improving memory function in behaving rats was recently presented (Berger et al., 2011, and the same group has also tested a device in primate prefrontal cortex aimed at restoring impaired cognitive functions Opris et al., 2012).
The realization of such prostheses implies that we know how to interact with neuronal cell assemblies, taking into account the intrinsic spontaneous activation of neuronal networks and understanding how to drive them into a desired state in order to produce a specific behavior. The longterm goal of replacing damaged brain areas with artificial devices requires neural network-like prosthetics or models that could be fed with recorded electrophysiological patterns and that could provide a substitute output to recover the desired functions. While ultimately this approach must be tested and applied in vivo, important insights could be gained using in vitro systems of increasing architectural complexity, which can be more easily and thoroughly accessed, monitored, manipulated, and modeled than in vivo systems (at least at present).
The final goal of the studies presented in this paper is to develop a test-bed for the development of a new generation of neuro-prostheses capable of restoring lost communication between neuronal circuits. These studies constitute the object of the European project BRAIN BOW (www.brainbowproject.eu). Healthy and lesioned in vitro neuronal circuits are characterized in parallel to the development of in silico neuronal networks, with the goal of establishing bi-directional communication to mimic or bypass an injured neuronal network. In order to develop an experimental and computational platform for the prototyping of neuro-prostheses, we followed a bottom-up approach using in vitro biological neuronal systems with increasing structural complexity. Our approach takes advantage of the unique features of in vitro neuronal cultures, which represent a powerful experimental model to investigate the inherent properties of neuronal cell assemblies as a complement to artificial computational models. We use engineered networks of increasing structural complexity, from isolated finite-size networks up to interacting assemblies, as a model of intercommunicating neuronal circuitries. Moreover, we scaled our studies up to the isolated whole guinea-pig brain (IWB), to translate to an in vivo model.
In this paper we present the overall multidisciplinary strategy and preliminary results on the different building blocks of the project. The structure-function relationship of "finite size circuits" was characterized with single cell resolution by combining calcium imaging and immunocytochemistry. Similarly to what previously observed in isolated neuronal clusters (Shein-Idelson et al., 2010), we found that the frequency of synchronous network events increased with circuit size. This result was reproduced by in silico neural network models based on Izhikevich neurons with scale-free connectivity. The feasibility of controlled network lesions was explored by optically transecting cell processes and monitoring the subsequent change in functional network connectivity. In addition, in a whole brain preparation, a focal ischemic lesion in the hippocampus was demonstrated to cause an interruption of the limbic olfactory pathway. Finally, a neural network hardware model with arbitrary connectivity based on Izhikevich neurons, working at nanosecond time scale, is presented. These experimental and computational platforms represent a starting point for restoring functional closed-loop communication in a neuronal network with lesioned circuitries.

EXPERIMENTAL MODELS
The repertoire of activity patterns exhibited by an in vitro neural network is strongly dependent on the complexity of its geometry . While homogeneous networks ( Figure 1A) tend to display highly stereotyped bursts which spread to most of the connected cells (Kamioka et al., 1996;Van Pelt et al., 2004;Chiappalone et al., 2006;Eytan and Marom, 2006), networks composed of smaller sub-networks with sparse connections (Figure 1C) usually present non-repetitive patterns of sparse spiking and local bursts (Macis et al., 2007;Shein-Idelson et al., 2010). The first cellular model proposed in this work is that of finite size network (Figure 1B), namely an isolated neuronal circuit consisting of a small number of neurons (dozens to a few hundreds) that is still able to spontaneously produce bursts similar to those observed in larger homogeneous networks (cf. section "Results"). Characterization of activity within these assemblies could allow their use as building blocks for larger, more complex structures of interconnecting sub-networks. At the other end of the complexity spectrum we set the isolated whole brain of a guinea pig ( Figure 1D). This model is used to investigate the properties of one complex functional neuronal assembly (the olfactory tract, see below) embedded in an intact brain (cf. section "Results").

Finite size networks: patterning, cell culture, and calcium imaging
The procedure adopted for the preparation of "finite size networks" is in accordance with the NIH standards for care and use of laboratory animals and was approved by the Tel-Aviv University Animal Care and Use Committee. Cultures were prepared as described in Herzog et al. (2011). After the fourth day in vitro, the growth medium was enriched with 0.5% Pen-Strep (Biological Industries Beit Haemek), 2% B-27 (Gibco), and 0.75% glutamax (Biological Industries Beit Haemek). Cells were plated at a density of 750 cells/mm 2 on a 23 mm square glass coverslip previously glued on a 35 mm petri dish. Coverslips were coated with spots of poly-D-lysine (PDL, Sigma), and petri dishes were homogenously coated with PDL. The cells attaching homogeneously on the free surface of the petri dish (i.e., not covered by the glass coverslip) functioned as a "supporting network" (Kleinfeld et al., 1988). PDL spots were created using either manual drop deposition or polydimethylsiloxane (PDMS) stencils. For manual drop deposition, an Eppendorf pipette with a tip of 10 μl capacity was used. The spots were created by touching the tip filled with 2 μl PDL on the coverslip surface and then drying the coverslips at 37 • C for 30 min.
When PDMS stencils were used, the procedure to create PDL spots was based on a soft lithography process, as described in Sorkin et al. (2006). Briefly, an SU8-2075 (Micro Chem) mould on a silicon wafer with a feature thickness of approximately 200 μm was used to shape the PDMS. The feature was composed of squares of 700 μm × 700 μm separated by at least 1 mm, in order to obtain isolated neuronal islands. The size of the square was chosen to fit the field of view of a 10× objective in the calcium imaging setup described in detail below and in Herzog et al. (2011). Once the PDMS substrate was shaped and dried on the silicon wafer, the PDMS stencils were detached and placed directly on the glass coverslips. Drops of the PDL solution were dripped onto the PDMS stencil until the features were completely covered. After mild vacuum degassing for 15 min, the excess PDL solution was removed and the sample was dried at 37 • C for 30 min. The PDMS stencil was removed before cell plating.
Calcium imaging of the patterned neuronal networks grown on coverslips was performed in buffered-ACSF solution (containing, in mM, 10 HEPES, 4 KCl, 1.5 CaCl 2 , 0.75 MgCl 2 , 139 NaCl, 10 D-glucose, adjusted with sucrose to an osmolarity of 325 mOsm, and with NaOH to a pH of 7.4). In order to load the cells with the calcium-sensitive dye, cultures were incubated for 30 min in 1 ml ACSF supplemented with 1 μl of 10% pluronic acid F-127 (Biotium 59000) and 1 μl Oregon-Green BAPTA-I AM (Invitrogen 06807) previously diluted with 7.6 μl anhydrous-DMSO. Following incubation, cultures were washed with ACSF and recorded at 37 • C. In order to avoid artifacts due to evaporation and pH change, the ACSF was replaced every 20 min during the recording session.
Calcium-fluorescence images were acquired with an EMCCD camera (Andor Ixon-885) mounted on an upright Olympus microscope (BX51WI) using a 10× water-immersion objective (Olympus, NA 0.4). Fluorescent excitation was provided via a 120 W mercury lamp (EXFO X-Cite 120PC) coupled to the microscope optical axis with a dichroic mirror, and equipped with an emission filter matching the dye spectrum (Chroma T495LP). Images were acquired at 59 fps in 2 × 2 binning mode using Andor software data-acquisition card (SOLIS) installed on a personal computer.

Immunocytochemical staining
At the end of calcium-imaging experiments, cultures were washed twice with PBS, then fixed with 4% PFA (15 min) and left in PBS for not more than 5 days before staining. For immunocytochemical staining, fixed cultures were washed three times with PBS (10 min each) and then incubated with 1% Triton ×100 in PBS for 30 min. Cultures were blocked with 2% BSA, 10% normal serum and 0.5% Triton × 100 in PBS for 2 h at room temperature. The cultures were incubated overnight with the first primary antibody (GAD67, 1:250, Millipore, MAB5406) in blocking solution at 4 • C. The next day cultures were incubated with the second primary antibody (MAP2, 1:500, Chemicon, AB5622) overnight at 4 • C. Cultures were then washed three times with TBS and incubated with the secondary antibodies in 2% BSA, 2 mM CaCl 2 in TBS for 1 h at room temperature. After being washed three times with TBS the cultures were mounted with aqueous mounting medium containing DAPI (vector).

In vitro whole brain
Young adult Hartley guinea pigs (150-300 g, Charles River) were used for IWB recordings. All procedures were approved by the Italian Department of Health and were conducted in accordance to FELASA guidelines and Italian and European directives (DL 116/92 and 2010/63/EU). Animals were anesthetized with sodium thiopental (125 mg/kg, i.p.) and transcardially perfused with a cold (4 • C), oxygenated (95% O 2 , 5% CO 2 ) saline solution composed of 126 mM NaCl, 3 mM KCl, 1.2 mM KH 2 PO 4 , 1.3 mM MgSO 4 , 2.4 mM CaCl 2 , 26 mM NaHCO 3 , 15 mM glucose, 2.1 mM HEPES, and 3% dextran (MW 70,000). The pH of the solution was corrected to 7.1 with 1N HCl. After assessing the absence of nociceptive and ocular reflexes, the brain was gently dissected out of the skull, transferred to a recording chamber, and perfused at 7 ml/min with the above solution (pH = 7.3, 15 • C) via a peristaltic pump (Minipulse II, Gilson, France) through a cannula inserted in the basilar artery ( Figure 5). Prior to recording, the temperature of the preparation was gradually increased to 32 • C (0.2 • C/min) (Llinas et al., 1981;Muhlethaler et al., 1993;De Curtis et al., 1998). In order to induce an ischemic insult in the hippocampal formation, a silk thread was positioned under the left rostral and caudal posterior cerebral arteries [r-and c-PCA, see Librizzi et al. (1999)] and a loose knot was prepared around the vessels. The flow was interrupted by pulling the thread ends to tighten the knot (Figure 5) (Pastori et al., 2007).

Optical manipulation and recording system for in vitro neural networks
The optical system combined a laser dissector with a microscope for simultaneous fluorescence and bright field imaging during electrophysiological recording of neural network activity, as previously described (Difato et al., 2011a). The light source used to perform calcium fluorescence imaging was composed of TTL modulable laser diodes (TECBL-15 G-473-TTL-FC, World Star Tech. Inc., USA) coupled to the microscope (BX51, Olympus, Italy) through a circle top-hat engineered diffuser (ED1-C20-MD, Thorlabs, Optoprim, Italy) to remove laser speckles. A pair of UV doublets (Thorlabs, Optoprim, Italy) coupled the laser light to the microscope objective (60×, 0.9 NA water dipping). The laser light was focused on the back focal plane of the microscope objective to produce a homogenous wide field illumination on the sample. A light emitting diode at 590 nm wavelength served as the bright field illumination source (M590L2, Thorlabs, Optoprim, Italy). The wavelength of the diode was chosen to avoid interference with the emission spectra of the fluorochrome (Fluo4-AM, Invitrogen) used to label the sample. A dichroic mirror separated the light coming from the sample (green and red portion of light spectra) onto two cameras. Green emission light was deviated on CCD1 (V887ECSUVB EMCCD, Andor, Lot Oriel, Italy) acquiring the calcium fluctuations due to network activity, and the red portion of the light spectra was deviated on CCD2 (Pilot PIA1000-48GM, Basler, Advanced Technologies, Italy) to perform bright-field imaging. The CCDs image acquisitions and light sources were synchronized with a TTL signal coming from a D/A board (PCI-6529, National Instruments, Italy). The use of TTL-modulable light sources for fluorescence and bright field imaging allowed a precisely timed illumination of the sample, thereby reducing phototoxicity and facilitating long term calcium imaging of neural networks. Bright-field images were acquired at 1 Hz to detect network topography before and after laser dissection of network connections. Cells were previously incubated for 10 min with 5 μm Fluo-4 AM (Invitrogen, Italy). To monitor the neural network activity before and after laser induced network lesions, calcium imaging was performed at 60 Hz (light exposure of 3 ms each frame, at an average power at the sample of 60 μW).
Cells were kept under the microscope at 35 • C using a Peltier device (QE1 resistive heating with TC-344B dual channel heater controller, Warner Instruments, Italy). For neuronal cultures plated on Petri dishes, pH and humidity were controlled by aerating a custom-designed polydimethylsiloxane (PDMS) sleeve, which integrated the objective for optical access, with humidified carbogen (95% O 2 , 5% CO 2 ).
A pulsed, sub-nanosecond UV Nd:YAG laser at 355 nm (PowerChip nano-Pulse UV laser PNV-001525-040, Teem Photonics, Italy) served as the source for performing laser microdissection experiments. The diaphragm of the epi-illuminator was substituted by a narrow-band laser mirror, which reflects 355 nm laser light while passing all other wavelengths coming from the laser diodes used for fluorescence microscopy (DM6, TLM1-350-45-P, CVI, Italy), thus allowing fluorescence imaging and laser dissection to be performed simultaneously. Damage to neural network was inflicted with laser pulse repetition rate settled at 100 Hz, and an average power at the sample of about 4 μW.

Electrophysiological system for the in vitro whole brain
Extra-and intracellular recordings were performed simultaneously in piriform and medial entorhinal cortex (PC and m-ERC). To test the viability of the preparation throughout the experiment, we monitored evoked local field potentials (LFPs) in PC and m-ERC in response to the electrical stimulation (0.5-3 mA, 0.3 ms) of the lateral olfactory tract (LOT) using custom-made bipolar electrodes made of twisted, insulated silver wires. Intracellular recordings were performed with sharp micropipettes filled with 3M potassium acetate (input resistance 70-110 M ) and attached to an electronically controlled micromanipulator (Sutter Instruments, Novato, CA, USA). Signals were amplified by an intracellular amplifier (IR-283A Cygnus Technology, PA, USA). Field potentials were recorded using glass pipette filled with 0.9% NaCl (resistance 2-5 m ) or microwire arrays (Tucker-Davis Technologies, Alachua, FL, USA) featuring 16 tungsten planar recording wires (filament diameter 50 μm, tip angle 45 • ), each separated by 250 μm (impedance 30-40 K ). The extracellular signals were acquired using a PBX3 preamplifier (Plexon, Dallas, TX, USA) configured to separately process spikes (150 Hz-8 KHz bandwidth) and local field potentials (0.7-300 Hz).
Data were digitized at 25 kHz using a PCI-6071E A/D board (National Instruments, Austin, TX, USA) and stored on the hard drive of a personal computer. Recordings were performed using ELPHO software developed by Dr. Vadym Gnatkovsky at the C. Besta Neurological Institute (Milan, Italy).

COMPUTATIONAL MODEL
In the following sections we will present the computational model used to mimic the dynamics expressed by finite size networks (cf. section "Experimental Models").

Neuron model
The neuron model used for the finite size networks is based on the Izhikevich equations (Izhikevich, 2003). The dynamics of this model depend on four parameters that, correctly chosen, reproduce the spiking behavior and voltage traces of specific types of cortical neurons. From a mathematical point of view, the model is described by a two-dimensional system of ordinary differential equations.
with the after-spike resetting conditions: In Equations (1-3), v is the membrane potential of the neuron, u is a membrane recovery variable which takes into account the activation of K + and inactivation of Na + channels; I syn describes the synaptic input from other neurons; I noise is a current source Frontiers in Neural Circuits www.frontiersin.org generator introduced to model the spontaneous subthreshold electrophysiological activity of the neurons. Practically, we introduced a stochastic source of noise (modeled according to an Ornstein-Uhlenbeck process) to each neuron described as follows: In Equation (4) the quantity ξ t is a white noise with zero mean and unitary variance. In this way, I noise is Gauss-distributed at any time t and, after a transient of magnitude τ I (correlation length), converges to a process with a mean equals to m I and standard deviation s I . For the simulation, we set τ I = 1 ms, m I = 25 pA, and s I = 9 pA. Among the possible firing patterns generated by the neuron model of Equations (1, 2), we implemented the family of regular spiking (RS) and the family of fast spiking neurons (FS) in percentage of 75% and 25%, respectively, in agreement with the experimental findings (cf. section "Finite Size Network Dynamics"). Mathematically, the four aforementioned parameters were set as follows: In Equation (5), the first row is relative to the excitatory, while the second one to the inhibitory neurons. r i is a random variable which spans from 0 to 1, and i the neuron index. r i was added in order to introduce a further variability in the neuron dynamics: for example, a neuron exhibits classic RS behavior if r i = 0, and bursting behavior if r i = 1.

Finite size network model
Graph theory was used to represent the network connectivity. All graphs are defined by nodes which represent the neurons, and edges which model the morphological connections among the neurons. The structure of the graph is described by the adjacency matrix, a square matrix of size equal to the number of nodes N with binary entries. If the element a ij = 1, a connection between the node j to i is present, otherwise a ij = 0 means no connection. All the auto-connections are avoided (a ii = 0, ∀ i). Then, the value 1 of the non-zero a ij elements has been replaced to mimic different synaptic strengths. Synaptic weights were chosen randomly from a normal distribution with a mean value and standard deviation equal to 10 and 3.5, respectively. To model the synaptic transmission we chose the approach of the pulse-coupled neural networks: practically, the firing of the j-th neuron causes an instantaneous change in the membrane potential of the neuron i-th by means of the weight s ij .
Among the possible graphs, following the experimental findings regarding the functional connectivity of such confined neuronal assemblies (cf. section "Finite Size Network Dynamics"), we implemented neuronal networks with a scale-free (SF) connectivity (Barabasi and Albert, 1999). Briefly, in SF networks the degree distribution follows a power law: if m is the number of edges incident to a node, i.e., the degree, the power law distribution is given by Dorogovtsev and Mendes (2002): where γ is the characteristic exponent. This law suggests that most nodes have just a few connections and other, named hubs, have a very high number of links. To build a SF network, we made use of the algorithm described in Batagelj and Brandes (2005), particularly efficient in terms of computation when dealing with large-scale networks. Nodes are added successively. For each node, m edges are generated. The endpoints are selected from the nodes whose edges have already been created, with a bias toward high degree nodes.
In order to mimic the experimental conditions of the confined assemblies described in section "Finite Size Network Dynamics," in section "Simulation Results" we presented the results regarding the ongoing activity of networks made up of 90, 100, 120, 150, 240, 320, and 520 neurons.

Analysis of network dynamics based on calcium fluorescence imaging
Custom software running in MATLAB (Crépel et al., 2007;Bonifazi et al., 2009) was used for the automatic identification of the cells loaded with the calcium indicator and for the extraction of their fluorescence signals as a function of time (time resolution 59 Hz). To detect the calcium events (i.e., the onset and offset of neuronal firing) from the fluorescent trace F ij of the neurons (1 ≤ i ≤ M, M number of neurons; 1 ≤ j ≤ N, N number of frames) we calculated the first derivative of the fluorescent signal ( F ij = F ij+1 − F ij ) and we integrated F ij in overlapping sliding time windows of 1 s (I ij = j ≤n≤j +59 F nj ; 1 ≤ j ≤ N − 59). A Gaussian fit centered at zero was used to extract the standard deviation σ i of the noise of the processed signal I ij . Signal transients exceeding the threshold of 3σ i for at least 5 consecutive points were considered as calcium events. The onset and the offset of these calcium events were determined using a four-parameter sigmoidal equation as described in Takano et al. (2012). The estimated onset and offset times were fixed respectively to the 5% and the 95% of the sigmoidal plateau.
The reconstruction of the functional connectivity of the network was based on pair-wise correlation analysis of the onset time series extracted from the calcium imaging data, as described in Bonifazi et al. (2009). Briefly, when the firing onset of cell j preceded in a repetitive way the firing onset of cell k, a functional connection directed from j to k was established. In order to reveal these temporal correlations, the post-stimulus time histogram of cell k centered on the firing onsets of cell j was calculated within a maximal time lag of 500 ms. Both the Student's t-test and the Kolmogorov-Smirnov test with a level of confidence of 5% were used to exclude the possibility that the poststimulus time distribution could be a Gaussian distribution with zero mean or a uniform distribution, respectively. In this way, we excluded cases where the activation of two neurons was completely uncorrelated Frontiers in Neural Circuits www.frontiersin.org March 2013 | Volume 7 | Article 40 | 5 (uniform distribution) or synchronous (Gaussian centered at zero). The cross correlation between firing onsets time series of individual neurons was used to estimate the average correlation and average time of activation of each neuron relative to all others, similarly to what described in Bonifazi et al. (2009) and Marissal et al. (2012). Briefly, the cross correlation between the time series of neurons a and b was calculated as follows: where σ a and σ b are the standard deviation of the time series, t is the sampling time, T the duration of the entire movie and |τ| ≤ 1 s. The maximum cross-correlation value (CC max ab ) and the time lag of its occurrence (τ max ab ) were used to calculate, respectively, the average correlation and average time of activation of neuron i to the following formulas CC max where n is the number of neurons displaying a positive cross-correlation with neuron I.

Processing of electrophysiological signals from the IWB
Raw data acquired by the ELPHO software were loaded into MATLAB (Mathworks Inc., Natick, MA, USA) for off-line processing. First, raw traces were band-pass filtered to select either multi-unit activity (MUA, 800 Hz-3 KHz) or local field potentials (LFP, 1-300 Hz). Stimulation artifacts were suppressed using an off-line MATLAB implementation of the SALPA algorithm (Wagenaar and Potter, 2002). Highly noisy channels were visually excluded from the analysis. Then, MUA raw data were spikedetected by means of the PTSD algorithm (Maccione et al., 2009) (peak lifetime period = 2 ms; refractory period = 1 ms; threshold = ±8 times the estimated noise standard deviation). The result of the spike-detection procedure consists of a series of point processes (i.e., spike trains), one for each recording channel (Bologna et al., 2010).
We evaluated the network-wide evoked response by computing the Peri Stimulus Time Histogram (PSTH; Perkel et al., 1967) for each recording channel of the array and for the full array [time bin = 4 ms, time window = (−100 ms, +400 ms) relative to the stimulus onset]. We also measured the intensity of the response as the average number of evoked spikes in a 200-ms time window following each stimulus. The final dataset comprised 4 recordings in control brains (duration ∼300 s, 10-20 paired pulses delivered to the LOT at 0.05 Hz, inter-pulse interval 200 ms) and 3 recordings before and after the induction of focal ischemia (same stimulation protocol).

Spontaneous synchronizations in finite size networks
To build an experimental model for the study of physiological and impaired communications between neuronal assemblies we grew finite size neuronal networks, i.e., networks composed of neuronal assemblies spatially separated by hundreds of micrometers and interconnected through long neuritis. As a first step, we focused on the properties of single modules, i.e., the structural and dynamical properties of isolated and spatially confined neuronal circuits (Figure 2). Isolated neuronal circuits located within an 800 × 800 μ m spot were obtained by plating the cells on glass cover slips previously coated with a geometrically defined molecular adhesive layer (PDL). The individual cell populations varied between a few dozen up to a few hundred neurons. Similar to homogenous and clustered cultures (Chiappalone et al., 2006;Shein-Idelson et al., 2010), finite size circuits displayed spontaneous synchronized events after 2 weeks in culture (Figure 2, panel B1) occurring with a frequency linearly correlated with the number of cells present in the circuit (Pearson correlation 0.88, Figure 2C1). Likewise, depending on the density of the plating and on the vicinity to the supporting network, finite size circuits organized into monolayers or in three-dimensional clusters, with a higher propensity of clustering at increased plating density or at larger distances from the supporting network (data not shown). We used calcium imaging of monolayer neuronal circuits (performed with a 10× objective) in combination with immunocytochemical staining to map the functional and structural properties of all the neurons in the circuits with singlecell resolution. GABAergic cells could be specifically identified ( Figure 2A3), allowing us to investigate their specific involvement in spontaneous synchronization processes, similar to the work of Bonifazi et al. (2009) in developing hippocampal networks.
A pair-wise analysis based on the cross-correlation between the firing onsets time series of pairs of neurons (see section "Materials and Methods") was used to estimate the average correlation and average time of activation of each neuron relative to all others (Bonifazi et al., 2009;Marissal et al., 2012). In all the circuits analyzed (n = 4) the time correlation graph presented a bimodal distribution (Figure 2C2), indicating that network events synchronized first the population of neurons plotted on the left side of the graph (i.e., with a time lag < 0), whereas neurons on the right (i.e., with a time lag > 0) were activated next. In addition, the presence of highly correlated early activated GABAergic neurons was observed (red points within the violet circle in Figure 2C2). Interestingly, the existence of a characteristic, earlyactivated neuronal population within the network synchronizations has been already documented in developing hippocampal circuits (Bonifazi et al., 2009) even in absence of GABAergic transmission (Marissal et al., 2012). Notably, in the presence of GABAergic transmission it has been shown that early-activated GABAergic neurons can play the role of hub cells in orchestrating network dynamics (Bonifazi et al., 2009). The similarity between these previous observations and the results presented here suggest that cortical circuits share common innate features in their functional organization.

Effect of laser ablation on functional connectivity
To monitor the synaptic re-organization of lesioned neuronal circuits with single cell resolution, we reconstructed the functional connectivity of a neuronal subset of a larger neuronal network 20 min before and after laser-induced ablations (see section "Materials and Methods").
Two micro-lesions (lesion 1 and lesion 2) were induced next to the center of the field of view, using an average laser power at the sample of 4 μW and 5 μW, respectively. The second lesion was performed at higher power to obtain a more pronounced alteration of the network. Indeed, this lesion produced a strong intracellular calcium increase in several cells, and a calcium "shockwave" started to propagate through the network. After a few minutes, only directly ablated cells displayed a saturated calcium fluorescence signal, while the other neurons recovered a relatively low basal calcium level and presented spontaneous activity (cf. Figure A1). The frequency of occurrence of spontaneous network synchronizations was not affected by the lesions (Figure 3, 4th and 5th rows) with no significantly statistical difference between the inter-burst interval distribution before and after lesion (student t-test, p > 0.05). However, the number of cells recruited within the network events in the imaged field (i.e., close to the location of the lesion) decreased by 31 ± 10% (student t-test, p < 0.05). Based on the calcium dynamics of the cells imaged in a circular field of 244 μm diameter (Figure 3), we reconstructed the functional connectivity of the neuronal population through a pair-wise analysis of the onset of firing (see section "Materials and Methods"). Briefly, if the activation of cell i reliably preceded the activation of cell j (i.e., over several repetitions with statistical significance, see section "Materials and Methods"), we inferred a functional connection directed from i to j. Cell pairs that were synchronously activated or not displaying any activation order were not included in the directed functional connectivity reconstruction (see section "Materials and Methods"). Figure 3 (1st row) shows the location of ten neurons with the highest number of functional INPUT (violet) and OUTPUT (yellow) connections before and after the lesions. Interestingly, after the lesions, top rank INPUT and OUTPUT neurons segregated into spatially distinct regions. Top rank OUTPUT neurons relocated in the bottom right region while top rank INPUT neurons remained in the rest of the circuit. In addition, just one out of the ten neurons for each group belonged to the top rank group before and after the lesion. The relocation of the functional connections (drawn for clarity just for the five best ranked neurons) can additionally be observed in Figure 3 (2nd and 3rd row).

In vitro WHOLE BRAIN
We also characterized the activity of an ex vivo experimental model (i.e., the isolated brain of a guinea pig, Figure 4) before and after a lesion induced by a focal ischemia.

Network response to LOT stimulation in the m-ERC
Electrical stimulation of the LOT induced a polysynaptic response in the m-ERC mediated by the interposed activation of the hippocampus (Biella and De Curtis, 2000;Gnatkovsky and De Curtis, 2006) (Figure 4). The intracellular correlate of the LOTevoked delayed response in neurons of m-ERC superficial layers was characterized by an early GABA A receptor-mediated inhibitory postsynaptic potential (IPSP; latency from LOT stimulation: 51 ± 1 ms, n = 12), followed by a relatively slow (duration 409 ± 36 ms) NMDA-dependent depolarizing component which often reached threshold for spike firing. Conversely, pyramidal cells in deeper layers responded to LOT stimulation with an early excitatory postsynaptic potential (EPSP) occurring 15 ± 1 ms after the population spike recorded in the dentate gyrus (DG, Figure 4). The EPSP often crossed the threshold for action potential firing and was followed by a relatively slow inhibitory potential mediated by GABA B receptors (Gnatkovsky and De Curtis, 2006). The early inhibition of the superficial principal cells is presumably due to a feed-forward mechanism sustained by interneurons recorded in layers II/III (i.e., basket and chandelier cells; Canto et al., 2008). In Figure 4 the firing of an interneuron  Figure A1. The field of view is a circular region of 244 μm diameter. The raster plot (representing the firing onsets) and the fraction of activated cells are shown respectively in the 4th and 5th row.

FIGURE 4 | The guinea pig isolated whole brain (IWB). (A)
Schematic view of IWB observed from its ventral surface. The circle of Willis with its principal branching arteries is highlighted in black. The whole brain is perfused by means of a peristaltic pump that delivers ACSF to the brain through a polyethylene cannula inserted into the basilar artery. The two vessels that are occluded to induce the hippocampal ischemia are marked by red crosses. In the same hemisphere a microelectrode array ( corresponds to the early IPSP measured in the pyramidal cells in the same layer. Spiking responses to paired-pulse LOT stimulation (interpulse interval 200 ms) were recorded by 16-channel MEAs implanted in the superficial layers of the m-ERC (200-500 μm from pial surface; Figure 5A). Figure 5B shows the peri-stimulus raster plots of two selected channels (19 and 24, experiment #1) in response to each of the two LOT stimulations for a selected experiment. An earlier phase, which we observed in almost all active recording channels, was characterized by two relatively sharp peaks: the first corresponding to the far-field response originating in the hippocampus and the other corresponding to the initial phase of the m-ERC response (Figure 4). This was followed by a late, long-lasting but less reliable component (cf. channels 19 and 24). The histogram in Figure 5C displays the number of spikes (mean ± S.E.M.) evoked by the 1st and the 2nd pulse for all experiments (control condition) as a measure of response intensity. In 2 out of 4 experiments (#1, #4) we observed a stronger activation after the 1st rather than 2nd pulse, whereas in the other 2 experiments (#2, #3) responses to the 2nd pulse were slightly stronger than to the 1st pulse (no significant statistical difference). However, one must consider that first evoked responses in experiments #1 and #4 were on average more intense, probably reflecting a relatively high probability of excitatory neurotransmitter release upon the first pulse. This would nearly deplete the available pool of synaptic glutamatergic vesicles, leading to a paired-pulse depression of the postsynaptic response.

Cutting the olfactory pathway: hippocampal focal ischemia
Occlusion of the posterior left cerebral arteries abruptly reduced ACSF perfusion of the hippocampus, resulting in a block of the propagation of the synaptic activity toward the entorhinal cortex (Figure 4). About 5 min after the ischemic insult, LOT stimulation failed to evoke any response (Figure 5). Stimulus-triggered raster plots and the corresponding pre-and post-lesion PSTH are shown in Figure 5E. The bar graph in Figure 5F summarizes the total number of spikes evoked by a paired-pulse stimulus before and after the ischemic lesion. A significant reduction of the response intensity caused by the lesion was observed in all three analyzed experiments.

SIMULATION RESULTS
In this section, we report the results of simulations in which we modeled the effects changing the number of neurons in confined networks. Each simulation lasted 10 min, sampled at 10 kHz. Networks were simulated in MATLAB (The Mathworks, Natik, US). Peak trains were stored and then processed by using SpyCode software (Bologna et al., 2010), conveniently adapted to manage large-scale networks.

Dynamics of finite size networks
We simulated the ongoing activity of neuronal networks made up a 90, 100, 120, 150, 240, 320, and 520 neurons. The choice of these networks sizes followed from the experimental findings described in section "Finite Size Network Dynamics" (assuming a neuron/glia ratio equal to 2:1). In addition, 25% of such neurons were considered inhibitory (Isaacson and Scanziani, 2011) and were modeled as FS neurons (cf. section "Computational Model"). Model neurons were connected following a scale-free (SF) topology. Figure 6A shows the degree distribution of the simulated SF networks. For all SF networks, the degree distribution was fitted by a power law and the corresponding exponent lay between −1.04 (network made up of 90 neurons) and −1.34 (networks made up of 520 neurons).
The simulated networks displayed spontaneous synchronized events (network bursts) independently of their size ( Figure 6B). However, the frequency of occurrence of those synchronized events varied in a linear manner with respect to the number of cells present in the circuit (Pearson correlation 0.96, Figure 6C). To facilitate comparison with Figure 2C1, the xaxis of Figure 6C reports the total cell number (neurons + glia), although the number of neurons effectively simulated is indicated near the blue dots. The results of the simulation were fit well with the experimental data, as confirmed by the slope of the linear fit (0.00015 vs. 0.00016). An interesting finding was that the simulated networks tended to show a higher proportion of random spiking activity and less bursting than normally observed in actual finite-size neuronal networks. This is consistent with other experimental results of interconnected finite-size networks previously reported in the literature (Macis et al., 2007).

HARDWARE SET-UP FOR A BRAIN PROSTHESIS
The hardware set-up that will be used to interface the biological component (either the neuronal culture or the in vitro whole brain) is a Spiking Neural Network (SNN) system. This SNN implements biologically realistic neural network models, spanning from the electrophysiological properties of one single neuron up to network plasticity rules. As already discussed in the modeling section, the choice of Izhikevich neuron model is relevant because (1) it is biologically realistic, and (2) it operates in biological real time. By real-time, we mean that computation results are provided within a firmly controlled delay (10 ns precision), which is lower than the sampling period (100 μs to 1 ms). Among these modules, the computation-critical task is the implementation of a SNN model, which represents the prosthesis itself, and the analysis of biological signals to produce events from the recorded activity.
The digital Izhikevich neurons and detection system are implemented as a configurable digital integrated circuit (fieldprogrammable gate array, FPGA) using the VHDL language. We implement Regular Spiking (RS) neurons (excitatory) and Fast Spiking (FS) interneurons (inhibitory) similar to those found in cell culture (Figures 2, 3 and 6). The hardware models follow the Izhikevich equations with parameters corresponding to RS activity (a = 0.02, b = 0.2, c = −65, and d = 8). In Figure 7A1 Frontiers in Neural Circuits www.frontiersin.org we describe the choice of the topology (Cassidy and Andreou, 2008) to implement the Izhikevich equations. We implement a neuron on FPGA board Xilinx Virtex 5 XC5VLX50. This neuron uses really few resources (only 2% of the FPGA) and works in real-time. In Figure 7A2 we compare the behavior f(I) of biological RS neurons and one RS neuron implemented into the FPGA.
Concerning the SNN, our goal was to implement a model using 80 neurons (FS and RS) with high connectivity capacity (e.g., 6400 synapses). Network structure is fully configurable, and synapses are excitatory or inhibitory conductances which provide current depending on the postsynaptic membrane voltage. Delays are also implemented to provide good accuracy on timing. The network is defined into the RAM of the digital board where lie all characteristics of all neurons and synaptic connections in the network. A synaptic connection is defined by a synaptic weight and the address of the neuron linked by this synapse. Added with complementary functions like loopback stimulation and monitoring, this system will be able to perform cross-platform neural computation.
The detection of neural electrophysiological activity is done by a reconfigurable acquisition based on wavelet detection circuit for in vitro biological signals. Our strategy for real-time spike detection is to implement a pre-processor, which emphasizes spikes shapes and attenuates out-of-band noise. This pre-processor provides two outputs corresponding to different wavelet detail levels. The first one is essentially composed of out-of-band noise used to determine a threshold level adapted to the signal amplitude. The second output is compared to the threshold to discriminate spike events. The pre-processing algorithm is the Stationary Wavelet Transform (SWT). The detection system computes in real-time the SWT, the adaptive threshold and the Frontiers in Neural Circuits www.frontiersin.org  Second row-(b). The same signal with added Gaussian white noise to reduce Signal to Noise Ratio. This step was added to stress the capability of the system to detect action potentials in difficult conditions. Third row-(c). Output of the stationary wavelet decomposition preprocessing module. We used a Haar mother wavelet with 16 bits fixed point computation. The output signal is the sixth level detail output of the decomposition tree. Fourth row-(d). Binary output of the detection module. This output is the result of a threshold applied to the signal in (c). The threshold is computed from the standard deviation of the first level detail output of the wavelet decomposition tree. The emphasized detected spike is a false positive. This shows that the signal (b) represents the limit of signals that can be reliably processed by our system. These signals were first recorded then input to the system with a waveform generator.
comparison. This method proved to be very efficient to extract action potential of excitable cells from very noisy signals (Raoux et al., 2012). Figure 7B shows the performance of the method on a single channel setup. Action potentials are emphasized by arrows on the signal A. We added significant noise [signal (b)] and then sent the signal to the detector that provides outputs (c) and (d).
To summarize, all modules (i.e., Izhikevich neuron, neural network specifications, detection and stimulation modules) will be implemented into the FPGA. This modular system will be used as a cross-platform neural computation unit. Microelectrode arrays will be used to record and electrically stimulate living neural networks, with a specific emphasis on stimulation localization. Dedicated integrated electronics will be designed to implement the communication channels between the living and the artificial networks. The biological signals (from living to artificial) will be processed by using on-line spike detection algorithms and a rate-based decoding (Rieke et al., 1997;Novellino et al., 2007;Tessadori et al., 2012), while the firing rate of an artificial neuronal sub-network will be translated into the stimulation frequency for the biological network (from artificial to living), thus following a similar rate-based strategy. The system including the artificial and living neural networks will form a closed loop with a regulated feedback (cf. next Section).

A BI-DIRECTIONAL NEURO-PROSTHESIS
The knowledge that we gained through the various studies presented here will contribute to the final realization of a bidirectional communication between in vitro and in silico models of interconnected cell assemblies. By studying the dynamics of in vitro networks (see Figures 2, 3), we will create a computational model (see Figure 6) exhibiting the same I/O function of its biological counterpart (Figure 8, panel A). Through this approach we also plan to further our knowledge about the  ). Second, we will use more complex modular networks and replace one sub-network module with our computational model of the finite-size network, in order to replicate the function of the intact system (A, right). Finally, the same conceptual approach will be adopted to recover the function of the olfactory-limbic circuit after an ischemic lesion of the hippocampus (B). The bidirectional interaction with a model reproducing the function of the damaged area will allow restoring the original I/O pattern. s(t): stimulus function; y live (t): response function of a healthy preparation; y sim (t): response function of the neuronal network model; y damaged (t): response function after lesion in the IWB; y hybrid (t): response function of the hybrid system resulting from the combination of biological and artificial components. In panel (B), the hippocampal areas targeted by the ischemic lesion are marked in red.
interplay between structural connectivity and dynamics in neuronal networks. Once we have realized and tested our model, we will bi-directionally integrate it into a biological network made up of few interconnected sub-networks in replacement of one of these that has been previously lesioned ( Figure 8A, right panel). The same conceptual approach will be applied to the olfactorylimbic pathway in the IWB (Figure 8, panel B). After a thorough characterization of spontaneous activity patterns (e.g., spontaneous periodic events, which strikingly resemble the ones shown by primary cortical cultures; see Figure 1) and LOT-evoked responses generated in the m-ERC (see Figure 8), we will include such information into a realistic computational model. We will then induce an ischemic lesion of the hippocampus and realize a functional model able to reproduce the same transfer function of the damaged part in order to restore the original pathway. Figure 8 summarizes this approach, both for in vitro interconnected finite-size networks and for the guinea pig IWB. The final step foreseen in the BRAIN BOW project is the hardware implementation of the signal processing algorithms and computational models to achieve our proof-of-concept neuroprosthesis based on a neuromorphic chip. Figure 9 illustrates the closed-loop architecture that we plan to develop. Raw traces recorded by means of either planar or implanted MEAs (depending on the biological sample) will be fed into the artificial element and pre-processed online to extract multi-unit activity patterns (MUA). Spatio-temporal spiking patterns will then be translated by the "decoding" block into signals delivered to the neuronal network model. After elaboration, output patterns produced by the model will be finally translated by the "coding" block into a stimulation delivered to the neural element (Figure 9).

DISCUSSION
This paper presents a bottom-up, multidisciplinary approach toward the realization of a neural prosthesis capable of replacing lesioned neuronal circuitries. The final goal of the studies consists of developing a neuromorphic chip reproducing the function of a lesioned circuit without replicating its specific architecture or structural organization. As a general model of a self-organized neuronal circuit, finite size neuronal circuits in culture are produced and studied in an isolated configuration to reveal innate (and therefore most general) features of intra-circuit organization (cf. Figure 1). Since finite size networks can spontaneously interconnect in a "multimodular" network organization, they also represent an optimal experimental model to reveal innate inter-circuit communication properties (cf. Figure 1), as shown in previous studies (Macis et al., 2007;Raichman and Ben-Jacob, 2008;Shein-Idelson et al., 2010. The structural-functional configuration of the finite size circuits can be replicated by an in silico neuronal network and then implemented on a neuromorphic prosthetic chip. The capability of the neuromorphic chip to replace the function of a lesioned circuit will be tested at increasing levels of network complexity from an in vitro modular network to an isolated whole brain system (IWB). In the attempt to present the overall scientific approach of the BRAIN BOW project (cf. section "Introduction" and "A Bi-Directional Neuro-Prosthesis"), this paper shows first results from the different level of investigation grounding the overall strategy.

FINITE SIZE CIRCUITS AND INNATE FUNCTIONAL ORGANIZATION OF CORTICAL CIRCUITS
As previously shown by Shein-Idelson et al. (2010), cultured cortical neuronal networks composed of at least a few dozen neurons are able to produce spontaneous collective dynamics known as network bursts, characterized by oscillatory activity in the gamma-theta range, and with the frequency of the bursts increasing with the number of neurons in the network. We confirmed these findings here using optical measurements on monolayer circuits (cf. Figure 2). By combining calcium imaging with immunocytochemistry, we have found that network events first recruit a characteristic population of neurons which includes GABAergic neurons. In particular, the time-lag correlation of the finite size cortical circuits is similar to what observed in developing hippocampal circuits (Bonifazi et al., 2009), in which a scale-free functional connectivity distribution was accompanied by the existence of GABAergic hub cells able to play a key role in the orchestration of the spontaneous network events. All together, these observations suggest that cortical neuronal circuits share a common innate functional organization which might include the existence of GABAergic hub cells.

MONITORING EFFECTS OF LESIONED NEURONAL CIRCUITS IN FINITE SIZE NETWORKS
After characterizing the spontaneous dynamics of the finite size networks we monitored how a focalized lesion can trigger functional reorganization in the neuronal circuit. We made controlled laser ablations of different intensities on our networks (e.g., targeting single modules, inter-connections between modules, single neuritis/cell bodies/cell assembly). After the lesions, the neuronal circuits continued to produce spontaneous networks events with no significant changes in the frequency of occurrence (Figures 3 and A1). These were presumably generated out of the imaged field where the lesions were performed. The number of cells recruited during network events decreased either because they were directly lesioned by the laser ablation or because of a change in the local functional organization of the circuits (see the functional connectivity graphs of Figure 4). In a previous study by Maeda et al. (1995) the authors made a lesion in a homogeneous network over a MEA to study the origin of spontaneous network bursting. More recently, Difato et al. (2011a) reported controlled sequential ablation of single connections in a neuronal network, causing modulation of its activity without irreversibly damaging it. By combining MEA recording and calcium imaging the authors found changes in electrophysiological patterns in the network and identified the contribution of neuronal sub-populations to the network activity (Difato et al., 2011b). To the best of our knowledge, our study is the first to make a spatially defined micro-lesion at the single cell scale and to analyze the neuronal dynamics and connectivity by means of optical-only tools. This methodology, which can be extended to the use of genetically encoded calcium sensors, allows a more detailed and prolonged monitoring of the functional reorganization of the circuit over hours or days with the advantage, when compared to electrophysiological recordings, that the high spatial resolution (i.e., single cell) can be linked to morphological/structural cellular properties through post-hoc immunocytochemical characterizations. This could also facilitate testing of methods to promote functional circuit repair, such as pharmacological approaches.

SIMULATION RESULTS (SOFTWARE AND HARDWARE)
Given the similarity between the synchronization dynamics observed in developing hippocampal networks (Bonifazi et al., 2009) and in the finite circuits (Figure 2) with early activated GABAergic cells forecasting synchrony, we hypothesized a common innate structural-functional organization in neocortical and paleocortical circuits. Therefore, we used a scale-free topology (Barabasi and Albert, 1999) to model a neuronal network based on Izhikevich neurons (Izhikevich, 2003) (Figure 6). The proposed model was able to reproduce the empirical dependence between bursting rate and circuit size. However, the model predicted a richer repertoire of firing patterns (e.g., Figure 6B). Indeed, such patterns can be found in biological networks (Segev et al., 2002;Macis et al., 2007;Marconi et al., 2012). Thus our synthetic models (conveniently tuned and adapted) are able to reproduce the dynamics found in in vitro networks. Our results also demonstrate that the hardware element of the prosthesis (cf. section "Hardware Set-up for a Brain Prosthesis" and "A Bi-Directional Neuro-Prosthesis") can be constituted by a neuromorphic model (SNN) built on the same equations as the computational model (Izhikevich, 2003), since it reproduces similar firing rate distributions (Figure 7). Thus, the computational (software) model serves as a bridge between the biological networks and the hardware implementation.

COMPARISON TO PREVIOUS WORK AND PROSPECTIVE RESULTS
In the last decades, great efforts have been made to develop neuroprostheses to restore lost sensory or motor functions (Taylor et al., 2002;Chader et al., 2009;Collinger et al., 2013), but very few groups have focused on neuro-prostheses targeting lesions at the level of the CNS and aimed at recovering lost cognitive capabilities (Berger et al., 2011;Prueckl et al., 2011;Bamford et al., 2012;Hampson et al., 2012;Opris et al., 2012). Although our studies are limited to simplified in vitro models of cell assemblies, their final aim is to provide useful insights for the design of future cognitive prostheses. We believe that our approach would help us understand how we can influence/drive the dynamics of a neuronal assembly by interfacing it to an artificial network, implemented either in software or hardware. This is not the first attempt to realize an in vitro closed-loop system: previous studies have used a robotic actuator or a control algorithm aimed at clamping network activity to a desired level (Demarse et al., 2001;Martinoia et al., 2004;Wagenaar et al., 2005;Wallach et al., 2011). However, we seek to extend these approaches by replacing a real biological network with a simulated network, and hence by implementing bi-directional communication between biological and simulated networks. This research project builds on previously published results in the field of in vitro closed-loop electrophysiology (Arsiero et al., 2007). It can also be generalized to a more structured experimental model like the in vitro whole brain of a guinea pig, which lies between in vivo (as it retains the original tridimensional architecture) and in vitro (as it is disconnected from any sensory input/motor output). In contrast to other groups which have exclusively investigated in vivo brain prostheses (Prueckl et al., 2011;Bamford et al., 2012;Berger et al., 2012;Hampson et al., 2012;Opris et al., 2012), we are trying to exploit the unique advantages of in vitro electrophysiology-accessibility, visibility and control of physical and chemical conditions-to study neural information processing in neuronal assemblies, and to understand which parameters are relevant for effectively interfacing biological and artificial networks. In addition to informing the design of future in vivo approaches, our approach could also illuminate how network structure constrains and drives network dynamics.