Large-Scale Mapping of Axonal Arbors Using High-Density Microelectrode Arrays

Abstract

Understanding the role of axons in neuronal information processing is a fundamental task in neuroscience. Over the last years, sophisticated patch-clamp investigations have provided unexpected and exciting data on axonal phenomena and functioning, but there is still a need for methods to investigate full axonal arbors at sufficient throughput. Here, we present a new method for the simultaneous mapping of the axonal arbors of a large number of individual neurons, which relies on their extracellular signals that have been recorded with high-density microelectrode arrays (HD-MEAs). The segmentation of axons was performed based on the local correlation of extracellular signals. Comparison of the results with both, ground truth and receiver operator characteristics, shows that the new segmentation method outperforms previously used methods. Using a standard HD-MEA, we mapped the axonal arbors of 68 neurons in less than six hours. The fully automated method can be extended to new generations of HD-MEAs with larger data output and is estimated to provide data of axonal arbors of thousands of neurons within recording sessions of a few hours.

Introduction

The classical view of axons is that of mere transmission cables (Hodgkin and Huxley, 1952), while dendrites integrate distinct synaptic inputs (Spruston, 2008), and learning and memory are perceived to be a consequence of synaptic plasticity (Redondo and Morris, 2011). Recent data, however, indicate that the functional capacity of axons may be much more complex (Debanne, 2004; Ohura and Kamiya, 2016; Rama et al., 2018).

According to the classical theory describing giant squid axons (Hodgkin and Huxley, 1952), voltage-gated sodium and potassium channels support self-sustained action potentials that propagate along the axon. That is also true for mammalian axons, but at least two additional types of cationic channels have been described. These channels are activated by G-protein dependent receptors or hyperpolarization and modify the shape and propagation of the action potentials (Elgueta et al., 2015; Ko et al., 2016).

Axons can be surrounded by a myelin sheet that leads to saltatory conduction of action potentials (Tasaki, 1939), instead of the continuous conduction that has been observed in non-myelinated axons, such as the giant squid axons (Hodgkin and Huxley, 1952). Mammalian axons also show considerable variations in their diameter (Gaussian distribution around a peak at 200 nm diameter) and their number of branch points and varicosities, which affect the conduction velocity of action potentials. Differences in axon lengths have been found to produce defined temporal delays for spike trains for coincidence detection in the auditory system (Seidl et al., 2010; Stange-Marten et al., 2017). Furthermore, axonal delays in the neocortex show sub-millisecond precision (Swadlow, 1994), which is compatible with mechanisms of spike-timing-dependent synaptic plasticity (Dan and Poo, 2004). These characteristics lead, according to computational studies (Izhikevich, 2006), to the emergence of precisely timed firing patterns.

Early computational studies have suggested that some axonal arbors can act as filters for spike patterns by the selective failure of action potential propagation (Lüscher and Shiner, 1990). Whereas cerebellar mossy fibers in the cerebellum can reliably transmit spike trains up to 1.6 kHz (Ritzau-Jost et al., 2014), occasional failures have been observed in the axons of CA3 pyramidal neuron at firing frequencies of 30–40 Hz (Meeks and Mennerick, 2007). Failure depends on the diameter and the branching morphology as well as on the activation of presynaptic A-type potassium channels. The morphology and molecular composition of the individual axonal arbors has been found to influence the timing and shape of the action potential (Bischofberger et al., 2002; Alle and Geiger, 2006; Cho et al., 2017) arriving at the pre-synapse, which affects synaptic transmission and plasticity. For example, broadening of action potentials due to slow inactivation of voltage-gated potassium channels during high-frequency spike trains was found to facilitate synaptic release (Geiger and Jonas, 2000). Furthermore, axonal signaling and synaptic connectivity have been found to constitute important parameters in induced-pluripotent-stem-cell (iPSC) models of Parkinson's (Kouroupi et al., 2017) and amyotrophic lateral sclerosis (Wainger et al., 2014).

The understanding of information processing in axons is a fundamental question in neuroscience; however, the availability of experimental data is severely limited due to the small axon diameter. Most data sets have been gathered by performing patch clamp recordings of axonal membranes, at axon terminals and boutons, at the intact axon shaft or at the axonal bleb that is formed upon cutting axons (Ohura and Kamiya, 2016). However, these patch-based methods are not very well-suited to track the propagation of action potentials at more than two sites or even across the full axonal arbor. Moreover, the patch recording time is limited to a few hours while axonal recording is a serial process and has to be done axon after axon. Finally, patch-based methods constitute endpoint measurements and following the development of axonal arbors over extended periods is not possible. Another possibility to study axons includes the use of imaging tools and optogenetics (Chen et al., 2013), and it has been shown that imaging of single axon terminals is possible (Hoppa et al., 2014). However these methods are still limited in signal-to-noise ratio as well as in temporal resolution (Emmenegger et al., 2019).

A viable alternative to the methods described above includes the use of high-density microelectrode arrays (HD-MEAs), based on complementary-metal-oxide-semiconductor (CMOS) technology (Eversmann et al., 2003; Berdondini et al., 2009; Frey et al., 2010; Ballini et al., 2014; Bertotti et al., 2014; Viswam et al., 2016; Tsai et al., 2017). HD-MEAs can be used to capture neuronal activity at high temporal resolution across spatial scales, including networks, dendrites and most importantly, axons (Obien et al., 2015). CMOS-based HD-MEAs have been used to study the action potential propagation along axonal arbors (Bakkum et al., 2013) and the initiation of action potentials at the axon initial segment (Bakkum et al., 2019), as well as to track single action potentials along axons (Radivojevic et al., 2017) and stimulate single axon initial segments (Ronchi et al., 2019). Furthermore, HD-MEA recordings can be combined with classical patch clamp recording to study postsynaptic currents in response to pre-synaptic stimulation (Jäckel et al., 2017).

The objective of this work was to provide a method for high-throughput scanning of axonal arbors and mapping their axonal delays with minimal or no need to adjust parameters for the detection of axonal signals. Using a standard HD-MEA (Frey et al., 2010), we mapped the axonal arbors of more than 68 neurons in <6 h.

Materials and Methods

Animal Use

Timed pregnant rats (Wistar) were obtained from a commercial vendor (Nihon SLC, Japan). Animals were sacrificed on the day of arrival to obtain embryos for primary neuron cultures. All experimental procedures on animals were carried out in accordance with the European Council Directive of 22 of September 2010 (2010/63/EU) and have been approved by the local authorities in Japan (Animal Care and Use Committee of RIKEN; QAH24-01).

High-Density Microelectrode Array (HD-MEA)

