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 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:

$$dI\_{\rm noise} = -\frac{I\_{\rm noise}}{\pi\_I}dt + \frac{m\_I}{\pi\_I}dt + s\_1\sqrt{\frac{2dt}{\pi\_I}}\xi\_t\tag{4}$$

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 *mI* and standard deviation *sI*. For the simulation, we set τ*I* = 1 ms, *mI* = 25 pA, and *sI* = 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:

$$a = \begin{bmatrix} 0.02 \\ 0.02 + 0.08r\_i \end{bmatrix} \quad b = \begin{bmatrix} 0.2 \\ 0.25 - 0.05r\_i \end{bmatrix}$$

$$c = \begin{bmatrix} -65 + 15r\_i^2 \\ -65 \end{bmatrix} \quad d = \begin{bmatrix} 8 - 6r\_i^2 \\ 2 \end{bmatrix} \tag{5}$$

In Equation (5), the first row is relative to the excitatory, while the second one to the inhibitory neurons. *ri* is a random variable which spans from 0 to 1, and *i* the neuron index. *ri* was added in order to introduce a further variability in the neuron dynamics: for example, a neuron exhibits classic RS behavior if *ri* = 0, and bursting behavior if *ri* = 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 *aij* = 1, a connection between the node *j* to *i* is present, otherwise *aij* = 0 means no connection. All the auto-connections are avoided (*aii* = 0*,* ∀ i). Then, the value 1 of the non-zero *aij* 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 *sij*.

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):

$$P(m) = m^{-\gamma} \tag{6}$$

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.

#### **DATA ANALYSIS AND STATISTICS**

## *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 *Fij* 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 (*Fij* = *Fij*+1 − *Fij*) and we integrated *Fij* in overlapping sliding time windows of 1 s (*Iij* = *j* ≤*n*≤*j* +59 *Fnj*; 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 *Iij*. 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 (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:

$$\text{CC}\_{ab}(\mathbf{r}) = \frac{\sum\_{0 \le t \le T} (a\_t + \mathbf{r} - \mathbf{c} \cdot \mathbf{a} \succ) \cdot (b\_t - \mathbf{c} \cdot \mathbf{b} \succ)}{\sigma\_a \cdot \sigma\_b} \tag{7}$$

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 *i* = 1 *n j* =*i CC*max *ij* and τmax *i* = 1 *n j* =*i* τmax *ij* 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).

#### **RESULTS**

#### **FINITE SIZE NETWORK DYNAMICS**

#### *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

indicator OGB are shown in the panel **(B2)** (objective magnification 10×, field

"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

circuits (Bonifazi et al., 2009).

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 GABAA 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 GABAB 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 3 | Directed functional connectivity before (left) and after (right) lesion.** The number of OUTPUT and INPUT functional connections has been calculated for all the imaged neurons based on the temporal correlation between the firing onsets of the neurons (see section "Materials and Methods"). The ten top ranked cells, i.e., the cells with the largest number of functional OUTPUT (yellow) and INPUT connections (pink), are represented in the top row. For graphic clarity, the connectivity graphs shown in the 2nd and 3rd rows (respectively INPUT and OUTPUT connections) include only the five top ranked cells. The data refers to a homogenous neuronal network where functional hub cells (i.e., neurons with a very large number of functional connections) were not identified. The fluorescent images show the cells loaded with the calcium indicator Fluo4 (see section "Materials and Methods"). The locations of the two lesions (L1 and L2) are marked by the white arrows. The green rectangle highlights the region shown in **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.

positioned in the center of the m-ERC (delimited by dotted line). S, stimulating electrode; LOT, lateral olfactory tract; PC, piriform cortex; l-ERC, lateral entorhinal cortex; DG, dentate gyrus; m-ERC, medial entorhinal cortex. **(B)** Stereomicroscope photograph of the isolated brain positioned in the perfusion chamber. **(C)** Electrical responses to LOT stimulation recorded in the m-ERC. Left, intracellularly recorded voltage traces from a superficial pyramidal cell lying at 200–300 μm from pial surface (black trace), a

pathway (PC and l-ERC) and subsequently invading the hippocampal structure (DG and CA1, dark and light gray spots, respectively). The left margin of the gray area is aligned to the first component of the m-ERC LFP. Right, simplified scheme of the polysynaptic neuronal circuitry within the m-ERC, based on the evoked response pattern and delay analysis of the neuronal response to LOT stimulation. The gray cell represents a putative interneuron mediating a feedback GABAergic inhibition onto a deep pyramidal cell and a feed-forward inhibition onto another interneuron.

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

**ischemic lesioning of the hippocampus. (A)** Local field potentials (LFP) and multi-unit activity (MUA) raw traces from two selected electrodes (19 and 24, experiment #1) recorded in response to an individual paired-pulse stimulus

(ISI 200 ms) delivered to the LOT. The volume-conducted components originating in DG and CA1 are indicated by the dark gray and light gray dots, respectively. **(B)** Peri-stimulus raster plots for the same two representative *(Continued)*

#### **FIGURE 5 | Continued**

electrodes. The corresponding PSTHs are superimposed (bin size = 4 ms). **(C)** Summary plot of mean number of evoked spikes (mean ± S*.*E*.*M*.*) after 1st and 2nd pulse for all four experiments. ∗*p <* 0*.*05, Mann–Whitney *U*-test. **(D)** LFP and MUA raw traces of one selected electrode recorded in response to a paired-pulse stimulus either before (black trace) or after (gray trace) an

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) ischemic lesion of the hippocampus. **(E)** Peri-stimulus raster plot for the same representative electrode, before and after the lesion. The corresponding PSTHs are superimposed (bin size = 4 ms). **(F)** Summary plot of mean number of evoked spikes by a paired pulse stimulus delivered to the LOT (mean ± S*.*E*.*M*.*) either before or after the ischemic lesion of the hippocampus for all analyzed experiments. ∗*p <* 0*.*05, Mann–Whitney *U*-test.

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 *x*axis 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**

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

term of frequency of the neuron *vs*. the stimulation current. **(B)** Outputs of the detection system to be implemented in the closed-loop set-up of the brain prosthesis. *First row*—*(a)*. Raw electrophysiological signal. *Second row—(b)*. The same signal with added Gaussian white noise to reduce Signal 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

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

function of the olfactory-limbic circuit after an ischemic lesion of the

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

lesion are marked in red.

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, 2011). The structural—functional configuration of the finite size circuits can be replicated by an *in silico* neuronal

complexity ("Neural element") will be interfaced to a hardware neuromorphic chip ("Artificial element—Hardware"), implementing both

"Decoding" blocks.

and its biological counterpart is accomplished by the "Coding" and

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

#### **REFERENCES**

illuminate how network structure constrains and drives network dynamics.

## **ACKNOWLEDGMENTS**

The research leading to these results has received funding from the European Union's Seventh Framework Programme (ICT-FET FP7/2007-2013, FET Young Explorers scheme) under grant agreement n. 284772 BRAIN BOW (www.brainbowproject.eu).

The authors wish to thank Dr. Marina Nanni and Claudia Chiabrera from NBT-IIT for technical assistance in cell culture preparation; Dr. Pieter Laurens Baljon for the MATLAB implementation of the SALPA algorithm; Prof. Marco de Curtis for helpful discussion about the *in vitro* guinea pig whole brain preparation and recording; Prof. Sergio Martinoia, Prof. Sylvie Renaud, Prof. Eshel Ben-Jacob, Prof. Ari Barzilai, and Prof. Yael Hanein for their support and for useful discussion.

The authors at TAU thank the help of Nitzan Herzog for calcium imaging set up, Mark Shein-Idelson for network patterning, Inna Brainis for cell culture and Moshe David-Pur for silicon wafer processing. Dr. Paolo Bonifazi is additionally supported by the Joint Italy—Israel Laboratory on Neuroscience.

A. (1996). Spontaneous periodic synchronized bursting during formation of mature patterns of connections in cortical cultures. *Neurosci. Lett.* 206, 109–112.

*Neural Circuits* 6:99. doi: 10.3389/ fncir.2012.00099

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

*Received: 01 November 2012; accepted: 25 February 2013; published online: 14 March 2013.*

*Citation: Bonifazi P, Difato F, Massobrio P, Breschi GL, Pasquale V, Levi T, Goldin M, Bornat Y, Tedesco M, Bisio M, Kanner S, Galron R, Tessadori J, Taverna S and Chiappalone M (2013) In vitro largescale experimental and theoretical studies for the realization of bi-directional brain-prostheses. Front. Neural Circuits 7:40. doi: 10.3389/fncir.2013.00040*

*Copyright © 2013 Bonifazi, Difato, Massobrio, Breschi, Pasquale, Levi, Goldin, Bornat, Tedesco, Bisio, Kanner, Galron, Tessadori, Taverna and Chiappalone. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in other forums, provided the original authors and source are credited and subject to any copyright notices concerning any third-party graphics etc.*

## **APPENDIX**

**FIGURE A1 | Subpopulation of a neuronal network before and after laser induced lesions.** The first panel, starting from the left upper corner, shows a subpopulation of a neuronal network loaded with Fuo4-AM. White arrows, depicted with L1 and L2, indicate the positions of the lesions inflicted to the network. The red arrows indicate the position of the UV laser focus spot. The average power delivered at the sample, during lesion one, is 4 μW, and during Lesion 2, is 5 μW. We delivered 300 UV light pulse for each lesion, at pulse-repetition rate of 100 Hz. At 25 s, after Lesion 1, the L2 position is centered onto the UV focus spot. The last panel shows the same field of view of the first one, after laser inflicted damages. The cells directly affected by the UV laser presented saturated calcium signal. Numbers indicate seconds. The field of view is 150 × 150 μm. Calcium imaging was acquired at 60 Hz.