Sub-cellular resolution extracellular recordings were obtained using a CMOS-based HD-MEA (Frey et al., 2010) with 11,011 electrodes, arranged in a hexagonal pattern and featuring an electrode density of 3,150 electrodes/mm². The culture chamber of the HD-MEA was prepared as described before (Heer et al., 2007) with minor modifications: after attaching the chamber ring (polycarbonate, 19 mm inner diameter, 8 mm high) using epoxy resin (EPO-TEK 301-2, Epoxy Technology Inc.), GlobTop (G8345D-37, Namics Inc.) was used to cover the bond wires while keeping the electrode area clean, and the remaining area was covered by a thin film of PDMS (Sylgard 184, Dow Corning). Platinum black was electrochemically deposited (Marrese, 1987) [with modifications of the original procedure (Heer et al., 2007)] on the electrodes to decrease their impedance in order to improve signal-to-noise characteristics (Viswam et al., 2014). Before plating, the surface of the HD-MEAs was rendered hydrophilic by oxygen-plasma treatment (40 s, 20 W), incubated for 4 h with of 50 μg/ml Poly-D-Lysine (Sigma-Aldrich, P7280) in PBS, washed twice with aqua dest and air-dried for 1 h.

Primary Neuron Cultures

Adult rats were anesthetized with isofluorane and killed using a guillotine. The embryos were removed from the uterus and decapitated. Their neocortex was dissected in ice-cold dissection medium (HBSS without Ca²⁺ and Mg²⁺; Gibco, NO.14175) and incubated for 20 min at 37°C in Trypsin/EDTA (Sigma-Aldrich). After washing twice with plating medium (Neurobasal A, supplemented with 10% Fetal bovine serum, 2% B27 Supplement, 1:100 GlutaMax, all from Gibco, Japan, and 10 μg/ml Gentamicin, Sigma-Aldrich, Japan), the tissue was mechanically dissociated, passed through a 40 μm nylon mesh, and centrifuged 6 min at 200 g. The supernatant was removed; the cells were suspended and counted. A 20 μl drop containing 10,000 cells was placed on the electrode area of the HD-MEA in the middle of the culture chamber. The cultures were covered with a membrane, permeable to gas but not to water vapor, and placed in a standard incubator (37°C, 5% CO2, 80% relative humidity). The neurons were allowed to settle and attach to the surface during 30 min. Thereafter, the culture chamber was filled with 600 μl serum-free, astrocyte-conditioned DMEM/Hams's F12 medium (Nerve Culture Medium, Sumitomo, Japan, #MB-X9501). Medium was completely exchanged with 600 μl conditioned medium after 4 days and then every 7 days until day 17 in vitro.

Recordings and Spike Event Detection

For recording, HD-MEAs were placed in a bench-top incubator (TOKAI HIT, Japan, INU-OTOR-RE) with temperature control, and 5% CO2 was supplied by a gas-mixer and humidified by a water bath. In order to avoid the evaporation of medium during prolonged recording intervals, the water bath and the lid temperature set point was set 1K and 3K above the sample temperature set point, which was 35°C. The HD-MEA recordings were performed using custom scripts, written in LabView (National Instruments, US), Matlab (Mathworks, US), C++ and Python running on a standard PC with a Linux operating system. Data underwent lossless data compression and were directly stored on a server on the local LAN.

Offline analysis of the recordings included filtering, event detection and averaging. First, a band-pass filter (2nd order Butterworth filter, 100–3,500 Hz) was used to remove slowly changing field potentials as well as high frequency noise. The remaining (background) noise was characterized by the median absolute deviation (MAD), which is resilient to outliers in the data but represents a consistent estimator of the standard deviation, sV = 1.4826 MAD(Vsig). Using a voltage-threshold method for event detection (Lewicki, 1998), negative signal peaks below a threshold of Vthr = 5sV (Vthr > 50 μV in all cases) were identified (Figure 1C). To avoid multiple detection of the same spike, successive events within <0.5 ms were discarded.

FIGURE 1

Figure 1. Outline of the main experimental methods. Neurons from embryonic-day-16 rat cortex were cultured on high-density microelectrode arrays (HD-MEAs) (A). The HD-MEA has a high number and density of electrodes in the recording area and a limited number of low-noise on-chip amplifiers placed at the periphery. A switch matrix allows for using different recording configurations by selecting recording electrodes and connecting them to the amplifiers (B). Extracellular recordings enabled the detection of spiking events by a simple thresholding procedure (C). The amplitude of the negative peak of these spikes was mapped across the entire HD-MEA by recording from all electrodes using non-overlapping block configurations (D). This map reveals the location of individual neurons, as neuronal axon initial segments (AISs) are strong contributors to a neuron's local extracellular field potential. Neuronal footprints (see, e.g., Figure 2B) were obtained by recording with multiple configurations consisting of both, fixed electrodes, as well as randomly selected electrodes (E). The fixed electrodes (indicated in black) were located near the AIS of the neuron and did not change between the configurations. The remaining electrodes of the HDMEA where sequentially recorded in batches (indicated by the same color) consisting of electrodes at randomly selected locations. Spike sorting of the spike events on the fixed electrodes was performed to obtain the activity of single neurons. Spike-triggered averages from all configurations were then combined into a single footprint for each neuron revealing the axonal delays (E). (B) was adapted with permission from Obien et al. (2015).

Mapping of Axonal Initial Segments

To initially identify the location of axonal initial segments (AISs), the whole array was scanned using configurations in which non-overlapping blocks of 6 × 17 electrodes (Figure 1D) were connected to the amplifiers through the switch matrix (Frey et al., 2010) (schematic in Figure 1B). After recording for 30 s from every electrode in the array, spike detection was performed for each electrode and taking the median value summarized the amplitude of the negative peaks. The median of the negative peak for each electrode was plotted as a map showing some areas with large negative peaks (Figure 2A). Such local minima were assumed to indicate the putative (proximal) AIS locations (Bakkum et al., 2019). We did not distinguish between inhibitory and excitatory neurons according to the AIS waveform (Mita et al., 2019).

FIGURE 2

Figure 2. Activity map and footprint of an example neuron. The marker in the activity map reveals the location of the AIS of the example neuron (A). The circle size indicates the square-root-scaled count of spiking events per electrode. The median negative amplitude of the spikes is color-coded with a cut-off at −200 [mu]V. Spike triggered averaging shows the axonal footprint (B). The circle diameter indicates the square-root-scaled amplitudes of the average APs. The axonal delay is color-coded. Close-ups of 3 regions (labeled I, II, III), showing the average AP waveforms, are presented in the lower panels. Gray axonal contours serve as guide to the eye and are estimated by observing the spatial movement of signal peaks in consecutive movie frames (Radivojevic et al., 2017). The asterisc (*) indicates an artifact probably from a large spike at the AIS of a presynaptic neuron.