## Modular neuronal assemblies embodied in a closed-loop environment: toward future integration of brains and machines

#### **Jacopo Tessadori 1 , Marta Bisio1 , Sergio Martinoia1,2 and Michela Chiappalone1\***

1 Department of Neuroscience and Brain Technologies, Istituto Italiano di Tecnologia, Genova, Italy

2 Department of Informatics, Bioengineering, Robotics and System Engineering, University of Genova, Genova, Italy

#### **Edited by:**

Ahmed El Hady, Max Planck Institute for Dynamics and Self Organization, Germany

#### **Reviewed by:**

Guo-Qiang Bi, University of Pittsburgh, USA Wolfgang Stein, Illinois State University, USA

#### **\*Correspondence:**

Michela Chiappalone, NeuroTech Unit, Department of Neuroscience and Brain Technologies, Istituto Italiano di Tecnologia, Via Morego 30, 16163 Genova, Italy. e-mail<EMAIL_ADDRESS>
Behaviors, from simple to most complex, require a two-way interaction with the environment and the contribution of different brain areas depending on the orchestrated activation of neuronal assemblies. In this work we present a new hybrid neuro-robotic architecture based on a neural controller bi-directionally connected to a virtual robot implementing a Braitenberg vehicle aimed at avoiding obstacles. The robot is characterized by proximity sensors and wheels, allowing it to navigate into a circular arena with obstacles of different sizes. As neural controller, we used hippocampal cultures dissociated from embryonic rats and kept alive over Micro Electrode Arrays (MEAs) for 3–8 weeks. The developed software architecture guarantees a bi-directional exchange of information between the natural and the artificial part by means of simple linear coding/decoding schemes. We used two different kinds of experimental preparation: "random" and "modular" populations. In the second case, the confinement was assured by a polydimethylsiloxane (PDMS) mask placed over the surface of the MEA device, thus defining two populations interconnected via specific microchannels. The main results of our study are: (i) neuronal cultures can be successfully interfaced to an artificial agent; (ii) modular networks show a different dynamics with respect to random culture, both in terms of spontaneous and evoked electrophysiological patterns; (iii) the robot performs better if a reinforcement learning paradigm (i.e., a tetanic stimulation delivered to the network following each collision) is activated, regardless of the modularity of the culture; (iv) the robot controlled by the modular network further enhances its capabilities in avoiding obstacles during the short-term plasticity trial. The developed paradigm offers a new framework for studying, in simplified model systems, neuro-artificial bi-directional interfaces for the development of new strategies for brain-machine interaction.

**Keywords: bi-directional, in vitro, hippocampal cultures, confinement, micro electrode array, robot**

## **INTRODUCTION**

Algorithms based on classical models of computation cannot compare with living beings capabilities in terms of dealing with unexpected situations. Different fields of study, such as developmental biology (West-Eberhard, 2003; Gilbert, 2009), embodied cognition (Clark, 1997), evolutionary robotics (Bongard, 2011), seem to indicate as a likely cause for this shortcoming the lack of a developmental phase in traditional silicon-based technology. This process is especially evident in the Central Nervous System (CNS), where morphological changes, both reversible and permanent, occur on a wide range of different time scales. One possible way to deal with this issue is the realization of hybrid systems, where biological components could be exploited for their plastic properties.

In the recent past, several different hybrid model systems have been developed (DeMarse et al., 2001; Martinoia et al., 2004; Mussa-Ivaldi et al., 2010; Warwick et al., 2010; Kudoh et al., 2011), consisting of living neurons coupled to a robotic system.

This solution allows the use of an artificial body whose dynamics can be easily and completely modeled, as opposed to the case of even the simplest animals. Furthermore, the exchange of information in a hybrid system can be limited to the desired level of complexity.

Following this "embodied neurophysiology" approach, we built a closed-loop electrophysiological system by interfacing a virtual mobile robot with a population of neurons, extracted from rat embryos and cultured over Micro Electrode Arrays (MEA; Novellino et al., 2007). The proposed paradigm represents an innovative, simplified, and controllable closed-loop system where it is possible to investigate the dynamic and adaptive properties of a neural population interacting with an external environment by means of an artificial body (i.e., the mobile robot). The main innovations of this experimental setup are: (i) the flexible software architecture at the base of the closed-loop experiments, here described in detail; (ii) the introduction of a modular network design. Starting from the observation of the high degree of modularity in the

brain, different studies point out how such a property is likely to have a profound impact on neural activity (Hubel et al., 1977; Sporns et al., 2000; Derdikman et al., 2003; Kumar et al., 2010; Pan et al., 2010; Boucsein et al., 2011). In this work, we took advantage of the modular structure of the network to obtain a better separation between interacting cell assemblies. A significant improvement to previous works would be the added capability of inducing plastic changes in a controlled fashion. A step in this direction is taken in this setup by the use of a tetanic stimulation to enhance interconnected pathways to improve robot behavior (Jimbo et al., 1999; Chiappalone et al., 2008), following a collision with an obstacle. It is worth pointing out that the final objective of this work is not to achieve the best possible control of the robot: excluding any biological component would, at this stage, easily provide better performance and more reliable results. What is being developed here is groundwork for the integration of electronic systems and neural networks, with the twofold long-term objectives of taking advantage of neural plasticity in more complex control systems and performing closed-loops experiment to gage the computational and learning properties of relatively simple neural preparations.

## **MATERIALS AND METHODS**

The setup developed for experiments of embodied electrophysiology is characterized by several different software, hardware and wetware components (**Figure 1**). The wetware part consists of hippocampal neurons cultured onto a standard 60-electrode MEA. The front-end electronics are constituted by a MEA1060-Inv-BC amplification system (Multichannel Systems, MCS, Reutlingen, Germany) and the computer used is a desktop machine (Dell Precision T5500, 2.66 GHz, 3.43 GB RAM) equipped with a DAQ E NI6255 (National Instruments, Austin, TX, USA) data acquisition board. An *ad hoc* adaptor was realized to interface the DAQ board with the amplification system. The software used for the

management and acquisition is HyBrain2, a specifically developed software based on what is described in a previous work (Mulas et al., 2010): it allows control of all the parameters of the neurorobotics experiments and performs the required data processing, such as the implementation of the coding, decoding and shortterm plasticity schemes. Information is sent to the culture as a series of electrical stimulations through a Stimulus Generator 4002 (Multichannel Systems). Three different robots can be used for the experiments: two physical ones (Khepera II and its successor Khepera III, from K-Team, Zi les Plains-Praz, Switzerland) and a virtual implementation within the HyBrain2 architecture. The relevant elements of the robot are a set of distance sensors and two independently controlled wheels. Both the physical and the virtual ones have a circular arena with obstacles to move in. In all of the experiments, the task the robot is trying to perform is obstacle avoidance. While both physical robots have been tested and are properly working within the setup, in the following, only experiments with the virtual one are reported. The main problems with the physical robot are the fact that it requires actual tracking from an image to compute its position (which is both machinetime consuming and occasionally fails) and the non-idealities of its sensors: among the other, ambient lighting conditions have an impact on the performance of the infrared distance sensor and it has been reputed unwise to add such a factor of unpredictability at this stage of the development.

## **NETWORK MODULE**

#### **Neuronal preparation: random and modular cultures**

Dissociated neuronal cultures were prepared from hippocampi of 18-day-old embryonic rats (pregnant female rats were obtained from Charles River Laboratories). Culture preparation was performed as previously described (Frega et al., 2012). Briefly, the hippocampi of 4–5 embryos were dissected out from the brain and dissociated first by enzymatic digestion in trypsin solution

a micro electrode array; (ii) a computer which hosts the developed software tool (i.e., HyBrain2) which manages the communication between the biological and the artificial part; (iii) the robotic module composed by a robot, either real or virtual, with sensors and actuators navigating into a circular arena with obstacles.

0.125% (30 min at 37˚C) and subsequently by mechanical dissociation with a fine-tipped Pasteur pipette. The resulting tissue was re-suspended in Neurobasal medium supplemented with 2% B-27, 1% Glutamax-I, 1% Pen-Strep solution, and 10% Fetal Bovine Serum (Invitrogen, Carlsbad, CA, USA), at the final concentration of 60 k cells/ml.

Cells were afterward plated onto standard 60-channel MEAs previously coated with poly-d-lysine and laminin to promote cell adhesion (final density around 1200 cells/mm2 ) and maintained with 1 ml of nutrient medium (**Figures 2A,B**). They were then placed in a humidified incubator having an atmosphere of 5% CO2 and 95% air at 37˚C. Half of the medium was changed weekly. Recordings were performed on cultures between 20 and 60 days *in vitro* (DIVs).

Considering the multitude of connections that usually forms in a random culture, a way to better control the network complexity consists of imposing a constraint to the neuronal cells growth along specific pathways (Chang et al., 2001;Boehler et al., 2012). To do this, a dual-compartment chamber with two interconnecting microchannels has been realized in polydimethylsiloxane (PDMS), a biocompatible, inert, and non-toxic polymer often used to this extent (Raichman and Ben-Jacob, 2008; Levy et al., 2012). The realization of the modular structures has been realized by replica molding using specific master with a previously developed technique (Berdondini et al., 2006). The obtained structures have been then placed on MEA substrates, in order to confine the growth of the neuronal cells that will be plated on it, as shown in **Figure 2B**.

#### **Micro electrode arrays**