Mapping of Electrical “Footprints”

To map the electrical “footprint” (spatial distribution of extracellularly measured electrical potentials obtained with the densely packed electrodes) of neuronal units and their axonal arbors over the entire array, we used a series of configurations in which so called “fixed electrodes” were always connected to the amplifiers through the switch matrix (Frey et al., 2010), whereas the remaining “variable electrodes” were connected in sequential configurations (schematic in Figure 1E). After recording the spontaneous activity for 115 s for each of these configurations, spike sorting was performed and the electrical “footprint” was calculated using the spike-triggered average as previously described Bakkum et al., 2013; Müller et al., 2015; Deligkaris et al., 2016; Radivojevic et al., 2016, 2017;Jäckel et al., 2017.

To perform spike-triggered averaging we need to reliably record spiking activity of the neurons. Therefore, we selected the fixed electrodes at the putative locations of the AISs but imposed a spatial restriction in order to not record the same neuron twice. For selection of fixed electrodes, all electrodes were ranked according to their median negative peak amplitude (see previous section). The electrode with the highest rank was selected, afterwards all electrodes in its proximity (within 100 μm distance) were discarded from the list, and the procedure was repeated.

Spike Sorting

Even in sparse cultures some electrodes will pick up spikes from multiple neurons in their vicinity. Therefore, spike sorting was performed for each fixed electrode after recording. Waveforms were extracted for a period of 1.5 ms before to 1.5 ms after the negative peak of each event comprising 61 samples for each event. In order to extract those features that best separate the different clusters of spikes, we performed a principal component analysis. We choose the first 10 principal components as spike features (Abeles and Goldstein, 1977), containing more than 85% of the energy of the signal. KlustaKwik (Kadir et al., 2014) was used to fit a mixture of Gaussians with unconstrained covariance matrices and automatically selected the number of mixture components. Clusters containing more than 2,500 spiking events in total were assigned to distinct neurons. This number corresponds to an average of 14 spikes in each of the 179 recording configurations. With 14 samples the spike triggered averaging reduces the noise by approximately a factor114=0.27 . This procedure omits neurons with a firing rate below 14 spikes/115 s ≈ 7 spikes/min for which we could not compute a reliable footprint. The uniqueness of the detected neurons was confirmed by assessing their different footprints (Figure 6). Clusters leading to identical footprints were manually merged to single neurons, and the footprints of the merged neurons were re-calculated.

Optimal Recording Configurations for High-Throughput Scanning of Electrical Footprints

The extracellular signals originating from axons and dendrites are very small with respect to the background electrical activity and noise, so that spike-triggered averaging was applied. We developed a set of recording configurations to map the electrical footprint of several neurons in parallel by utilizing the switch matrix of our HD-MEA. The switch matrix can be dynamically configured to connect a large number, e, of electrodes to a smaller number, a, of amplifiers. In a first-order approach, one could use one electrode as trigger and the remaining electrodes to scan the neuronal footprint, which would result in a large number of configurations, cw, needed to scan the whole array, cw = e/a. For recording axonal arbors of n neurons, electrodes near the AISs of these n neurons have to be always connected to amplifiers (n fixed electrodes). The remaining amplifiers can then be connected to the remaining electrodes in several successive configurations (variable electrodes). Following this procedure, the whole array can be scanned with

c(n)=e-na-n

configurations. For n neurons we need c(n) configurations, which means on average C(n) = c(n)/n configurations for one neuron. An optimal strategy means to choose n in such a way that C(n) → min for 0 < n < a. With n < a ≪ e, the number of neurons being much smaller than the number of electrodes, we can approximate:

C(n)=e-n(a-n)n≈e(a-n )n

The right-hand side has a minimum for n = a/2. Therefore, approximately half of the amplifiers should be connected to fixed electrodes. The other half of the amplifiers can then be used to scan the whole array in 2cw configurations, or on averageCa=C(a/2)=4e/a2 configurations per neuron. The axonal arbors of a single neuron may extend over the whole array, but by scanning the axonal arbors of many neurons in parallel, the average number of configurations, Ca, per neuron is much less than the number of configurations, cw, required for scanning the whole array for large a:

Ca=4ea2≪ea=cw

This is due to the fact that increasing the number of amplifiers quadratically decreases the average time to scan a single neuron, which allows for high-throughput acquisition of axonal delay maps.

Live Imaging