Micro electrode arrays (Multichannel Systems, MCS, Reutlingen, Germany) consist of 60 TiN/SiN planar round electrodes (30µm diameter; 200µm center-to-center inter-electrode distance, see

**Figure 2A**) arranged in an 8 × 8 square grid excluding corners. In some devices, one recording electrode is replaced by a larger ground electrode. Each electrode provides information on the activity of the neural network in its immediate area. A microwire connects each micro electrode of the MEA to a different channel of a dedicated amplifying system with a gain of 1100. The amplified 60-channel data is then conveyed to the data acquisition card which samples them at 10 kHz per channel and converts them into digital, 12 bit data (**Figures 2C,D**).

#### **HYBRAIN2 SOFTWARE**

The need for real-time access to data led to the adoption of a general-purpose acquisition card (NI6255, National Instruments, Austin, TX, USA) and required the development of a specific software: Hybrain2. The core of the program handles incoming data from the acquisition card and graphically displays them in a panel such as the one shown in **Figure 3A**. Spike detection options can be selected from this panel, such as threshold amplitudes or update times, as well as software blanking of stimulus artifacts. While a rather sophisticated algorithm (i.e., SALPA filtering; Wagenaar and Potter, 2002) for blanking has been included and validated, it has not been used in the described experiments, as it tends to compete for CPU-time with the rest of the system, leading to occasional resource starvation. In its current version, Hybrain2 does not make use of raw data other than for displaying. Instead, incoming data is processed by a spike detection algorithm (Maccione et al., 2009) whose output is a series of time stamps.

As explained later in more detail, both the coding and decoding algorithms for the closed-loop control of the robot are rate-based, therefore spike time stamps are a lossless representation of incoming data. **Figure 3B** shows the panels used for configuration of

**FIGURE 2 | Random and modular neuronal assemblies over micro electrode arrays. (A)** On the left, a random culture grown on a standard MEA device. On the right, the MEA layout is shown: a squared matrix of 59 micro electrodes (the missing one is the reference electrode), in which the inter-electrode distance is 200µm and the micro electrode diameter is 30µm. **(B)** On the left, a confined culture on a MEA substrate. On the right, the

bi-compartmental system realized in PDMS with two interconnection microchannels. Compartments height is 700µm, and width is 1500µm. Microchannels height is 100µm, and width is 50µm. **(C)** Spontaneous electrophysiological activity of a confined culture of hippocampal neurons, registered from all the micro electrodes. **(D)** A typical hippocampal burst waveform recorded from a single channel.

display panel, including options for data visualization, artifact filtering, and spike detection. **(B)** Several panels allow the configuration of coding and decoding algorithms and saving of data during experiments. **(C)** The robot

robot, this can also be used to draw the arena itself. **(D)** The physical robot inside the arena where two obstacles are placed. The dotted red line represents the trajectory of the robot inside the arena.

the parameters of these algorithms, such as selection of recording and stimulation electrodes, pulses amplitudes and lengths, and maximum and minimum allowed wheel speeds for the robot.

A module of the software is dedicated to managing the robot itself: in **Figure 3C**, a sample experiment with a virtual robot is shown. Here, the software is generating the robot environment as well as controlling all the relevant parameters of the robot itself, while, in the case of a physical robot (such as that in **Figure 3D**), the software provides a simple tracking feature on images provided by a webcam positioned over the arena and the required communication with the robot itself. All the data produced during experiments, including electrode readings, time stamps, and robot navigation data can be stored for later analysis both in text and/or binary format, while common parameters configurations can be saved and loaded in order to minimize experiment setup times and human errors.

#### **ROBOTIC MODULE**

The robot, either virtual or physical, is basically a two-wheeled sensor platform: six infrared sensors are mounted on the robot at different angles, providing information about the distance of surrounding objects in different directions,whereas the speed profile of each wheel determine the direction and velocity of the robot itself.

The arena consists of an enclosed space containing several different round obstacles in random positions and the robot.A typical experiment with the virtual robot is shown in **Figure 3C**: the robot is moving in a 400 × 400 pixels circular arena, where dark green pixels represent obstacles or arena walls, whereas light green pixels are free for the robot to move in. The robot (small pink circle) is collecting information about its environment through its six sensors: each black line departing from the robot represents the line of sight of a different sensor; their angles are fixed with respect to the robot heading (in this case, 30˚, 45˚, and 90˚ on both sides of the robot direction), while the length of each line is equal to the distance from the robot center to the closest obstacle in the sensor direction. This distance defines the reading of the sensor: the output is 0 if the robot is in direct contact with an obstacle, 1 if the closest obstacle is at the maximum distance possible (the diameter of the arena, in this case). The three sensor readings on each side are averaged to provide the neuronal network with a single value per side.

In the case shown in **Figure 3C**, the robot is performing an obstacle avoidance task, as can be inferred by the red trajectory. The speed of a wheel is inversely proportional to the average of the sensor readings on the same-side, therefore the robot turns away from close obstacles. The ideal behavior of the robot is that of a Braitenberg vehicle (Braitenberg, 1984) in the case of no loss of information and no significant delays between sensor data collection and motor command execution. Obtaining a behavior as close as possible to this one is the goal of the coding-decoding-short-term plasticity process implemented here.

During experiments, collisions with obstacles or walls are unavoidable: following such an event, the robot moves back to a previous position in its path, at a fixed distance from any obstacle.

#### **INTERFACING THE NETWORK AND THE ROBOTIC MODULE Decoding scheme**

Although many different decoding schemes are possible, so far the only one implemented has been a frequency rate-based algorithm (Adrian, 1928; Rieke et al., 1997; Martinoia et al., 2004). For this scheme, only a feature of the recorded signals is useful: the frequency of spikes at each location. A group of electrodes (i.e., a sub-population of neurons) on the MEA is selected and defined as the "output area" through the procedure described in the Section "Experimental protocol." The number of spikes occurring over that area in 100 ms, non-overlapping windows constitutes the basis for calculating the motor signal for the corresponding wheel. In the current architecture, a linear relation is implemented between wheel speed and motor signal: if no spikes are detected in a time window, the corresponding wheel turns at a set minimum speed, increasing linearly with the number of detected spikes, up to a defined maximum rate. A low-pass filtering effect is added by taking into account previous samples, in order to smooth robot movements.

Dissociated neural networks are especially prone to bursting (Chiappalone et al., 2006) and this pattern of activity has been shown to code different information than just the sum of its spikes (Cozzi et al., 2006). A module for the detection of bursts has been already added to the Hybrain2 software, but its output is not yet part of the control loop of the robot.

For each wheel, the speed is therefore defined as:

$$\boldsymbol{\omega}\_{i} = \begin{cases} \frac{f\_{i,t} + f\_{i,t-1}}{2f\_{i}^{\text{MAX}}} \left( \boldsymbol{\omega}\_{i}^{\text{MAX}} - \boldsymbol{\omega}\_{i}^{\text{min}} \right) + \boldsymbol{\omega}\_{i}^{\text{min}} & \text{for } f\_{i} < f\_{i}^{\text{MAX}}\\ \boldsymbol{\omega}\_{i}^{\text{MAX}} & \text{for } f\_{i} \ge f\_{i}^{\text{MAX}} \end{cases}$$

where subscript *i* denotes wheel side, ω is the wheel speed, and *fi ,t* is the averaged firing rate over all the electrodes corresponding to the *i*-th recording area at time sample *t*. ω MAX, ω min, and *f* MAX are parameters set by the experimenter before the start of the experiment.

#### **Coding scheme**

Likewise, the coding scheme is linear and rate-based: two groups of electrodes are defined as "input areas" and assigned to the sensors on the left and right side of the robot body. The details for area selection are fully explained in the Section "Experimental Protocol." Each sensor provides a reading, normalized to 1 for an object in direct contact with the robot and 0 for an object at the far end of the designed arena (while this behavior is nearly ideal for the virtual robot, it is far from so in the case of the physical robot, as already mentioned in the Section "Materials and Methods." The readings from the sensors on the same-side of the robot are then averaged and coded back to the corresponding sensory area. As mentioned before, the coding is linear and frequency based: a fixed stimulus is delivered at the sensory area at a frequency directly proportional to the averaged, same-side sensors readings. The stimulation rate for each input region is determined as:

$$s\_i = \left(s\_i^{\text{MAX}} - s\_i^{\text{min}}\right)r\_i + s\_i^{\text{min}}$$

where *si* is the stimulation rate of the *i*-th input area and *ri* the normalized average of all the sensor readings on the corresponding side of the robots, whereas *s* MAX *i s* min *i* are user-set parameters fixing the maximum and minimum stimulation rate.

#### **Short-term plasticity protocol**

In order to progress toward the desired behavior, it is necessary to define a learning rule that allows a modification of connectivity between input and output areas by rewarding "good behavior," while discouraging "bad behavior." The effect of tetanic stimulation in these networks was already demonstrated by our group and by others in the past, showing that a 20 Hz stimulation should strengthen the synaptic connections of receiving neurons (Jimbo et al., 1999; Tateno and Jimbo, 1999; Madhavan et al., 2007; Chiappalone et al., 2008; le Feber et al., 2010). In all these papers the effect of the tetanus on the change of firing rate was studied in a time frame comparable to that of our experiments (30 min to 1 h). Additionally, in a previous paper from our group (Chiappalone et al., 2008), we were able to demonstrate that a single tetanic shock to a neuronal network had an immediate effect in terms of increase in the Post Stimulus Time Histograms (PSTH) area (i.e., increase in the number of spikes evoked by a stimulus), a medium-term effect (i.e., few hours after the tetanus delivery), and a long-term effect (i.e., 1 day after the tetanus delivery).

The above observations have been used to define the learning rule in the current implementation of the software: following each robot collision, a 2-s-long, 20 Hz stimulation is delivered to the same-side input area. The rationale for this choice is that collisions are usually caused by poor correlation between stimulation in an input area and detected activity in the corresponding output area, thus making the network responses to stimulation insufficient to steer the robot in the correct direction. Our hypothesis is that tetanic stimulation strengthens all participating connections, thus correcting the problem, as demonstrated in the studies cited above. A tetanic stimulation induces short-term plasticity effects which allow the groups of neurons involved in the obstacle avoidance tasks to fire at a higher frequency, thus inducing the corresponding wheel to increase the angular velocity. Since input-output regions were selected according to connection strength (see Experimental Protocol below), this should increase responses detected from the desired electrodes upon delivery of a stimulus from the input electrodes. This bring to a generalized strengthening of connections in the network and to an improvement in the driving of the robot.

## **ON-LINE PROCESSING OF ELECTROPHYSIOLOGICAL SIGNALS**

#### **Spike detection**

The electrophysiological signals acquired from MEA electrodes must be preprocessed in order to remove the stimulus artifact and to isolate spikes from noise. The spike detection algorithm uses a differential peak-to-peak threshold to follow the variability of the signal and a set of controls are performed in order to make the algorithm as reliable as possible (Maccione et al., 2009). The threshold is proportional to noise SD and is calculated separately for each individual channel (typically as six or seven times SD) before the beginning of the actual experiment (i.e., during phase 1 of the protocol described below).

## **Blanking of stimulus artifact**

Stimulus artifacts are detected when the recorded signal exceeds a defined threshold much higher than the one used for spike detection. The artifact is then suppressed by canceling the first samples in the spike train occurring immediately after it, corresponding to a signal blanking of 4 ms after stimulus delivery.

## **EXPERIMENTAL PROTOCOL**

The typical experimental protocol followed in this work consists of a five-step procedure:

During the first step of the experimental session, spontaneous activity of the network is subject to observation, in order to determine, empirically,which electrodes are the most likely candidate as "input" sites (i.e., sites from which stimulation must be delivered). Typical features to look for in this phase are a sustained mean firing rate (i.e., sufficient number of spikes per second, usually higher than 0.1 spikes/s) and patterns of activity not synchronous with other regions. The best candidates (usually a set of 8–10 sites) are then selected for the second step of the experiment. From each of the candidate "input" channel, in turn, a 500-µs, 1.5V peak-topeak, bipolar square wave is delivered every 5 s, until a total of 40 stimuli per channels have been delivered, while spiking activity is detected from other electrodes.

At the end of this phase, for every stimulation electrode involved, 59 PSTH are generated (Chiappalone et al., 2007): these graphs report the average number of spikes detected from each electrode in the 600 ms following each stimulation and therefore provide information on the strength of the connections in the culture. Through a custom-made script developed in the Matlab environment (The Matworks, Natick, MA, USA), the generated PSTHs are then compared in order to look for areas that present a significant degree of specificity, i.e., where responses are not elicited by stimulation delivered from all the electrodes, but from some of them. In this way, it is possible to define an output (recording) area that will respond mostly to stimulation from the corresponding input area, while remaining silent during stimulation from the opposite input area (cf., see "Input and Output Sites of a Neuronal Population" of the Results).

During steps 3 and 4, the robot is left free to roam the arena with the rules described above, with a tetanic stimulus following each collision with an obstacle delivered during step 4. If the starting hypotheses hold true, this will progressively drive the network toward the desired condition of reliable and specific evoked responses.

Finally, we collect the data on the robot performances. In order to verify the neural-based behavior of the robot, we compared the results obtained (i) in a neuron-controlled experiment (a MEA with living neurons grown on, bi-directionally connected to the robot), (ii) in a open-loop experiment (a MEA with living neurons grown on, but without sensory feedback), and (iii) in an "empty" MEA experiment (a MEA with culturing medium only). In case (ii), the robot performs in a way imposed by the spontaneous firing rate of the neural network, usually in a random pattern, while in the case of the "empty" MEA (iii) the robot basically drives in a straight line (see the Supplementary Videos and Closed-Loop Robot Navigation of the Results).

#### **DATABASE OF EXPERIMENTS, DATA ANALYSIS, AND STATISTICS**

Experiments on a total of *N* = 17 different cultures, ranging from 20 to 60 DIV, have been conducted: 11 of those were random hippocampal cultures, while the other six experiments were conducted on hippocampal cultures, divided into sub-populations by a confinement mask, as described above. Those six cultures were also compared for spontaneous activity evaluation with a subset of six random cultures (age range of the subset: 21–42 DIV).

In order to highlight differences in term of synchronization between the two populations, a cross-correlation algorithm was applied to spike trains, a technique already introduced previously (Frega et al., 2012). Briefly, the cross-correlation function (i.e., cross-correlogram) is defined by the incidence of a spike at electrode *y* after that a spike was fired at electrode *x*. More specifically, given two spike trains (i.e., *x* and *y*) from two electrodes of a MEA, we count the number of spikes in the *y* train within a time frame around the spikes of the *x* train of ±T (in the order of tens of milliseconds), using bins of amplitude ∆τ (usually set at multiple of the sampling frequency). The correct *Cxy*(τ) is obtained by means of a normalization procedure, by dividing each element of the array by the square root of the product between the number of peaks in the *x* and the *y* train. If the obtained *Cxy*(τ) shows a distribution that clearly deviated from flat, electrodes x and y are considered correlated. For each cross-correlogram *Cxy*(τ) we then estimated the coefficient *C*peak. *C*peak represents the value of the cross-correlogram in an area around the maximum detected peak and it is usually evaluated in order to quantify the *correlation level* among two recording channels. The statistical distribution of all *C*peak values was computed for the two experimental groups during spontaneous activity (i.e., random vs. modular cultures). For each robot run, two different parameters have been computed in order to evaluate the performance of the robot, namely the average distance traveled by the robot between hits (measured in pixels) and the average number of hits per second. The virtual robot is implemented so that following a collision against an obstacle, it

is immediately moved to the last location where its center was at least 20 pixels away from any other object. Since the robot radius is 5 pixels, the lower limit for the average distance traveled by the robot during each robot run is that of 15 pixels.

Statistical tests were employed to assess the significant difference among diverse experimental conditions. The normal distribution of experimental data was assessed using the Kolmogorov-Smirnov normality test. According to the distribution of the data, we performed either parametric (e.g., ANOVA, **Figure 7**) or nonparametric (e.g., Mann–Whitney U test, **Figures 4**–**6** and **8**) tests and *p* values < 0.05 were considered significant. Statistical analysis was carried out by using OriginPro (OriginLab Corporation, Northampton, MA, USA).

#### **RESULTS**

#### **NETWORK DYNAMICS: SPONTANEOUS ACTIVITY IN RANDOM AND MODULAR NETWORKS**

Hippocampal cultures grown *in vitro* over MEAs show a spontaneous (i.e., ongoing) activity, similar to that exhibited by *in vivo* systems during their development (Ben-Ari, 2001) or during deep sleep (Corner, 2008). Their electrophysiological behavior is characterized by spontaneous spiking which becomes synchronized with the maturation of the network, giving rise to phenomena called "bursts," network bursts (Pasquale et al., 2010) or network spikes (Eytan and Marom, 2006). These network bursts are the fingerprints of a steady-state in which the network dynamic found a balance between excitation and inhibition (on average 70–80% of neurons are excitatory ones and the remaining 20–30% is constituted by inhibitory interneurons). Such state can be easily pharmacologically disrupted by acting on the glutamatergic as well as on the gabaergic receptors or by adding neuromodulators (Keefer et al., 2001; Eytan et al., 2004; Frega et al., 2012). Another possibility to alter such stereotyped behavior is to introduce modularity (i.e., interconnected populations) instead of having a single uniform and random culture (Raichman and Ben-Jacob, 2008; Shein Idelson et al., 2010; Kanagasabapathi et al., 2012).

**Figure 4** shows the spontaneous activity from a representative random (**Figure 4A**, top) and a modular culture (**Figure 4A**, bottom) during the fourth week of development.While in the random culture the activity is highly synchronized and packed in the form of "network bursts" (van Pelt et al., 2004; Pasquale et al., 2010), in the modular culture we can identify two different temporal patterns of activity with moments of synchronized bursts interleaved with sparse spiking periods. Synchronized network bursts spread to the whole culture also in the modular networks, even if, globally, modular cultures are much less correlated than the random ones (**Figure 4B**).

#### **NETWORK DYNAMICS: EVOKED ACTIVITY IN RANDOM AND MODULAR NETWORKS**

It is possible to electrically modulate the activity of the network by means of electrical stimulation. The typical response of a

same (S) or other (O) compartment with respect to stimulating electrodes. No statistical differences can be noted in a random culture. N = 11 random cultures. **(D)** Box-plot of the latency from the first evoked spikes in the same (S) or other (O) compartment with respect to stimulating electrodes. In a modular network, the latency between the stimulus and the first evoked spike is statistically lower for the electrodes belonging to the same cluster of the stimulating electrodes. N = 6 modular cultures. Box range: percentile 25–75; Box whiskers: percentile 5–95; line: median; square: mean. Mann–Whitney test for not-normal data, significance level = \*p < 0.05.

network can be evaluated through the Post Stimulus Time Histogram (PSTH, cf., see Materials and Methods). In **Figure 5A** the maps of the PSTH obtained as a consequence of the stimulation from site 13 (top) and site 72 (bottom) are reported in a non-confined culture. Typically, the PSTH is characterized by an "early response," lasting 20–40 ms, and by a late response, lasting more than 100–200 ms, usually due to the generation of an evoked burst synchronized over the whole network (Gal et al., 2010). The integral calculated over the PSTH profile represents the average number of evoked spikes at a specific site and it is used for quantifying the strength of the connection between a specific stimulation site and all the recording ones (Chiappalone et al., 2008). This parameter is at the base of the choice of the input-output connections for our neuro-robotic studies (cf., see Input and Output

Bottom. PSTHs obtained by stimulating electrode 28 (black square) in the bottom compartment of the same confined network. Shaded area indicates the top compartment. X -axis: time (0, 400) ms, bin 4 ms; Y -axis: probability of

> Sites of a Neuronal Population). **Figure 5B** reports the maps of the PSTH obtained in a modular network. When stimulation is delivered from site 21 (top compartment, **Figure 5B** top), mainly the electrodes of the top compartment respond to the stimulation. Few activations can be observed also in the bottom compartment, but with a dominant late response and an almost absent early one. In the same network, when stimulation comes from one electrode of the bottom compartment (electrode 28, **Figure 5B** bottom) practically only that compartment responds to the stimulus.

> To further test the actual confinement of the evoked responses, we also analyzed the distribution of the mean latencies (i.e., the distance between the stimulus and the first evoked spike) obtained for each couple of stimulation-recording electrodes (Mainen and Sejnowski, 1995; Tateno and Jimbo, 1999): simply by eye, it is

**FIGURE 6 | Input-output selection. (A)** Map obtained in a representative "random" culture for the selection of the output sites, given two inputs sites (e.g., 26 and 47): red, left recording area; blue, right recording area. **(B)** Schematic representation of the input (yellow and light blue) and respective recording (red and blue) areas for the same experiment reported in A ("random" culture): note that the selected electrodes are quite spread over the entire recording area. **(C)** Map obtained in a representative "confined" culture for the selection of the output sites, given two inputs sites (e.g., 27 and 62): red, left recording area; blue, right recording area. **(D)** Schematic representation of the input (yellow and light blue) and respective recording (red and blue) areas for the same experiment reported in **(B)** ("confined" culture): note that the selected

recording electrodes are close to the stimulating electrode and they follow the structure of the underlying network. **(E)** A box-plot representing the distances from bisector of the selected recording electrodes in the set of random and confined cultures used within this study (N = 11 random and N = 6 modular cultures). The distribution of the distances in the modular case is significantly higher than in the random case. Box range: percentile 25–75; box whiskers: percentile 5–95; line: median; square: mean. Mann–Whitney test for not-normal data, significance level = \*p < 0.05. **(F)** Pie chart representing the percentage of networks in which at least 50% of the recording electrodes were selected in the same compartment of the stimulating electrode. The percentage is higher for the modular networks (N = 11 random and N = 6 modular cultures).

clear that the evoked response is (mostly) limited to the compartment hosting the stimulation electrode. **Figures 5C,D** reports the distribution of the latencies from the electrodes in the same compartment (i.e., top or bottom) of the stimulating electrode (S) compared to those from the electrodes in the other compartment (O). Only in the case of confined networks (**Figure 5D**) the two distributions are statistically different, being the latencies evaluated in the electrodes belonging to the same compartment of the

stimulation significantly lower than those of the electrodes in the other compartment. This proves that dividing the neural network in two sub-populations has indeed an effect on stimulus response.

#### **INPUT AND OUTPUT SITES OF A NEURONAL POPULATION**

The simplest architecture that can be adopted for the proposed task includes two electrodes to deliver coded sensory information, one for each set of sensors. While the same could be said for output sites, the point of interest in this work was the response of the network, therefore a set of 8–10 electrodes is chosen to act as output sites for each wheel.

The main disadvantage in dealing with dissociated cultures instead of experimental models with a preserved neural structure is the lack of predefined architecture. For this reason, before starting an experiment, a procedure has been performed to define the stimulation (sensory input) and recording (motor output) areas of the network. During this procedure (i.e., phase 2 of our experimental protocol, cf., see Materials and Methods), we stimulated the cultures by delivering trains of 40 electrical stimuli (1.5V peak-to-peak, biphasic pulses, 500µs total duration) from 8 to 10 sites in a serial way. Then, the PSTH area (i.e., the number of spikes in the 600 ms following each stimulation) between each pair of stimulation-recording electrodes is computed and the related maps, like the one reported in **Figures 6A,C**, are

produced. The coordinates of each square in that map represent the PSTH areas at a specific recording site relative to stimulation from the two stimulating sites reported on the axis (Stim[26] and Stim[47] in **Figure 6A**, for example). All the possible inputoutput combinations are explored and only the pathways producing "selective" responses are retained. These "selective" pathways are identified by pool of recording sites with respect to a couple of stimulating sites for which the responses measured fall far away from the bisector (i.e., pool of recording site closer to the axis).

Those specific pathways of sensory-motor activations can be then conveniently utilized for driving the robot and for implementing simple reactive behaviors (e.g., obstacle avoidance). **Figures 6B,D** report the selected inputs (i.e., two electrodes, one for the left and one for the right area) and output regions, characterized by eight electrodes each, corresponding with maps **Figures 6B,D**, respectively for two representative cultures (i.e., random and modular).

The presence of a confinement structure tends to generate networks showing a higher degree of functional separation (i.e., selectivity), as well as a physical one, when compared to totally random networks: as can be seen in **Figure 6E**, the average distance from the bisector of the evoked response pair is significantly increased in the case of the modular network. The geometry of the

stimulation-recording pairs is also affected, as they are more likely to be clustered together on the same half of the culture (**Figure 6F**).

#### **CLOSED-LOOP ROBOT NAVIGATION**

All the parameters relevant to the movement of the robot are recorded during the experiment. In **Figures 7A,B**, more than 1000 s of signal recordings are plotted (**Figure 7A** for the left side and **Figure 7B** for the right side). The top panels are showing sensory information, with the blue trace representing the average value of proximity sensors on the left side of the robot and the red one the average value of those on the right. In the second graph, a measure of stimulation is shown, expressed as the mean stimulation rate. The third line of graphs reports the firing rates, measured in spikes per second; wheel speeds (shown in the lower graph, expressed in pixels per seconds), closely follow neural activity.

The results of the behavior described so far can be observed in **Figures 8A–C**, where a virtual arena is shown along with the path drawn by the robot (in red) in a 20-min long robot run, respectively in an "empty" experiment (**Figure 8A**), an open-loop experiment

The values are obtained in N = 5 experiments for the empty and the open-loop case and in the N = 17 experiments reported in the text (light blue = empty MEA; blue = open-loop MEA; cyan = closed-loop MEA). The closed-loop experiments give the best results. Statistical analysis was carried out by using one-way ANOVA (\*p < 0.05) for normal distributions (Kolmogorov–Smirnov test of normality), while for mean comparison both the Tukey and the Bonferroni tests were used.

false positives in the spike detection algorithm on background noise or stimulation artifacts. As can be inferred from the image, though, their total impact is almost null and the robot moves almost precisely in a straight line. **(B)** Reconstruction of a 20-min long robot trajectory in open-loop. During this robot run the control loop has been opened by stopping stimulation to the neural culture. As a result, the robot is, similarly to the

(**Figure 8B**) and finally a closed-loop experiment (**Figure 8C**). While collisions are fairly frequent even in the latter case the behavior of the robot is still much closer to the desired one rather than in an open-loop configuration, or (obviously) in the absence of a biological substrate. As can be observed from the graph in **Figure 8D**, the average path traveled between hits is significantly higher in the case of a close-loop.

#### **IMPACT OF MODULARITY AND TETANIC STIMULATION ON ROBOT NAVIGATION**

Despite the improvement in performance of the closed-loop scenario compared to the control cases, robot collisions against obstacles are still a frequent occurrence in random networks. Observation of PSTHs reveals that random networks show a very high degree of connectivity, with evoked responses showing a strong overlap regardless of the stimulating electrodes position (**Figure 5A**). The introduction of a confinement mask shows a marked separation in the responses obtained from stimulation, as can be observed from **Figure 5B**. This, in turn, leads to a reduction in the amount of "cross talk" between input and output channels, with a consequent increase in the navigation performance of the robot. **Figures 9A,C** compare the improvement in performance between the random network structure and the modular one. Specifically, **Figure 9A** shows the comparison between performances evaluated as the average distance between consecutive collisions in different conditions (without and with tetanic stimulation, respectively on the left and right graphs), while **Figure 9C** displays the same performances evaluated through a different parameter, average number of hits per second. The tetanic stimulation leads to a further improvement in the performance, especially when performed on a network with a modular geometry, as can be observed in **Figures 9B,D**: the first couple of graphs show the increase in performance following the introduction of the tetanic stimulation routine (in a random network, left, and in a modular one, right) evaluated as distance between collisions, while the graphs in **Figure 9D** show the performance obtained in the same experiment as average number of hits per second. Examples of changes in effective connectivity obtained in modular and random networks can be observed in **Figure A1** in the Appendix. Even if quantification will be necessary, preliminary analyses of changes in connectivity show that tetanic stimulation does affect the network response, by strengthening the connections on one side and weakening or not affecting the connections on the opposite side.

While all of the described comparisons yield statistically significant results in the case of the average distance parameter, it is not the case for the average number of collisions: the only condition that causes a large enough change to be significant is the introduction of a tetanic stimulation on a modular network.

#### **DISCUSSION AND FUTURE DEVELOPMENT**

In this paper we successfully interfaced, in a bi-directional way, a network of neurons coming from the hippocampus of embryonic rats with a virtual robot. The robot, which has sensors and wheels,is forced to move in a static arena with obstacles and its task consists in avoiding collisions. Looking at the spontaneous electrophysiological activity of the network, we first select a set of possible "inputs," then we evaluate the evoked response of the entire culture

by delivering patterns of electrical stimulation. This procedure allows us to select the "outputs" of our network. Then, by applying a linear rate-based decoding strategy,we were able to transform the spikefrequency into velocity and the sensory information collected by the robot"eyes"into stimulationfrequencyfor our neurons. The behavior of the robot during the closed-loop experiments resulted significantly better than that in open-loop (i.e., without any sensory feedback) or the "empty" MEA condition, proving that the activity driving the robot is actually neural-based (cf. **Figure 8**). In general, these results prove that an *in vitro* network of biological neurons can control an external agent. While ours is not the first setup to achieve this goal, in our knowledge, no previous work reports an extensive set of experiments like the ones we performed (DeMarse et al., 2001; Martinoia et al., 2004; Novellino et al., 2007; Bakkum et al., 2008; Kudoh et al., 2011), but, rather they focus on a single thesis supported by data obtained from a limited number of analogous preparations. Here, we introduce for the first time statistical comparisons obtained on a sizable number of different preparations with highly different spiking behaviors, such as those observed on random and modular networks. Furthermore, bi-modularity of cultures is introduced here for the first time in the context of closed-loop interfaces and its impact is shown to be relevant for the performance of the embodied agent.

Early experiments on random networks showed the tendency of these cultures to evolve toward a degenerate state where mostly network-wide synchronous activity can be observed. The addition of a confinement mask and the consequent modularity qualitatively changed the behavior of the network, preventing or at least strongly reducing the appearance of synchronized network bursts (cf. **Figure 4**). This change alone was enough to provide a significant increase in the performance of the robot (cf. **Figures 5**, **6**, and **9**). These results lead to two possible investigation lines on the same experimental setup: increasing the modularity of the network might allow more complex behavior to emerge, while chronic stimulation since the day of plating might be used in future experiments to define functionally but not physically distinct sub-populations of neurons within the same culture.

Another point of novelty in our approach has been the systematic use of tetanic stimulation on hippocampal cultures over MEA. Previous approaches aiming at demonstrating plasticity in neuronal assemblies by using stimulation protocols from embedded extracellular electrodes were always applied to cortical cultures (Jimbo et al., 1999; Madhavan et al., 2007; Chiappalone et al., 2008; Stegenga et al., 2010). Here we used hippocampal cells and we proved that tetanic stimulation worked successfully, providing an increase in performance both in random and modular networks (cf. **Figure 9**). A further analysis on data is being conducted to determine whether it is possible to define a clear relationship between spontaneous activity of the network and its impact on the observed changes in connectivity strength, since the patterns of induced change proved to be more complex than expected (see **Figure A1** in the Appendix for a preliminary example of effective connection changes induced by tetanic stimulation). This could allow the design of a more successful learning scheme. The exact biological mechanisms linking performance increase and tetanic stimulation are still unclear and further investigations and targeted experiments are needed. Along this direction, the use of

networks, in the absence or presence of tetanic stimulation (respectively, left and right graph), evaluated as average distance (in pixels) between consecutive hits. **(B)** Comparison of robot performances between different conditions of tetanic stimulation, in random (left graph) and modular networks (right graph), evaluated as average distance between hits. **(C)** Comparison between robot performance in random and modular networks, in the absence or presence of tetanic stimulation (respectively, left and right graph),

(left graph) and modular networks (right graph), evaluated in hits per second. All the values are obtained in the experiments described in text (N = 11 experiments for the random condition, N = 6 for the modular), with a tetanic stimulation session following each standard robot run. Box range: percentile 25–75; box whiskers: percentile 5–95; line: median; square: mean. Statistical analysis was carried out by using Mann–Whitney test for not-normal data, significance level = \*p < 0.05.

pharmacological manipulation could allow to change the state of the network and thus to investigate roles of synaptic transmission and receptors involved in the process of adaptation and learning depending on specific stimulation protocols.

As expected, the final performance of the robot is worse than what was possible to achieve without including biological components in the closed-loop (data not shown): for the task of obstacle avoidance, it would be possible to program the robot so that it can perform the navigation task with no risk of hitting obstacles. However, our neuro-robotic framework proved to be a valid tool for the study of mechanisms of neural coding and the computational and adaptive properties of neuronal assemblies with the final goal to facilitate progress in understanding neural pathologies, designing neural prosthetics, and creating fundamentally different types of artificial or hybrid intelligence.

## **ACKNOWLEDGMENTS**

The authors wish to thanks Dr. Marina Nanni for the technical assistance in cell culture preparation. The authors would like to thank Alessandro Bosca and Luca Berdondini for their help and technical assistance in the realization of the PDMS masks used for part of the experiments reported in the paper.

#### **SUPPLEMENTARY MATERIAL**

The Supplementary Material for this article can be found online at http://www.frontiersin.org/Neural\_Circuits/10.3389/ fncir.2012.00099/abstract

#### **REFERENCES**

*(GECCO '11)*, ed. N. Krasnogor (New York: ACM).

**Video S1 | Video of a closed-loop robot run.** This video of a virtual robot run is running at 40× real speed. The arena is composed of dark green solid obstacles and light green "floor" which the robot can move upon. The magenta circle is the virtual robot itself, the red dots highlight the path followed by the robot center over time, while black circles represent hits against obstacles. While the amount of obstacles hit by the robot shows that control is not perfect, the robot is able to take advantage of sensory information to extricate itself from all the situations encountered in a limited amount of time and hits.

**Video S2 | Video of an "empty MEA" robot run.** This video of a virtual robot run is running at 40× real speed. The arena is composed of dark green solid obstacles and light green "floor" which the robot can move upon. The magenta circle is the virtual robot itself; the red dots highlight the path followed by the robot center over time, while black circles represent hits against obstacles. The starting direction of the robot in this trial is rotated 90˚ clockwise with respect to the other two shown videos. Total lack of biological material on the MEA prevents a closing of the sensory-motor loop. As a consequence, the robot shows a total inability to navigate its environment. The small changes in robot heading are likely false positives in the spike detection algorithm on background noise or stimulation artifacts. As can be inferred from the video, though, their total impact is almost null and the robot moves almost precisely in a straight line.

**Video S3 | Video of an open-loop robot run.** This video of a virtual robot run is running at 40× real speed. The arena is composed of dark green solid obstacles and light green "floor" which the robot can move upon. The magenta circle is the virtual robot itself, the red dots highlight the path followed by the robot center over time, while black circles represent hits against obstacles. During this robot run the control loop has been opened by stopping stimulation to the neural culture. As a result, the robot is, similarly to the previous case, lacking any capability of navigating its environment. Changes in robot direction are, in this case, provoked by the spontaneous activity of the neural network.

neuronal signals. *J. Neurosci. Methods* 177, 241–249.

in dissociated rat cerebral cortex cell cultures on multi-electrode arrays. *Neurosci. Lett.* 361, 86–89.

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

*Received: 06 August 2012; accepted: 18 November 2012; published online: 12 December 2012.*

*Citation: Tessadori J, Bisio M, Martinoia S and Chiappalone M (2012) Modular neuronal assemblies embodied in a closed-loop environment: toward future integration of brains and machines. Front. Neural Circuits 6:99. doi: 10.3389/fncir.2012.00099*

*Copyright © 2012 Tessadori, Bisio, Martinoia and Chiappalone. This is an openaccess article distributed under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in other forums, provided the original authors and source are credited and subject to any copyright notices concerning any third-party graphics etc.*

## **APPENDIX**

**FIGURE A1 | Maps of changes in effective connectivity. (A)** Changes in effective connectivity occurring during a tetanic stimulation experiment. The large dots, in yellow and light blue, represent electrodes used for delivery of both tetanic stimulation and sensory information for the left and right inputs. The smaller dots in blue and red indicate the position of electrodes used for recording from the two "motor" areas. Change in effective connectivity is defined as the difference in the area of PSTHs measured after and before the short-term plasticity experiment, divided by the average of these two values. Variations greater than 20% are represented as lines on the maps, with gray and black lines indicating, respectively, a decrease and an increase in functional connectivity. Only connections involving either stimulating electrode have been represented for clarity, with thicker lines highlighting the connections used in the closed-loop control of the robot (left input-left output and right input-right output areas). This map, in particular, is displaying the change in connectivity observed on a random culture during a 30-min short-term plasticity experiment. In this culture, tetanic stimulation led to a widespread increase of connection strengths involving the electrode represented in yellow, while those involving the one in light blue underwent a mixed change, with about half of them resulting strengthened and half of them weakened. **(B)** Same map as **(A)** obtained from recordings on a modular culture before and after a 30-min short-term plasticity experiment. While tetanic stimulation was delivered to both the yellow and light blue electrodes (respectively in the "lower" and "upper" halves of the culture), only one of the sub-populations was significantly affected, with a diffuse increase in connectivity.

## Toward a self-wired active reconstruction of the hippocampal trisynaptic loop: DG-CA3

## *Gregory J. Brewer 1,2\*†, Michael D. Boehler 1, Stathis Leondopulos 1,3, Liangbin Pan3, Sankaraleengam Alagapan3, Thomas B. DeMarse3 and Bruce C. Wheeler 3*

*1 Department of Medical Microbiology, Immunology and Cell Biology, Southern Illinois University School of Medicine, Springfield, IL, USA*

*2 Department of Neurology, Southern Illinois University School of Medicine, Springfield, IL, USA*

*3 Pruitt Family Department of Biomedical Engineering, University of Florida, Gainesville, FL, USA*

#### *Edited by:*

*Ahmed El Hady, Max Planck Institute for Dynamics and Self Organization, Germany*

#### *Reviewed by:*

*Ulrich Egert, University of Freiburg, Germany Inah Lee, Seoul National University, South Korea*

#### *\*Correspondence:*

*Gregory J. Brewer, Department of Medical Microbiology, Immunology and Cell Biology, Southern Illinois University School of Medicine, PO Box 19626, Springfield, IL 62794-9626, USA e-mail<EMAIL_ADDRESS>
#### *†Present address:*

*Gregory J. Brewer, Natural Sciences II, Department of Bioengineering, University of California, Irvine, USA*

#### **INTRODUCTION**

The mammalian hippocampus crucially encodes the formation of long-term episodic memories and spatial navigation, yet the staged encoding mechanisms remain elusive. While we know molecular details of many types of synapses in the major regions of the hippocampus important to learning and memory at the single neuron level, we don't know if these regions self-wire into the anatomically accurate network or require external electrical or chemical inputs. Further, brain functional studies are saddled with a tradeoff between high-spatial resolution (e.g., MRI, fMRI, EEG, EcoG) vs. high-temporal resolution (e.g., *in vivo* electrode arrays, single cell patch clamp). To bridge this gap, our strategy employs *in vitro* culture in an attempt to recapitulate entire *in vivo* brain regions in culture. Today a wide variety of cell types from various areas of the brain can easily be explanted, cultured *in vitro*, and studied in detail. While *in vitro* technology does provide exquisite temporal and spatial access it too has significant shortcomings. A key hurdle toward reconstructing brain areas *in vitro* has been the difficulty controlling the structural connectivity among cells to begin to recapitulate the *in vivo* architecture. The connections in the hippocampal formation of the brain uniquely propagate forward excitatory communication from one region to the next with the CA3 region distinctive for recurrent collateral excitation. We begin to create a functional tri-synaptic network of the hippocampal formation from the entorhinal cortex (EC) to the dentate gyrus (DG) to the CA3 to the CA1 (Cajal,

The mammalian hippocampus functions to encode and retrieve memories by transiently changing synaptic strengths, yet encoding in individual subregions for transmission between regions remains poorly understood. Toward the goal of better understanding the coding in the trisynaptic pathway from the dentate gyrus (DG) to the CA3 and CA1, we report a novel microfabricated device that divides a micro-electrode array into two compartments of separate hippocampal network subregions connected by axons that grow through 3 × 10 × 400µm tunnels. Gene expression by qPCR demonstrated selective enrichment of separate DG, CA3, and CA1 subregions. Reconnection of DG to CA3 altered burst dynamics associated with marked enrichment of GAD67 in DG and GFAP in CA3. Surprisingly, DG axon spike propagation was preferentially unidirectional to the CA3 region at 0.5 m/s with little reverse transmission. Therefore, select hippocampal subregions intrinsically self-wire in anatomically appropriate patterns and maintain their distinct subregion phenotype without external inputs.

**Keywords: multielectrode array, dentate gyrus, GAD67, burst, GFAP**

1968; Amaral and Lavenex, 2006). Other aspects of hippocampal anatomy not modeled here are the connections from the EC through the perforant path to the CA3 in a feed-forward fashion. In addition, the hippocampus receives modulatory inputs from the amygdala and basal forebrain. Output from the CA1 proceeds through the subiculum and returns to the EC to complete the loop. With smaller numbers of electrodes placed in the rat brain, others have monitored activity from each hippocampal region in behaving animals to describe specific patterns of activity for each region (Rolls and Kesner, 2006; Leutgeb et al., 2007), suggesting staged encoding, but we lack information about the inputs necessary to evoke these patterns and their network relationships.

To achieve these staged connections, we combine microfabrication (MEMS) technology to channel connections between cultured subregions of the hippocampus on a multi-electrode array to simultaneously monitor activity (**Figure 1**). Inspired by Campenot (1977, 1987), the MEMS device creates compartments in which we separately place cells from each major area of the hippocampus (EC, DG, CA3, CA1) connected by microscale tunnels through which axons can pass between the wells to define neuronal communication pathways between each well (Taylor et al., 2005; Dworak and Wheeler, 2009; Pan et al., 2011; Kanagasabapathi et al., 2012). Due to the defined geometry of the hippocampus, the cells can be dissociated from micro-dissected subregions of DG, CA3, CA1, and EC (Mattson et al., 1989; Baranes et al., 1996; Zhao et al., 2001; Lein et al., 2004). Cells

how the tunnels promoted selective growth of axons from one compartment into another.

are loaded into the compartments after the device is placed over an array of extracellular electrodes (microelectrode array or MEA) to measure neural activity in each subnetwork as well as communication between the subnetworks through the tunnels (Morefield et al., 2000; Czarnecki et al., 2012; Downes et al., 2012; Kanagasabapathi et al., 2012; Dranias et al., 2013).

Here we reconstructed paired components of the tri-synaptic pathway, with a focus on the DG to CA3 connection to determine: (1) whether specific subregions of the hippocampus be reproducibly dissected as evidenced by region-restricted gene expression? (2) Will these regions maintain and establish their original identity in a uniform culture environment when removed from external hormonal gradients and input activity? (3) Given that CA3 development precedes DG *in vivo* (Bayer, 1980), is the natural axonal polarity of DG*>*CA3 intrinsically controlled by the neurons or does *in vivo* recapitulation of connectivity require external cues? (4) Do the dynamics of neural activity differentiate between each area *in vitro* and to what extent are they similar to activity patterns seen *in vivo*? We addressed the above issues by quantitative PCR of region-restricted gene expression, by evaluation of distinct spike and burst dynamics in each sub-region compartment and by establishing the polarity of directional communication between sub-regions, whether random or anatomically accurate from the DG to the CA3. Surprisingly, intrinsic capabilities of the DG neurons promote axon extension toward the CA3 neurons, with limited back propagation.

## **RESULTS**

In order to reconstruct subregions of the rat hippocampus, we microscopically dissected these regions from postnatal day 4 rats. At this time, the CA3 and CA1 are well-developed and the dentate granule and hilar region (DG) have nearly completed neurogenesis (Bayer, 1980). Single cell suspensions from each region were plated into separate compartments of a microfabricated PDMS device, positioned over a 60 electrode microarray (**Figure 1**). Cells were plated at physiological density ratios of 100 for DG to 33 for CA3 or 41 for CA1. Between the two compartments was a 400 um barrier perforated by a series of 51 narrow channels that excluded cell somata, but promoted growth of axons along the length of the tunnels (Taylor et al., 2005; Dworak and Wheeler, 2009; Berdichevsky et al., 2010; Pan et al., 2011; Kanagasabapathi et al., 2012; Wang et al., 2012). For quality control of the dissection and to determine whether the microdissected DG, CA1, and CA3 regions maintain their *in vivo* identity in culture, we performed qPCR on the neurons that developed in the compartments for 3 weeks. We selected several genes based on their demonstrated enrichment in specific regions of the adult hippocampus (Lein et al., 2004). When standard, uniform cultures from a single subregion were assessed from glass slips without tunnel devices, the adult animal enriched gene expression was replicated with detection of specific transcripts for each region DG, CA1, and CA3 cultured in the common culture medium (**Figure 2A**). But would the subregion types of expression be maintained across the microtunnel devices? **Figure 2B** shows that the same subregionenriched gene expression is maintained between compartments with the same subregion in each compartment as well as in devices with different subregions on each side (**Figures 2C–E**). These results indicate the fidelity of the dissection and culture process as well as the ability of these hippocampal subregions to maintain their specific identities in the absence of external vascular, hormonal, or electrical instruction.

We examined network spike and burst dynamics with the goal of decoding communication between hippocampal subregions. We recorded from paired compartments of DG and CA3 neurons as a model of this part of the brain anatomy. As controls, we recorded activity from networks comprised of either DG on both sides or CA3 neurons on both sides of tunnel-connected compartments or in random single compartment models. Regardless of configuration, 80% of electrodes were active with an average spike rate around 12 Hz in the NbAct4 medium (data not shown). However, burst dynamics differed between hippocampal subregions. The spike rate outside of bursts for DG increased with anatomically correct tunnel connection to CA3 neurons, while CA3 neurons showed the opposite trend (**Figure 3A**). **Figure 3B** shows that about 60% of spikes occurred within bursts in DG networks, regardless of configuration, while only 45% of spikes occurred in bursts in the CA3 networks apposed in tunnels. The average duration of each burst (**Figure 3C**) showed trends that

**phenocopies selective expression in the adult hippocampus as measured by qPCR.** Specific genes probed for expression after 3 weeks in **(A)** the indicated homogeneous random cultures without tunnel devices. Note enrichment of Trpc6 in the DG cultures, Prkcd in the CA3 cultures, and Nov in the CA1 cultures. **(B)** The same hippocampal subregions plated into each of two compartments of the tunnel device. Note same enrichment profile in the tunnel device as in the random cultures. **(C)** Enriched expression of the DG gene Trpc6 whenever DG neurons are present in heterologous combinations of hippocampal subregions cultured between tunnels, normalized to DG cultured on both sides of the tunnels. **(D)** Enriched expression of the CA3 gene Prkcd when CA3 neurons are present in heterologous combinations of hippocampal sub-regions cultured between tunnels, normalized to CA3 cultured on both sides of the tunnels. **(E)** Enriched expression of the CA1 gene Nov when CA1 neurons are present in heterologous combinations of hippocampal sub-regions cultured between tunnels, normalized to CA1 cultured on both sides of the tunnels. Note the similar expression of each region-specific gene to neurons of that region in combination with heterologous regions, while the other 4 combinations without this region express lower levels of this marker mRNA (*n* = 3 separate cultures).

mirrored the extra-burst spike rate, with DG apposed to CA3 showing longer burst durations while CA3 bursts were shorter. The larger changes in inter-burst interval (**Figure 3D**) showed a co-modulation upon anatomical connection with increased intervals for the DG and CA3 apposed configuration. Intraburst spike rate (**Figure 3E**) and spikes per burst (**Figure 3F**) also differed with configuration.

*In vivo*, granule cells in DG uni-directionally synapse with pyramidal cells in CA3 with no back-propagation of connections from CA3 to DG. In our preparation, cells from DG and CA3 were plated simultaneously in apposing compartments that would permit connectivity in either direction. If axonal polarity of DG → CA3 is intrinsically controlled by the neurons (i.e., self-wire) in the absence of other external cues found *in vivo*, then polarity

corresponding random cultures and (ii) is higher for DG than CA3 apposed across tunnels. **(B)** Percent spikes in bursts are generally higher for any DG culture than any CA3 culture. Log normal distribution statistics apply to **(C–F)**. **(C)** Burst duration was longer for DG than CA3 when they were apposed. For random cultures without tunnel devices, DG burst durations are much lower than CA3 random cultures. **(D)** Inter-burst intervals are lengthened by 50% in DG apposed to CA3 compared to DG self-apposed across tunnels and 300% compared to DG in random networks. Similarly, inter-burst times are longer for CA3 apposed to DG than CA3 apposed to itself or random CA3 cultures. **(E)** Intra-burst spike rates are shortened by 20% in DG apposed to CA3 compared to DG self-apposed across tunnels but longer in the reverse direction. **(F)** Spikes per burst decreased by 14% in DG apposed to CA3 compared to DG self-apposed across tunnels and even less in the reverse direction. In all cases *n* displayed is total degrees of freedom from burst or non-burst segments from 3 min recordings of networks of 4 random DG, 4 random CA3, 8 DG(DG), 8 CA3(CA3), 5 DG(CA3), and 5 CA3(DG). Different letters above bars indicate significant differences (a shared letter indicates a non-significant comparison) by *post-hoc* Tukey multiple-comparison analysis after significant ANOVA, *p <* 0*.*05, normal distribution statistics.

of connectivity from DG-CA3 should be maintained *in vitro* and could also account for the distinct burst dynamics of DG apposed to CA3. Previously, we showed that selective axon polarity from one side to the other, as opposed to bi-directional axon crossing, could be achieved in cortical neurons across tunnels by plating and culture of one side of a device followed after 1 week by plating and culture on the opposite side (Dworak et al., 2010; Pan et al., 2011). Here we tested whether intrinsic DG neuron properties would mimic *in vivo* conditions and preferentially cross the tunnels to innervate the CA3 neurons in the other compartment when both compartments were plated on the same day, with fewer axons crossing from CA3 to DG (reverse direction). This axon polarity could be measured in our devices by the direction of the time delay between spikes detected on the two microelectrodes embedded in the tunnels as shown in **Figure 4Ai** where a 0.48 ms delay from the DG side to the CA3 side over the 0.2 mm distance indicates a forward conduction velocity of 0.42 m/s. In contrast, when CA3 was connected to CA3 (or DG-DG), **Figure 4Aii** shows a 0.40 ms delay in the opposite direction, implying a reverse speed of 0.50 m/s that was evident in about half the spikes in these configurations. Individual tunnels were examined to determine whether axon conduction velocity was directionally polarized or whether evidence of multiple spike heights and shapes suggested several axons in a tunnel. Statistics based on all spike pairs above

**FIGURE 4 | Native polarity established from DG to CA3.** Delay times of spikes traveling in axons in tunnels were determined from the difference in spike times at two tunnel electrodes separated by 200µm. **(Ai)** Example of spike travelling from DG to CA3 with a 480µs delay indicating a velocity of 0.42 m/s. **(ii)** Example of spike propagation from the top to the bottom compartment (arbitrarily designated reverse direction for CA3-CA3). **(B)** Statistical analysis of directional propagation indicates 62% of tunnels spontaneously connect axons with anatomical accuracy from DG-CA3, while homologous regions across tunnels fail to show polarity (Wilcoxin non-parametric test).

a high positive threshold indicated 81% unidirectional axon conduction from DG to CA3 at a mean velocity of 0*.*54 ± 0*.*02 m/s (SD., *n* = 9167 spike pairs). Other positive conduction velocities were DG:DG 0*.*47 ± 0*.*02 (*n* = 2615) and CA3:CA3 0*.*51 ± 0*.*01 (*n* = 1343) m/s, all within the range for the hippocampus *in vivo* (Patolsky et al., 2006).

Closer examination on a tunnel by tunnel basis revealed a more complicated situation. In 3 tunnels on one array, *>*99% of the spikes propagated from DG to CA3. In other tunnels, examination of the waveforms for negative directions from CA3 back to DG indicated multiple roughly simultaneous spikes, likely the spikes from two or more axons, whose sum was detected by the electrode as a shift in the peak within the 2 ms detection window. We never observed a tunnel with only back propagation, suggesting that these events are rare and supporting the conclusion that axons from DG neurons preferentially connect to CA3. For all the DG-CA3 tunnels with measurable spike pairs (*n* = 26 tunnels), more than *>*60% of the spikes propagated from the DG to CA3 direction (**Figure 4B)**. This polarity contrasts with the nearly equal directional distributions of DG-DG and CA3-CA3.

Differences in GABAergic neuron and astroglia content could also affect burst dynamics. Transfection of networks with a Lenti virus carrying a GFP reporter driven by the inhibitory neuron GAD67 promoter was used to evaluate neurons expressing the GABA synthetic enzyme GAD67. **Figure 5** shows higher GAD+ inhibitory neuron density in heterologous sub-region connections in DG compared to CA3 (**Figures 5A,B**) and higher GAD+ inhibitory neurons/nucleus in homologous sub-region connection in DG compared to CA3 (**Figures 5C,D**). Quantitation of GAD67 neurons relative to nuclei (**Figure 5G**) showed a 5-fold increase in GAD+ neuron density in DG over CA3 for heterologous and a 3-fold increase in the homologous configuration. The dissected DG region which includes the hilus contains a higher percentage of GAD67 expressing neurons than CA3 (Harvey and Boksa, 2011) that could contribute to a stronger inhibitory drive to enable a higher percentage of spikes in bursts (**Figure 3B**) and a longer burst duration (**Figure 3C**) for DG than CA3. The same mechanism could also contribute to the longer inter-burst intervals in DG apposed to CA3 (**Figure 3D**) and offer more opportunity for higher extra-burst spike rates. Some of these GAD67 neurons sent axons across the tunnels (not shown), better seen with GAD65 immunoreactivity (**Figure 5I**), as a feed-forward inhibitory component (Cabezas et al., 2012).

By increasing glutamate uptake, higher astroglial density could also affect burst dynamics (Boehler et al., 2007). A nuclear count of more than 2-fold above the plating density for CA3 apposed to DG (**Figure 5H**) suggested proliferation of astroglia. **Figures 5E,F** show that GFAP stain for astroglia are indeed more activated in the CA3 side than the DG side.

## **DISCUSSION**

Here we reconstructed rat hippocampal sub-regions in pairs connected by axon-conducting tunnels to demonstrate intrinsic retention of *in vivo* behavior in the absence of external electrical and hormonal stimuli. The subregions maintained their physiological distinctions based on qPCR expression of

**FIGURE 5 | DG networks contain 3–5× more GABAergic (GAD67-GFP) neurons than CA3 networks while CA3 has more astroglia.** Some GAD+ neurons traverse the tunnels. **(A–D)** Green is GAD67-GFP expression 11 days after infection with Lenti-virus with GAD-67 promoter fused to GFP. Red is pseudocolored for blue bisbenzamide labeled nuclei. **(E,F)** GFAP immunostain in DG or CA3 compartments. **(G)** Nuclei per somata from bisbenzamide stain for DNA. Note red vertical striped CA3 apposed to DG is 50% higher than CA3 by itself. **(H)** Percent GAD67 labeled neurons per nuclei. (*N* = 6 20× fields from each of 2 networks). **(I)** GAD65 immunolabeled axons traverse tunnels.

subregion-enriched genes, distinct spike dynamics, GABAergic neurons, astroglia, and preferential wiring direction.

Our observations on spike and burst dynamics are consistent with a synchronously connected network; neither DG nor CA3 operates as an independent oscillating/bursting center. Isolated DG neurons in the network burst at a higher rate than CA3 neurons, suggestive of regions preprogrammed to drive faster bursting DG onto slower bursting CA3, but operating in constrained fashion with predominantly forward connectivity from DG to CA3, but sufficient recurrent GABA-ergic innervation in DG and astroglia (Boehler et al., 2007) in CA3 to modulate the dynamic behavior. Further, the higher extra burst spike rate, slightly longer burst duration and GABAergic neuron density in DG provide extra drive from DG to CA3, establishing the background readiness of CA3 to be receptive to other inputs for learning or information fusion.

A key hurdle toward reconstructing brain areas *in vitro* has been the difficulty controlling the structural connectivity among cells to reflect, or even begin to adequately recapitulate the *in vivo* architecture in an *in vitro* model. A variety of novel *in vitro* technologies address this difficult problem. These technologies have been targeted toward modifying the surface chemistry to provide guidance cues that promote preferential attachment and growth (Boehler et al., 2012), using microfluidics (Morin et al., 2006), or alternatively, to capitalize on the intrinsic neuronal property to follow topographical features in their environment such as pillars (Dowell-Mesfin et al., 2004), ridges (Curtis and Wilkinson, 1997), or gradients (Hattori et al., 2010). The tunnel approach used in this paper confirms the intrinsic ability of *ex vivo* neurons to reconnect in an *in vivo* order (Czarnecki et al., 2012; Downes et al., 2012; Kanagasabapathi et al., 2012; Dranias et al., 2013).

Our model is a greatly reduced analogue of the *in vivo* circuit due to its anatomical incompleteness as a two-dimensional network. It fails to include modulatory cholinergic, noradrenergic, seratonergic, or dopaminergic inputs. It certainly lacks hormonal fluctuations and efficient removal of waste metabolites. With extracellular electrodes, we only monitor the net effects of thousands of synapses as individual action potentials per neuron.