Live-cell visualization of whole neurons was performed by transfection (Bakkum et al., 2013; Radivojevic et al., 2016). Transfection was performed using a pLV-hSyn-RFP plasmid from Edward Callaway (The Salk Institute, US; Addgene, #22909) and Lipofectamine 2000 (Life Technologies) in accordance with the manufacturer's protocol. A Leica DM6000 FS microscope, a Leica DFC 345 FX camera, and the Leica Application Suite software were used to produce the micrographs.

Results

To initially identify the location of axonal initial segments (AISs), the whole array was scanned by using configurations in which non-overlapping blocks (Figure 1D). The median of the negative peak for each electrode was plotted as a map showing some areas with large negative peaks (Figure 2A). Such local minima were assumed to indicate the putative (proximal) AIS locations (Bakkum et al., 2019). In low-density cultures, several AISs could be distinguished (Figure 1A). So called “fixed electrodes” were then selected as trigger electrodes at the putative locations of the AISs, while the remaining “variable electrodes” were selected in sequential configurations to map the electrical footprint of neuronal units and their axonal arbors (Bakkum et al., 2013; Müller et al., 2015; Deligkaris et al., 2016; Radivojevic et al., 2016, 2017; Jäckel et al., 2017) over the entire array (Figure 1E). After spike-triggered averaging, the amplitude and delay of axonal signals could be mapped across the whole recording area for each individual neuron (Figure 2B).

High-Throughput Scanning of Electrical “Footprints”

We used HD-MEA with 11,011 electrodes and 126 amplifiers (Frey et al., 2010). We connected 62 amplifiers to fixed electrodes and used the remaining amplifiers to scan the whole array selecting 64 electrodes at random positions, which would yield a total ofc(n)=e-na-n=11011-62126-62=1094964=172 configurations (see section Spike Sorting). However, due to the design of the switch matrix in our HD-MEA, some electrodes cannot be selected in the same configuration, which increased the number of recording configurations to 175 (see Table 1).

TABLE 1

Table 1. Comparison of high-throughput mapping of axonal delays using different types of CMOS-based HD-MEA.

In our culture, 53 of the 62 fixed electrodes recorded single-unit activity. However, only 28 axonal arbors could be identified, because 25 electrodes did not record enough spikes to reliably determine arbor structures (see Methods). The remaining 9 fixed electrodes recorded multi-unit activity, so that spike sorting was used to identify single-neuron activity and to obtain additional 40 axonal arbors. In total, we mapped the axonal arbors of 68 neurons. For further analysis we only used the n = 46 neurons with axonal arbors extending over more than 50 electrodes (Figure 6).

Identification of Axonal Arbors

Previously (see Figure 6 in Bakkum et al., 2013), axonal arbors were traced according to the occurrence of negative peaks in signal amplitudes in the spike-triggered averages of extracellularly recorded electrical signals that exceeded 5 times the background noise (sV). This method (hereafter termed “method I”) evaluates the spike-triggered averages for each electrode separately, without considering their spatial arrangement and signal correlations between neighboring electrodes. We used another method (hereafter termed “method II”) to benefit from the fact that the delay of the negative signal amplitude peak in the extracellularly measured potential originating from a single axon is very similar at neighboring electrodes it passes by. When we mapped the delay of these negative signal peaks in the spike-triggered averages (Figure 3A), the resulting map showed a distinct region of similar delays against a background of random delays (Figure 3B). This finding is due to the fact that signals originating from a common source, e.g., from an axon of the same neuron, are very similar across neighboring electrodes, whereas in the case of random signals the negative peak could occur anywhere in the interval [Δtpre, Δtpost] for which the spike triggered averaging was performed (see histogram in Figure 3C). Δtpre and Δtpost represent the boundaries relative to the spike triggering, so that the total spike-triggered average has the length T = Δtpost − Δtpre. To quantify the smoothness of the delay map, the delay was sampled in a neighborhood (compare Figure 4H) of hexagonal shape around each electrode, and the sample standard deviation Sτn was calculated. Note, that in the case of a uniform distribution over the interval [0, 1], the standard deviation of a sample is bounded by 0 ≤ s ≤ 0.5 and shows a characteristic distribution, depending on the number of observations, N, in each sample.

FIGURE 3

Figure 3. Segmentation of an axonal arbor based on the spatially correlated spontaneous activity of a single neuron. Spike-triggered averages (A) of signals from electrodes located close to the (proximal) AIS (red trace), close to axons (black) and for electrodes recording background activity and noise (gray). The negative peak at the AIS appears slightly earlier than at the trigger electrode. Mapping (B) and histogram (C) of the delay of the negative peak, [tau], showing an irregularly shaped area with a “smooth” gray value outlining the axonal arbor, which is surrounded by a “salt-and-pepper” patterned background area. Axonal signals appear at 0 ms < [tau] < 2 ms. Spike-triggered averages for N = 7 neighboring electrodes, located in the “salt-and-pepper” region (D), feature a large sample standard deviation for the delays, s[tau]n, as compared to those located in the “smooth” region (E). Mapping (F) and histogram (G) of s[tau]n. The small irregularly shaped area outlining the axonal arbor is dark, whereas the surrounding area is displayed in lighter tones. Segmentation is done by placing the threshold sthr ≈ 0.5 ms in the valley between the sharper peak (black), close to 0 ms, and the broad peak (gray) around the expectedsbackground=8/12ms (open triangle) for random delays. Mapping of electrodes, where the negative peak appears after the negative peak of the AIS (H), with s[tau]n < smin (I), which record presumably axonal signals (J). The crosshair symbol shows the location of the (proximal) AIS, the green and blue dots represent a patch of 7 neighboring electrodes located in the “salt-and-pepper” and “smooth” areas, respectively. Corresponding negative peaks are indicated by triangles of the same color.

FIGURE 4

Figure 4. Evaluation of axon segmentation based on ground truth. Mappings and Haussdorff distance,H,are shown for method I (A,D,F) and II (B,E,G). The high threshold, employed by method I, leads to a higher false negative rate and a largerH≈150[mu]m compared withH≈100 [mu]m for method II. Lowering the threshold from 5sV (A) to 3sV (D) for method I leads to more false positive electrodes, far away from the axon (black outlines) andH≈ 300 [mu]m. In contrast, method II is robust (G) with respect to an increased electrode distance (H): increasing the distance from r ≈ 18 [mu]m (B) to r ≈ 36 [mu]m (E) yields a higher true positive rate andH≈ 200[mu]m. Axons were manually traced from fluorescence images (DsRed fluorescence displayed using an inverted grayscale) (C).

This distribution shows a sharp peak around the mean standard deviation:

srandom¯=112.

There is no analytical expression for this distribution, but after a coordinate transformation of the interval, it can be approximated by a beta distributionB(α,β).

In case that axonal signals are present, the delays in the neighborhood of an electrode with a negative peak at t are distributed in the interval[t -rc;t +rc], depending on the velocity, c, of the action potential and the distance, r, between the electrodes. Therefore, the mean of the sample standard deviation for axonal delays is:

saxon¯=2rc12=rc3.

For an HD-MEA with r = 18 μm and a typical conduction velocity for short-range-projecting axons in the rat neocortex of 0.3−0.44 m/s (Lohmann and Rörig, 1994; Telfeian and Connors, 2003), a standard deviation of around 30 μs can be expected (compare with Figure 3E). At the boundary, more and more neighboring electrodes will no more pick up the axonal signal, so that the sample standard deviation shifts toward:

srandom¯=T12.

In our case, with T = 8 ms, a standard deviation of 2.3 ms for background signals would be expected (compare with Figure 3D). Empirically, the distribution of saxon can be approximated by an (truncated) exponential distributionE(λ) (see below). Therefore, axons have a distribution of sτn with a peak close to zero, which is clearly distinguishable from the distribution for the background. A threshold smin, placed at the local minimum in the sτn distribution (Figure 3G), can be used to separate both populations. The electrodes with a sτn below this threshold represent negative peaks that are consistent across neighboring electrodes (Figure 3I). If these peaks appear after the negative peak at the AIS (Figure 3H), they are assumed to originate from the axonal arbor of the neuron (Figure 3J).

For a limited number of neurons, the ground truth in the form of fluorescence images was available, so that we could compare the performance of the new method (method II; Figure 4B) with the old method (method I; Figure 4A) employing a fixed threshold at 5sV (Bakkum et al., 2013). An example of a transfected neuron is shown in Figure 4C. More electrodes were selected by method II than by method I, which could be compensated for by lowering the threshold, e.g., to 3sV (Figure 4D). In order to compare the electrode selection (a set of electrode coordinates, E) with the ground-truth axon information (a set of pixels in the image, A), we used the Haussdorff distanceH(A,E), which is commonly used in computer vision to measure how strongly two shapes differ from each other. After registering the fluorescence image to the electrode coordinate system, we calculated

H=max{supa∈Ainfe∈Ed(a,e),supe∈Einfa∈Ad(a,e)}

using the Euclidean distance d(a, e) between the coordinatee=(xe,ye)Tof an electrode recording an axonal signal and the coordinatea=(xa,ya)T of the pixel representing an axon in the fluorescence image. For smaller threshold for method I lead to a larger deviation from the ground truth than method II. This was mainly due to the fact that, more electrodes far away from the axonal arbors were selected. The new method inherently relies on adjacency and rejected these “outliers” and produced more compact maps that more closely followed the ground truth. In other words, the distributions of the feature used to classify the electrical activity as either axonal signals or as background, showed a larger overlap for method I than method II. We tested the robustness of method II against the spatial distance of the electrodes by increasing the spatial extension of the neighborhood while keeping the number of electrodes in each neighborhood constant (N = 7). When hexagonal patterns with r = 2 × 18 μm (Figure 4E) or r = 3 × 18 μm distance between electrodes were selected, the distance to the ground truth only slightly increased (Figure 4G). However, it seemed that for a carefully chosen threshold (e.g., around 4.5sV, H≈100 μm, Figure 4F), the original method performed as well as the new method.

Due to the limits of the Haussdorff distance in estimating the quality of the extracted mappings, we also calculated the receiver-operator characteristics (ROC) (Fawcett, 2006), which is possible, as both segmentation methods are binary classifiers. We fitted the respective empirical distributions of their scores with a mixture of two partially overlapping distributions representing axonal signals (“positive” class, P) and background activity (“negative” class, N) with:

1. two normal distributionsN (Figure 5A) for the segmentation with method I (Figure 5D)

P(x)≈pNN(x;μN,σN2)+pPN(x;μP,σP2) withx=logVnsn , μN < μP

2. a beta distributionB and a truncated exponential distributionE (Figure 5B) for the segmentation with method II (Figure 5E):

P(x)≈pNB(x;αN,βN)+pPE(x;λP)withx=sτnT/2,E(x;λ)=λe-λx1-e-λxfor 0 ≦ x ≦ 1

FIGURE 5

Figure 5. Comparison of the axon segmentation methods, based on the receiver-operator characteristics (ROC). Distributions and mappings are shown for the segmentation of an individual neuron, segmented by methods I (A,D) and II (B,E). The empirical distribution (NP) of amplitudes Vn (log normalized by signal noise, [sigma]V) and the sample standard deviations of the delays, s[tau]n (normalized by T/2) of the negative peaks were fitted (fit NP) to obtain the distributions of axonal signals (positive class, P) and background activity (negative class, N). The corresponding true-positive rate (TPR) and false-positive rate (FPR) were calculated for each possible threshold and plotted as ROC curve (C). The cross depicts the position of the (fixed) threshold of method I (FPR = 0.00009, TPR = 0.7), whereas the circle indicates the (adaptive) threshold of method II (FPR = 0.011, TPR = 0.85). Method II (gray shading) performs better than method I (blue shading) as shown by the larger area under the curve (AUC). This held true for all n = 46 neurons (F), and, although method I has a lower FPR (H), its TPR (G) was much lower than that of method II, as it missed out on more than 50% of the axonal signals.

For both methods, we used the fitted distributions to estimate the probability observing axonal signals (P⁺) or background activity (P⁻).

1.P-(x)= pNN(x;μN,σN2)andP+(x)=pPN(x;μP,σP2) for method I,

2.P-(x)=pNB(x;αN,βN)andP+(x)=pPE(x;λP)for method II.

We then numerically calculated the cumulative distributions in order to calculate the true positive rate (TPR) and false positive rate (FPR) for each threshold, x, and plotted them as an ROC curve (Figure 5C). The area under the curve (AUC) showed that the new method consistently had a better performance than the original method for a total of n = 46 neurons (Figure 5F). Furthermore, the automatic threshold procedure yielded a much better TPR (Figure 5G) at the expense of a slightly increased false-positive rate FPR in comparison to the original method (Figure 5H).

Discussion

HD-MEAs with more than 3,000 electrodes per square millimeter and dedicated low-noise on-chip amplifiers are suitable tools to record the electrical activity of individual axonal arbors. We first optimized a recording scheme for the switch matrix HD-MEA that relied on combinations of fixed and variable recording sites for high-throughput parallel mapping of as many neurons as possible per total recording time. The method described here shows a promising way to obtain axonal arbors at large scale from potentially all neurons during a single recording session. As an example, 68 neurons were mapped in parallel, 48 neurons of which featured large axonal arbors (Figure 6). This yield can be further improved using an HD-MEA design with an increased number of simultaneously active recording channels, decreasing the number of necessary measurement configurations and, hence, the on average the required measurement time per neuron to a few seconds (Table 1).

FIGURE 6

Figure 6. High-throughput mapping of axonal arbors. The activity map reveals the locations of the AISs of several neurons (A). The circle size indicates the square-root-scaled count of spiking events per electrode. The median negative amplitude of the spikes is color-coded with a cut-off at −200 [mu]V. Spike-triggered averaging shows the axonal footprint of 46 neurons (C). Only neurons with axonal arbors extending over more than 50 electrodes are shown. The circle size indicates the square-root-scaled amplitudes of the average APs. The axonal delay is color-coded. Gray axonal contours serve as guide to the eye and have been estimated by observing the spatial movement of signal peaks in consecutive movie frames. The axonal contours of all neurons were color-coded and combined, showing the axonal arbors in the recorded neuronal network (B).

We then developed a method to distinguish axonal signals from the background noise. Axonal arbors reveal themselves by typical waveforms of the extracellular electric field potentials (Bakkum et al., 2013, 2014; Petersen et al., 2015; Deligkaris et al., 2016). Previously (Bakkum et al., 2013), axonal arbors were traced according to the occurrence of negative peaks in signal amplitudes in the spike-triggered averages of extracellularly recorded electrical signals that exceeded 5 times the background noise. With this threshold more than 50% of the smaller axonal signals are lost (Figure 5G). Lowering the threshold can alleviate this problem, however, the likelihood that background signals are falsely assigned to axonal arbors also increases. Some electrodes measuring spurious signals can easily be detected, if they are located at comparably large distance from electrodes that record genuine axonal signals. Typically, these false positives are removed manually after visual inspection of the recorded waveforms and their isolated location in the recording area. However, this procedure is not possible for more than a few recorded axonal arbors. Therefore, we developed a new method that compares the axonal delays on neighboring electrodes. If these delays are similar, the corresponding signals most likely originated from the same axon. If there are large differences in the delay, the corresponding signals may represent background noise. Both methods, amplitude thresholding (Bakkum et al., 2013) and the new method, worked well when we validated them by comparing their results with the morphologies of axonal arbors obtained by sparse transfection (Bakkum et al., 2013; Radivojevic et al., 2016). Upon comparing both methods according to their ROCs on a larger data set of 46 axonal arbors, the new method based on local correlations outperformed the classic amplitude thresholding. Furthermore, the new method does not require manual intervention, e.g., for setting a threshold.

The example workflow demonstrated here enables high-throughput scanning of axonal arbors and mapping of their axonal delays without the need to adjust parameters for the detection of axonal signals. This method can be extended to new generations of HD-MEAs (Ballini et al., 2014; Viswam et al., 2016; Yuan et al., 2018) and can be used to obtain data from axonal arbors of thousands of neurons within a recording session of a few hours (Table 1). Our example workflow yielded axonal arbors of 68 neurons, 46 of which featured axonal arbors extending over more than 50 electrodes (Figure 6C). These axonal arbors showed considerable variation in total length, but also in local branching patterns and axonal delays (Supplemental Figure 1). Interestingly, some neurons featured footprints in the form of disconnected patches, which could be the signature of saltatory conduction of action potentials. This observation is consistent with the reported formation of nodal components in neuronal cultures (Freeman et al., 2015), but requires further experimental validation.

Our method can be used for automated selection of neurons with suitable axonal arbors for stimulation experiments (Jäckel et al., 2017), single action potential tracking (Radivojevic et al., 2017), automated patching (Obien et al., 2019), for mapping of ion receptors by local drug application (Sasaki et al., 2011) and even for automated single-cell phenotyping of axonal conduction in human iPSC-derived neuronal cultures.

Data Availability

The Hana (high density microelectrode array recording analysis) analysis pipeline is open source. All source code as well as example data to replicate the figures are available at: http://github.com/tbullmann/hdmea_axon. The example data consists of spike triggered-averages that were extracted from the raw recordings.

Ethics Statement

All experimental procedures on animals were carried out in accordance with the European Council Directive of 22 September 2010 (2010/63/EU) and had been approved by the local authorities (Animal Care and Use Committee of RIKEN; QAH24-01).

Author Contributions

TB designed the study, performed experiments, wrote the software, analyzed data, assembled figures, interpreted the results, prepared, and revised the manuscript. MR performed recording and imaging experiments. SH implemented and tested analysis algorithms. KD performed cell culture experiments. AH interpreted the results, prepared, and revised the manuscript. UF planned the study, supported the experiments, interpreted the results, and revised the manuscript.

Funding

The HD-MEA work at ETH Zurich was financially supported by the European Community through the European Research Council Advanced Grant 694829 neuroXscales (Horizon 2020) and the Swiss National Science Foundation Grant 205321_157092/1 (Axons). Financial support through the Swiss Commission for Technology and Innovation project 25933.2 PFLS-LS is also acknowledged.

Conflict of Interest Statement

UF is a co-founder of MaxWell Biosystems AG, Mattenstrasse 26, Basel, Switzerland.

The remaining 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.

Acknowledgments

We thank Alexander Stettler and Peter Rimpf for post-processing CMOS chips as well as Manuel Schröter, Roland Diggelmann and Felix Franke for helpful discussions about spike sorting.

Supplementary Material

The Supplementary Material for this article can be found online at: https://[see PDF]/articles/10.3389/fncel.2019.00404/full#supplementary-material

Supplementary Figure 1. Sholl analysis for number of electrodes (A) and axonal delay (B) for n = 46 neurons shown in Figure 6.

References

Abeles, M., and Goldstein, M. H. (1977). Multispike train analysis. Proc. IEEE 65, 762–773. doi: 10.1109/PROC.1977.10559

Alle, H., and Geiger, J. R. (2006). Combined analog and action potential coding in hippocampal mossy fibers. Science 311, 1290–1293. doi: 10.1126/science.1119055

Bakkum, D. J., Frey, U., Radivojevic, M., Russell, T. L., Müller, J., Fiscella, M., et al. (2013). Tracking axonal action potential propagation on a high-density microelectrode array across hundreds of sites. Nat. Commun. 4:2181. doi: 10.1038/ncomms3181

Bakkum, D. J., Obien, M. E. J., Radivojevic, M., Jäckel, D., Frey, U., Takahashi, H., et al. (2019). The axon initial segment is the dominant contributor to the neuron's extracellular electrical potential landscape. Adv. Biosyst. 3:1800308. doi: 10.1002/adbi.201800308

Bakkum, D. J., Radivojevic, M., Jaeckel, D., Franke, F., Russell, T. L., Frey, U., et al. (2014). “The axon initial segment is a strong contributor to a neuron's local extracellular field potential,” in Society for Neuroscience (SfN) Conference (Washington, DC).

Ballini, M., Müller, J., Livi, P., Chen, Y., Frey, U., Stettler, A., et al. (2014). A 1024-channel CMOS microelectrode array with 26,400 electrodes for recording and stimulation of electrogenic cells in vitro. IEEE J. Solid-State Circuits 49, 2705–2719. doi: 10.1109/JSSC.2014.2359219

Berdondini, L., Imfeld, K., Martinoia, S., Tedesco, M., Neukom, S., Koudelka-Hep, M., et al. (2009). Active pixel sensor array for high spatio-temporal resolution electrophysiological recordings from single cell to large scale neuronal networks. Lab Chip 9, 2644. doi: 10.1039/b907394a

Bertotti, G., Velychko, D., Dodel, N., Keil, S., Wolansky, D., Tillak, B., et al. (2014). “A CMOS-based sensor array for in-vitro neural tissue interfacing with 4225 recording sites and 1024 stimulation sites,” in Proceedings IEEE 2014 Biomedical Circuits and Systems Conference, BioCAS 2014 (Lausanne).

Bischofberger, J., Geiger, J. R., and Jonas, P. (2002). Timing and efficacy of Ca²⁺ channel activation in hippocampal mossy fiber boutons. J. Neurosci. 22, 10593–10602. doi: 10.1523/JNEUROSCI.22-24-10593.2002

Chen, T. W., Wardill, T. J., Sun, Y., Pulver, S. R., Renninger, S. L., Baohanm, A., et al. (2013). Ultrasensitive fluorescent proteins for imaging neuronal activity. Nature 499, 295–300. doi: 10.1038/nature12354

Cho, I. H., Panzera, L. C., Chin, M., and Hoppa, M. B. (2017). Sodium channel β2 subunits prevent action potential propagation failures at axonal branch points. J. Neurosci. 37, 9519–9533. doi: 10.1523/JNEUROSCI.0891-17.2017

Dan, Y., and Poo, M.-M. (2004). Spike timing-dependent plasticity of neural circuits. Neuron 44, 23–30. doi: 10.1016/j.neuron.2004.09.007

Debanne, D. (2004). Information processing in the axon. Nat. Rev. Neurosci. 5, 304–316. doi: 10.1038/nrn1397

Deligkaris, K., Bullmann, T., and Frey, U. (2016). Extracellularly recorded somatic and neuritic signal shapes and classification algorithms for high-density microelectrode array electrophysiology. Front. Neurosci. 10:421. doi: 10.3389/fnins.2016.00421

Dragas, J., Viswam, V., Shadmani, A., Chen, Y., Bounik, R., Stettler, A., et al. (2017). A multi-functional microelectrode array featuring 59760 electrodes, 2048 electrophysiology channels, stimulation, impedance measurement and neurotransmitter detection channels. IEEE J. Solid-State Circuits 52, 1576–1590. doi: 10.1109/JSSC.2017.2686580

Elgueta, C., Köhler, J., and Bartos, M. (2015). Persistent discharges in dentate gyrus perisoma-inhibiting interneurons require hyperpolarization-activated cyclic nucleotide-gated channel activation. J. Neurosci. 35, 4131–4139. doi: 10.1523/JNEUROSCI.3671-14.2015

Emmenegger, V., Obien, M. E. J., Franke, F., and Hierlemann, A. (2019). Technologies to study action potential propagation with a focus on HD-MEAs. Front. Cell Neurosci. 2019:159. doi: 10.3389/fncel.2019.00159

Eversmann, B., Jenkner, M., Hofmann, F., Paulus, C., Brederlow, R., Holzapfl, B., et al. (2003). A 128 × 128 CMOS biosensor array for extracellular recording of neural activity. IEEE J. Solid-State Circuits 38, 2306–2317. doi: 10.1109/JSSC.2003.819174

Fawcett, T. (2006). An introduction to ROC analysis. Pattern Recognit. Lett. 27, 861–874. doi: 10.1016/j.patrec.2005.10.010

Freeman, S. A., Desmazières, A., Simonnet, J., Gatta, M., Pfeiffer, F., Aigrot, M. S., et al. (2015). Acceleration of conduction velocity linked to clustering of nodal components precedes myelination. Proc. Natl. Acad. Sci. U.S.A. 112, E321–E328. doi: 10.1073/pnas.1419099112

Frey, U., Sedivy, J., Heer, F., Pedron, R., Ballini, M., Mueller, J., et al. (2010). Switch-matrix-based high-density microelectrode array in CMOS technology. Solid-State Circuits IEEE J. 45, 467–482. doi: 10.1109/JSSC.2009.2035196

Geiger, J. R. P., and Jonas, P. (2000). Dynamic control of presynaptic Ca²⁺ inflow by fast-inactivating K⁺ channels in hippocampal mossy fiber boutons. Neuron 28, 927–939. doi: 10.1016/S0896-6273(00)00164-1

Heer, F., Hafizovic, S., Ugniwenko, T., Frey, U., Franks, W., Perriard, E., et al. (2007). Single-chip microelectronic system to interface with living cells. Biosens. Bioelectron 22, 2546–2553. doi: 10.1016/j.bios.2006.10.003

Hodgkin, A. L., and Huxley, A. F. (1952). A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 117, 500–544. doi: 10.1113/jphysiol.1952.sp004764

Hoppa, M. B., Gouzer, G., Armbruster, M., and Ryan, T. A. (2014). Control and plasticity of the presynaptic action potential waveform at small CNS nerve terminals. Neuron 84, 778–789. doi: 10.1016/j.neuron.2014.09.038

Izhikevich, E. M. (2006). Polychronization: computation with spikes. Neural Comput. 18, 245–282. doi: 10.1162/089976606775093882

Jäckel, D., Bakkum, D. J., Russell, T. L., Müller, J., Radivojevic, M., Frey, U., et al. (2017). Combination of high-density microelectrode array and patch clamp recordings to enable studies of multisynaptic integration. Sci. Rep. 7:978. doi: 10.1038/s41598-017-00981-4

Kadir, S. N., Goodman, D. F., and Harris, K. D. (2014). High-dimensional cluster analysis with the masked EM algorithm. Neural Comput. 26, 2379–2394. doi: 10.1162/NECO_a_00661

Ko, K. W., Rasband, M. N., Meseguer, V., Kramer, R. H., and Golding, N. L. (2016). Serotonin modulates spike probability in the axon initial segment through HCN channels. Nat. Neurosci. 19, 826–834. doi: 10.1038/nn.4293

Kouroupi, G., Taoufik, E., Vlachos, I. S., Tsioras, K., Antoniou, N., Papastefanaki, F., et al. (2017). Defective synaptic connectivity and axonal neuropathology in a human iPSC-based model of familial Parkinson's disease. Proc. Natl. Acad. Sci. U.S.A. 114, E3679–E3688. doi: 10.1073/pnas.1617259114

Lewicki, M. S. (1998). A review of methods for spike sorting: the detection and classification of neural action potentials. Network 9, R53–R78. doi: 10.1088/0954-898X_9_4_001

Lohmann, H., and Rörig, B. (1994). Long-range horizontal connections between supragranular pyramidal cells in the extrastriate visual cortex of the rat. J. Comp. Neurol. 344, 543–558. doi: 10.1002/cne.903440405

Lüscher, H. R., and Shiner, J. S. (1990). Simulation of action potential propagation in complex terminal arborizations. Biophys. J. 58, 1389–1399. doi: 10.1016/S0006-3495(90)82485-1

Marrese, C. A. (1987). Preparation of strongly adherent platinum black coatings. Anal. Chem. 59, 217–218. doi: 10.1021/ac00128a049

Meeks, J. P., and Mennerick, S. (2007). Action potential initiation and propagation in CA3 pyramidal axons. J. Neurophysiol. 97, 3460–3472. doi: 10.1152/jn.01288.2006

Mita, T., Bakkum, D., Frey, U., Hierlemann, A., Kanzaki, R., and Takahashi, H. (2019). Classification of inhibitory and excitatory neurons of dissociated cultures based on action potential waveforms on high-density CMOS microelectrode arrays. IEEE J. Trans. Electron. Inf. Syst. 139, 615–624. doi: 10.1541/ieejeiss.139.615

Müller, J., Ballini, M., Livi, P., Chen, Y., Radivojevic, M., Shadmani, A., et al. (2015). High-resolution CMOS MEA platform to study neurons at subcellular, cellular, and network levels. Lab Chip 15, 2767–2780. doi: 10.1039/C5LC00133A

Obien, M. E. J., Deligkaris, K., Bullmann, T., Bakkum, D. J., and Frey, U. (2015). Revealing neuronal function through microelectrode array recordings. Front. Neurosci. 9:e423. doi: 10.3389/fnins.2014.00423

Obien, M. E. J., Hierlemann, A., and Frey, U. (2019). Accurate signal-source localization in brain slices by means of high-density microelectrode arrays. Sci. Rep. 9:788. doi: 10.1038/s41598-018-36895-y

Ohura, S., and Kamiya, H. (2016). Excitability tuning of axons in the central nervous system. J. Physiol. Sci. 66, 189–196. doi: 10.1007/s12576-015-0415-2

Petersen, A. V., Johansen, E. Ø., and Perrier, J.-F. (2015). Fast and reliable identification of axons, axon initial segments and dendrites with local field potential recording. Front. Cell Neurosci. 9:429. doi: 10.3389/fncel.2015.00429

Radivojevic, M., Franke, F., Altermatt, M., Müller, J., Hierlemann, A., and Bakkum, D. J. (2017). Tracking individual action potentials throughout mammalian axonal arbors. eLife 6:e30198. doi: 10.7554/eLife.30198

Radivojevic, M., Jäckel, D., Altermatt, M., Müller, J., Viswam, V., Hierlemann, A., et al. (2016). Electrical identification and selective microstimulation of neuronal compartments based on features of extracellular action potentials. Sci. Rep. 6:31332. doi: 10.1038/srep31332

Rama, S., Zbili, M., and Debanne, D. (2018). Signal propagation along the axon. Curr. Opin. Neurobiol. 51, 37–44. doi: 10.1016/j.conb.2018.02.017

Redondo, R. L., and Morris, R. G. (2011). Making memories last: the synaptic tagging and capture hypothesis. Nat. Rev. Neurosci. 12, 17–30. doi: 10.1038/nrn2963

Ritzau-Jost, A., Delvendahl, I., Rings, A., Byczkowicz, N., Harada, H., Shigemoto, R., et al. (2014). Ultrafast action potentials mediate Kilohertz signaling at a central synapse. Neuron 84, 152–163. doi: 10.1016/j.neuron.2014.08.036

Ronchi, S., Fiscella, M., Marchetti, C., Viswam, V., Müller, J., Frey, U., et al. (2019). Single-cell electrical stimulation with CMOS-based high-density microelectrode arrays. Front. Cell. Neurosci. 13:208. doi: 10.3389/fnins.2019.00208

Sasaki, T., Matsuki, N., and Ikegaya, Y. (2011). Action-potential modulation during axonal conduction. Science 331, 599–601. doi: 10.1126/science.1197598

Seidl, A. H., Rubel, E. W., and Harris, D. M. (2010). Mechanisms for adjusting interaural time differences to achieve binaural coincidence detection. J. Neurosci. 30, 70–80. doi: 10.1523/JNEUROSCI.3464-09.2010

Spruston, N. (2008). Pyramidal neurons: dendritic structure and synaptic integration. Nat. Rev. Neurosci. 9, 206–221. doi: 10.1038/nrn2286

Stange-Marten, A., Nabel, A. L., Sinclair, J. L., Fischl, M., Alexandrova, O., Wohlfrom, H., et al. (2017). Input timing for spatial processing is precisely tuned via constant synaptic delays and myelination patterns in the auditory brainstem. Proc. Natl. Acad. Sci. U.S.A. 114, E4851–E4858. doi: 10.1073/pnas.1702290114

Swadlow, H. A. (1994). Efferent neurons and suspected interneurons in motor cortex of the awake rabbit: axonal properties, sensory receptive fields, and subthreshold synaptic inputs. J. Neurophysiol. 71, 437–453. doi: 10.1152/jn.1994.71.2.437

Tasaki, I. (1939). The electro-saltatory transmission of the neve impulse and the effect of narcosis upon the nerve fiber. Am. J. Physiol. 127, 211–227. doi: 10.1152/ajplegacy.1939.127.2.211

Telfeian, A. E., and Connors, B. W. (2003). Widely integrative properties of layer 5 pyramidal cells support a role for processing of extralaminar synaptic inputs in rat neocortex. Neurosci. Lett. 343, 121–124. doi: 10.1016/S0304-3940(03)00379-3

Tsai, D., Sawyer, D., Bradd, A., Yuste, R., and Shepard, K. L. (2017). A very large-scale microelectrode array for cellular-resolution electrophysiology. Nat. Commun. 8:1802. doi: 10.1038/s41467-017-02009-x

Viswam, V., Chen, Y., Shadmani, A., Dragas, J., Bounik, R., Milos, R., et al. (2016). “2048 action potential recording channels with 2.4 μ Vrms noise and stimulation artifact suppression,” in 2016 IEEE Biomedical Circuits and Systems Conference (BioCAS) (Shanghai).

Viswam, V., Jäckel, D., Jones, I., Ballini, M., Muller, J., Frey, U., et al. (2014). “Effects of sub-10 μm electrode sizes on extracellular recording of neuronal cells,” in MicroTAS (San Antonio, TX).

Viswam, V., Obien, M. E. J., Franke, F., Frey, U., and Hierlemann, A. (2019). Optimal electrode size for multi-scale extracellular-potential recording from neuronal assemblies. Front. Neurosci. 13:00385. doi: 10.3389/fnins.2019.00385

Wainger, B. J., Kiskinis, E., Mellin, C., Wiskow, O., Han, S. S., Sandoe, J., et al. (2014). Intrinsic membrane hyperexcitability of amyotrophic lateral sclerosis patient-derived motor neurons. Cell Rep. 7, 1–11. doi: 10.1016/j.celrep.2014.03.019

Yuan, X., Emmenegger, V., Obien, M. E. J., Hierlemann, A., and Frey, U. (2018). “Dual-mode microelectrode array featuring 20k electrodes and high snr for extracellular recording of neural networks,” in Proceedings of 2018 IEEE Biomedical Circuits and Systems Conference, BioCAS 2018 (Cleveland, OH).