Microglial peri‐somatic abutments classify two novel types of GABAergic neuron in the inferior colliculus

Emerging evidence suggests functional roles for microglia in the healthy, mature nervous system. However, we know little of the cellular density and ramified morphology of microglia in sensory systems, and even less of their inter‐relationship with inhibitory neurons. We therefore conducted fluorescent multi‐channel immunohistochemistry and confocal microscopy in guinea pigs of both sexes for Iba1, GAD67, GFAP, calbindin, and calretinin. We explored these markers in the inferior colliculi (IC), which contain sub‐regions specialized for different aspects of auditory processing. First, we found that while the density of Iba1+ somata is similar throughout the IC parenchyma, Iba1+ microglia in dorsal cortex are significantly more ramified than those in the central nucleus or lateral cortex. Conversely, Iba1+ ramifications in ventral central nucleus, a region with the highest density of GAD67+ (putative GABAergic) neurons in IC, are longer with fewer ramifications. Second, we observed extensive abutments of ramified Iba1+ processes onto GAD67+ somata throughout the whole IC and developed novel measures to quantify these. Cluster analyses revealed two novel sub‐types of GAD67+ neuron that differ in the quantity of Iba1+ somatic abutments they receive. Unlike previous classification schemes for GAD67+ neurons in IC, these clusters are not related to GAD67+ soma size. Taken together, these data demonstrate that microglial ramifications vary between IC sub‐regions in the healthy, adult IC, possibly related to the ongoing demands of their niche. Furthermore, Iba1+ abutments onto neuronal somata are a novel means by which GAD67+ neurons can be classified.


| INTRODUCTION
Inhibition is an essential element of neural processing. GABAergic inhibition is prevalent in the auditory system, particularly in the principal auditory midbrain nuclei, the inferior colliculi (IC). The IC has a tonotopic topography that can be divided into distinct sub-regions. The central nucleus (CNIC) is dominated by neurons sharply tuned to simple auditory stimuli. The dorsal cortex (DCIC) has much broader frequency tuning and receives extensive corticofugal input and is specialized for synaptic plasticity (Bajo & Moore, 2005;Bajo et al., 2010;Herbert et al., 1991;Winer et al., 1998). The other major sub-region is the lateral cortex (LCIC) which exhibits polysensory tuning (Aitkin et al., 1978). Around a quarter of neurons in IC are GABAergic (Merchán et al., 2005), but we know little of how they interact with resident glia.
Microglia are integral cell types with varying densities and morphologies throughout the brain (Lawson et al., 1990;Matias et al., 2019), suggesting adaptation to the functional demands of their local milieu. Microglia have been shown to interact with neurons during "normal" processing, including sensing and responding to local signaling (Deemyad et al., 2018;Pocock & Kettenmann, 2007;Schafer et al., 2012;Wake et al., 2009). Despite the essential role the IC plays in auditory processing, little is known about microglia in the IC, nor how they may vary between sub-regions.
Understanding how GABAergic cell types vary throughout the brain, to specialize for distinct functions, is a key area of neuroscientific study (Freund & Buzsaki, 1996;Tremblay et al., 2016). Most investigations into sub-types of GABAergic neurons focus on the cells per se, including their morphology, electrophysiological firing characteristics, expression of cytoplasmic calcium binding proteins, and RNA transcriptome. However, another approach is to characterize and classify GABAergic neurons based on differences in the afferent peri-somatic inputs they receive (Beebe et al., 2016;Ito et al., 2009).
Here we take advantage of the functional organization of the IC to investigate anatomical inter-relationships of Iba1+ microglia and GAD67+ (putative GABAergic) neurons. We employed multi-channel fluorescence immunohistochemistry and confocal microscopy and tested our hypotheses that (a) microglia exhibit sub-region-specific morphological adaptations to their niche, and (b) abutments from Iba1+ processes onto GAD67+ somata can be used to classify GAD67+ neuron sub-types in IC.

| Regulation and ethics
All animals were housed and procedures performed in accordance with the terms and conditions of a license (PPL 60/3934) issued by the UK Home Office under the Animals (Scientific Procedures) Act 1986. Ethical approval was provided by the Local Ethical Review committee at Newcastle University and the School of Healthcare Science Ethics Committee at Manchester Metropolitan University.

| Animals
Results are described from four adult (one at 4 months old, three at 6 months old) outbred, tricolor guinea pigs (Cavia porcellus) of both sexes (three male, one female). Similar data was derived from for all four cases. Data from male and female guinea pigs are reported together throughout the main body text as there were no substantial differences found between sexes (Beery, 2018;Beery & Zucker, 2011). Data derived from the males (8,209,8,210,8,224) and the female (8,212) can be found in the Supplementary Material and images taken from each case, including the sex of the animal, are outlined in figure legends. Animal weights on the morning of each respective perfusion ranged from 675 to 867 g. We aimed to minimize the number of animals used and their suffering.

| Anesthesia and tissue processing
Animals were deeply anesthetized with sodium pentobarbital (i.p. injection; Euthanal, Merial; 200 mg/ml, 2 ml volume). After 5 min, the pedal withdrawal and blink reflexes were assessed to confirm a deep plane of anesthesia. This was followed by transcardial gravity perfusion with 500 ml of 0.1 M heparinized PBS followed by 500 ml of freshly made 4% paraformaldehyde in 0.1 M PBS. Both solutions were of pH 7.2 directly prior to use.
Brains were removed with rongeurs (Micro Friedman, 0.8 mm jaws, WPI) and post-fixed in 30% sucrose in 4% paraformaldehyde, for at least 3 days at 4°C. Once brains sank, they were cut in the coronal plane with a razor blade through the parietal and temporal lobes, at the rostro-caudal location of the medial geniculate. The tissue block was then placed in an embedding mold (Peel-a-way; Shandon), covered in embedding medium (OCT; Agar Scientific) and frozen at −80°C. The 60-µm sections were taken on a cryostat (HM560, Microm) and collected in 12-well plates in cryoprotectant (30% sucrose, 30% ethylene glycol, 1% polyvinyl pyrrolidone-40 in 0.1 M PBS) and stored at −20°C until use (Olthof et al., 2019;Watson et al., 1986).

| Antibody characterization
The following primary antibodies were used:
Mouse anti-GFAP (1:500; monoclonal; clone G-A-5; G3893; lot# 045M4889V; Sigma; RRID: AB_477010)-according to the manufacturer, this antibody is raised against an epitope from the C-terminus of GFAP in purified pig spinal cord (Debus et al., 1983;Latov et al., 1979). The antibody has been shown to recognize a single band of approximately 50kDa and reacts with homologous, conserved residues across mammals (Lorenz et al., 2005). The use of this antibody has been demonstrated in many species, including mouse (Komitova et al., 2005), rat (Lennerz et al., 2008;Sanchez et al., 2009), tree shrew (Knabe et al., 2008), guinea pig (Kelleher et al., 2011(Kelleher et al., , 2013 and human (Toro et al., 2006). Labeling observed in this study was consistent with these studies and the known morphology of astrocytes.
Rabbit anti-calbindin D-28k (1:1,000; polyclonal; AB1778; lot# 2895780; Millipore; RRID: AB_2068336)according to the manufacturer, this antibody recognizes a single band at 28kDa in human, mouse, and rat brain tissues. It does not bind to calretinin and pre-adsorption of diluted antiserum with calbindin removed all the labeling in human brain (Huynh et al., 2000). Previous labeling of mouse olfactory bulb (Kotani et al., 2010), rat piriform cortex (Gavrilovici et al., 2010) and guinea pig enteric nervous system (Liu et al., 2005) all showed highly selective cytoplasmic labeling of neurons. We observed labeling consistent with previous reports.

| Fluorescence immunohistochemistry
Sections through the superior colliculus and the rostral-most third along the rostro-caudal axis through the IC were first used to optimize labeling. Data are presented from two sections from each animal in the middle third of the IC along the rostro-caudal axis, which contained the CNIC, DCIC and LCIC using an approach similar to that reported by Beebe et al. (2016). The location of each section through the rostrocaudal axis was referenced to an atlas of the guinea pig brainstem (Voitenko & Marlinsky, 1993).
All steps in the labeling protocol involved continuous gentle agitation of sections. Free-floating sections were brought to room temperature and washed 3 × 5 min in PBS. Sections were blocked and permeabilized in 5% normal goat serum (Vector) and 0.05% Triton X-100 (Sigma) in PBS for 1 hr. Following blocking, a cocktail of primary antibodies was added to the blocking solution and applied to sections overnight at room temperature. The next day, sections were washed 3 × 5 min in PBS and incubated for 2 hr in appropriate secondary antibodies (Invitrogen; 1:250 in blocking solution). For double labeling of Iba1 and GAD67, goat anti-rabbit AlexaFluor 488 and goat anti-mouse AlexaFluor 568 were used. For double labeling of calbindin and GFAP, goat anti-rabbit AlexaFluor 488 and goat anti-mouse AlexaFluor 647 were used. For triple labeling of Iba1, GSL1 (pre-conjugated rhodamine fluorophore) and GFAP, goat anti-rabbit AlexaFluor 488 and goat anti-mouse AlexaFluor 647 were used. Sections were then mounted on slides and coverslipped using Vectashield (Vector Labs, H-1000) and kept at 4°C until imaged. All experiments had control slides where the primary, secondary or both the primary and secondary antibodies were excluded. This allowed detection of autofluoresence and any aspecific signal and ensured only labeling from primary and secondary binding to epitope targets was acquired.

| Image acquisition
Sequentially acquired micrographs were taken with a confocal microscope (Leica SP5) using a wide field stage and | 5817 zoom functions. Micrographs were acquired via a 40× objective (NA = 1.25) for images of the entire cross-section of the IC, and a 63× objective (NA = 1.4) for panoramas within IC sub-regions. Images of the IC section were taken using 5 µm equidistant slices in the Z-plane to produce maximum intensity tiled projections (40× objective, pixel size: x&y = 0.7583 µm; z-slice distance = 2.97µm). For GAD67 and Iba1 panoramas, 5-row × 6-column (432 × 552 µm) tiled images were taken using 1 µm z-slices and rendered as maximum intensity projections (63× objective, pixel size: x&y = 0.2406µm; z-slice distance = 1µm). A hybrid scanning detector (HyD 3) was used to optimize image quality, which negated offsetting properties. The overall laser power used for all images was set at 30%. For individual wavelengths, the 488-exciation power was set at 30% and the 555-exciation power was set at 54% for both 40× and 63× objective images. The gain setting for all images and objectives was set at its lowest (10%). The pinhole for 40x objective images was set at 121.33 µm for the 488-and 199.22 µm for the 555-exciation power. For 63× objective images these were set at 124.37 µm for the 488-exciation power and 150.11 µm for 555-exciation power. We defined sub-regions of IC by overlaying images of the whole IC on top of those presented in Coote and Rees (2008) and used their landmarks to demark the three sub-regions of IC. We generated panoramas clearly within each sub-region, with each positioned a marked distance from other panoramas in that section. All image acquisition settings were kept the same for panoramas imaged and analyzed across all animals and sub-regions.

| Image analyses
All analyses were undertaken using Fiji ImageJ. For Iba1+ soma density estimates, whole IC section images were subject to manual counts. The peripheral borders of the IC were delineated and a contour drawn, and each image cropped to its respective contour. To make fair comparisons between cases, 450 µm 2 grids were placed across each whole IC section image and centered on the middle pixels of each micrograph in ImageJ. Only those grids which were filled entirely by stained parenchymal tissue were subject to counts. Comparisons were then made between cases, such that only grids that were present in images from all four animals were included in calculation of group means and standard deviations per grid.
Panoramas were rendered as maximum intensity projections and analyzed for initial Iba1+ soma counts and percentage field of view analyses, which were calculated using the analyze particles plugin after thresholding. Parameters in the analyze particles plugin were set to exclude labeling under an area of 20 µm 2 . Only Iba1+ soma counts and percentage field of view analyses were rendered as maximum intensity projections and subjected to this analysis. Masks of all included counts were visually checked to ensure labeling was morphologically specific for ramified Iba1+ somata and processes before inclusion in statistical analyses. This provided a basis for comparison of sub-regional Iba1+ microglial labeling. Individual 1 μm z-slices were reconstructed across three-dimensions and were analyzed by (a) manual GAD67+ soma counts, (b) individual Iba1+ Sholl and (c) skeleton analyses. All panorama micrographs used in Sholl and skeleton analyses were first processed by filtering monochrome images using a median pixel (1.5) filter and then thresholded to binary by implementing the IsoData algorithm.
For 3-D volumetric Sholl analyses, individual Iba1+ microglia were cropped and a series of equidistant radiating 1 µm concentric spheres were plotted from the center of the soma to the furthest radiating extent of ramification. Each intersection was counted. The Skeletonize algorithm was used to display a one-pixel thick framework. The Analyze Skeleton plugin calculated number of branches and branch lengths across three dimensions.
Five GAD67+ somata were selected randomly from each sub-region panorama. A z-stack was collected throughout the z-plane of the section using a selection criterion that the entirety of the soma must be contained within the x, y and z dimensions of the panorama. GAD67+ soma perimeters were manually contoured in each z-slice after ensuring they met the inclusion criteria above. To determine Iba1+ labeling that abutted the contour, the Huang algorithm was implemented on both markers and if there were no pixels between the two thresholded signals, it was counted as an Iba1+ process abutting a GAD67+ soma. Iba1+ microglia processes and their respective somata were included in this analysis based on their proximity to the selected GAD67+ soma. Because of the ramified, non-overlapping territories of parenchymal Iba1+ cells, all Iba1+ cells that had processes long enough to contact the selected GAD67+ soma were included in this analysis.
To examine the nature of Iba1+ putative interactions with the selected GAD67+ somata, we calculated five different dependent quantitative measures of GAD67+ somata abutments (n = 20 per brain section, 40 per IC sub-region). (a) We measured GAD67+ soma diameter by quantifying the number of 1 µm slices of GAD67+ somatic labeling present through the z-axis-that is, the z-distance from the first to last z-planes containing GAD67+ somatic labeling (hereafter named Maximum Diameter of GAD67+ soma; MDG). (b) To quantify how many Iba1+ cells abutted each GAD67+ soma, the number of individual microglial cells that made contact with each measured soma were manually counted through all z-planes (hereafter named Number of Iba1+ cells contacting GAD67+ soma; NCG). (c) To determine the 3D extent of microglial abutments, the percentage of the x, y and z dimensions of each whole GAD67+ soma abutted by Iba1+ labeling was measured. This | 5818 was done by quantifying the percentage of the total distance of Iba1+ labeling abutting the perimeter of every GAD67+ z-slice across the z-axis, then averaging to produce a measure for each GAD67+ cell (hereafter named Percentage of GAD67+ Soma Perimeter abutted by Iba1+; PSP). (d) As a result of the orders of magnitude over which IC neuron somata vary in size, a process of normalization was undertaken to compensate for the Iba1+ abutments they receive. For every GAD67+ soma, the total distance of Iba1+ labeling abutting the perimeter of every GAD67+ z-slice across the z-axis was calculated. The total distance of all Iba1+ abutments onto GAD67+ somata perimeters was then normalized by dividing this value for each cell by the value of the cell with the greatest distance of all Iba1+ abutments. This process compensated for any potential bias toward larger GAD67+ somata having more Iba1+ abutments (hereafter named Normalized Length of GAD67+ abutments made; NLG). (e) Finally, we counted the number of distinct processes from each Iba1+ cell abutting each GAD67+ soma (hereafter named Number of Distinct Iba1+ Processes abutting GAD67+ soma; NDP).
While the above analyses considered the inter-relationships of Iba1+ and GAD67+ labeling from the perspective of the latter, to consider the former, we also selected 10 Iba1+ cells per panorama and counted the number of GAD67+ somata abutted by each Iba1+ microglia.

| Statistical analysis
Data were collected in Excel spreadsheets. Statistical hypothesis testing was performed in Prism 7 (GraphPad). Factorial analyses were conducted using the non-parametric Kruskall-Wallis ANOVA with sub-region as the factor in all cases. Where appropriate, post hoc tests with Dunn's method were conducted. For post hoc analyses the α was Šidák corrected for multiple comparisons. Spearman's rank correlations were used to investigate potential associations between dependent variables.
Principal component and two-step cluster analyses were conducted in SPSS v25 (IBM). The two-step cluster analysis employed Euclidean distance measures with Schwarz's Bayesian clustering criterion. All GAD67+ somata were successfully classified by this analysis. A chi-squared test was used to analyze the ratio of somata in IC sub-regions in each cluster. All reported p values are exact and two tailed.

| GFAP+ astrocytes and Iba1+ microglia form the glia limitans externa in IC
We first sought to identify the distribution of GFAP+ astrocytes and Iba1+ microglia in adult guinea pig IC. Coronal, 60 µm sections showed pronounced GFAP+ and Iba1+ labeling of the glia limitans externa lining the dorsal and lateral borders of the IC (Figure 1a). Extensive labeling was also distributed medially, lining the cerebral aqueduct, with ramified GFAP+ astrocytic processes radiating into the periaqueductal gray, as well as the commissure of the IC. Interestingly, we found no GFAP+ astrocytes throughout the IC parenchyma, save for sparse labeling in the outermost layers of the DCIC and LCIC. Conversely, ramified Iba1+ microglia tiled the parenchyma in non-overlapping domains with similar density throughout the IC, as quantified in Figure 1b.
Combining Iba1+ and GFAP+ labeling with the fluorescent-conjugated lectin GSL1 revealed extensive peri-vascular labeling along putative penetrating arteries and arterioles ( Figure 1c). Somata with neuronal morphology expressing cytoplasmic calbindin or calretinin were distributed in the outermost regions of the cortices of the IC, matching previous reports Zettel et al., 1997) and in close proximity to GFAP+ processes ( Figure 1d).
These findings demonstrate that many aspects of IC glial organization mirror those reported in other brain regions, with GFAP+ astrocytes forming the glia limitans externa and lining adjacent to blood vessels (Hafidi & Galifianakis, 2003), alongside Iba1+ microglia. The uniform spacing of Iba1+ microglial somata throughout IC parenchyma suggested potential interactions with neurons. As GAD67+ cell density is known to vary between sub-regions of IC, this led us to question if there are interactions between Iba1+ and GAD67+ somata, and furthermore, whether these differ between sub-regions of IC. We therefore focused our analyses on GAD67+ somata and the ramified processes of Iba1+ microglia.
We predicted that the stronger Iba1+ microglia labeling in DCIC neuropil was primarily because of a greater number and extent of ramifications compared to other sub-regions of IC. To test this, 3D volumetric Sholl analyses were implemented on a total of 64 microglia per sub-region. Labeling from each individual z-slice was quantified and these combined to generate volumetric 3D datasets. The slices were projected as 2D maximum intensity projections to enable a contour to be drawn around the maximum extent of the soma and ramified processes (Figure 5a). Labeling was then cropped and split into individual slices to remove background/non-cellular labeling from every slice (Figure 5b-i). Slices were then merged as three-dimensional projections (Figure 5j) and thresholded to generate binary images (Figure 5k). The number of F I G U R E 1 Iba1+ Microglia and GFAP+ astrocytes form the glia limitans externa and peri-vascular borders of IC but only Iba1+ microglia tile the parenchyma. (a) Tiled maximum intensity projection confocal micrograph of GFAP+ astrocytes (magenta) and Iba1+ microglia (green) in IC of animal 8,209. Scale bar 400 µm. Note that GFAP labeling is restricted to the peripheral borders and penetrating vessels while Iba1 is evenly distributed throughout. (b) Quantification of Iba1+ somata counts per 450 µm 2 grid (shown inset) within the parenchyma across all four cases. Numbers on top refer to individual animals. Brackets refer to sexes. We found no clear differences between data derived from different sexes so data were pooled. (c) Confocal micrograph showing GFAP, Iba1 and GSL1 (blue) labeling of a penetrating arteriole coursing into IC. Scale bar 50 µm. (d) Confocal micrograph showing a calbindin (green) expressing neuron in the outer layers of the LCIC surrounded by GFAP+ peri-somatic processes. Scale bar 10 µm | 5820 intersections at every micrometer distance radiating from the center of the soma was calculated (Figure 5l). Binary thresholds were also skeletonized to derive information about the shape and structure of ramifications, including the number of branches and maximum branch length (Figure 5m). Iba1+ microglia were more ramified in DCIC than mid-CNIC, LCIC, or VCNIC (Figure 6a). The total number of intersections of every ramification (independent of distance from the soma) had a median value of 324 in DCIC (IQR = 269-375). This was greater than in mid-CNIC (284; ±222-323), LCIC (247; ±210-312) and VCNIC (262; ±207-318) (Figure 6b). A Kruskal-Wallis one-way ANOVA on ranks, with sub-region as the factor found the differences between sub-regions of IC were likely a real effect (H(3) = 22.44; p < 0.0001). Post hoc analyses via Dunn's tests showed differences between DCIC and the other three sub-regions were likely to be a real effect (mid-CNIC p = 0.0222, LCIC p = 0.0001, VCNIC p = 0.0006).
Conversely, the average branch length of every process, defined as the distance covered by any ramification of skeletonized Iba1+ microglia without branching, followed the opposite trend. The longest processes were found in VCNIC (median = 4.87 µm; IQR = 4.35-5.19) (Figure 6d

| Iba1+ putative abutments reveal two novel types of GAD67+ neurons in IC
We hypothesized that GAD67+ neurons, which are known to receive a variety of types and densities of presynaptic contacts (Beebe et al., 2016;Ito et al., 2009), may also receive varying amounts of Iba1+ abutments. To examine the nature of Iba1+ putative interactions with GAD67+ somata, we quantified five dependent variables from each GAD67+ soma based on analyses performed as in Figure 7 (see Section 2). These variables (MDG, NCG, PSP, NLG, and NDP) were calculated from micrographs such as the representative example in Figure 8a, which shows a GAD67+ cell being abutted by two Iba1+ microglia. A correlation matrix revealed weak associations between MDG (i) and the other four dependent variables (Figure 8b). However, there were strong correlations between NCG, PSP, NLG, and NDP.
As these variables were only weakly correlated with GAD67+ soma diameter, we further investigated whether a multivariate analysis could better explain the observed distributions. We conducted a principal component analysis for the five variables in all sub-regions of the IC. The data showed a clear dissociation between MDG in one cluster and F I G U R E 2 GAD67+ somata vary in density between sub-regions of IC. (a) Tiled confocal micrograph showing Iba1+ microglia tiling IC parenchyma from animal 8,209. Scale bar 400µm. Same scale for panels b and c. (b) GAD67+ neuropil can be seen to demark the medial and ventral borders of the IC. GAD67+ somata are found throughout IC but vary in density. (c) Merge of a (green) and b (magenta). (d) Borders of IC sub-regions were delineated using those defined by Coote and Rees (2008). Panoramas were located within distinct sub-regions of IC that could be clearly distinguished from one another | 5821 the other variables, which clustered together separately from MDG (Figure 8c). Both clusters were categorized using a standard correlation coefficient of >0.5 as a cut-off value, which showed that one cluster was explained by only the MDG variable, while the other cluster had significant contributions from each of the other four variables. This was also true for all sub-region-specific analyses in IC (Figure 8d). We then conducted a two-step cluster analysis including all five variables. We employed Euclidean distance measures with Schwarz's Bayesian clustering criterion (Figure 8e). Contrasting iterations up to 15 clusters, the analysis found two clusters demonstrated the best explanatory power, displaying good (0.5) silhouette measures of cohesion and separation ( Figure 8f). All 160 GAD67+ somata were classified into one of the two clusters determined by the two-step analysis.
Representative examples of each cluster found in each sub-region are shown in Figure 9. Note that while Iba1+ processes are present in all eight example images, only F I G U R E 3 Representative panoramas show differences between Iba1+ and GAD67+ somata in sub-regions of IC in case 8,212. Scale bar in (a) 50 µm. Same scale for all panels. Maximum intensity projections of tiled confocal panoramas in all sub-regions show (left column; a,d,g,j) Iba1+ microglia tiling the parenchyma with similar density. Conversely, labeling of GAD67+ somata (middle column, b,e,h,k) reveals varied soma density between subregions. Merging both labels (right column; c,f,i,l) reveals intercalating of Iba1+ processes with GAD67+ labeling | 5822 those in cluster 1 exhibit extensive peri-somatic abutments. There were 23 cases in cluster 1, and 137 in cluster 2. To visualize the contribution of each of the variables, they were plotted as a function of MDG (Figure 10a-d). These scatterplots revealed a dissociation with little overlap between the two clusters using PSP (Figure 10a). There was almost perfect discrimination between the clusters using NLG (Figure 10b). Conversely, NCG had a significant degree of overlap with little difference between clusters (Figure 10c), while NDP had little overlap between distributions ( Figure 10d). To compare the ability of each of these variables to independently discriminate between the two clusters, we conducted ROC analyses (Figure 10e). These revealed that while each variable had an area under the curve >0.5, the three variables relating to the nature of Iba1+ processes abutting GAD67+ somata had the best discriminatory power, with NLG exhibiting an area under the curve of 0.99.

| Iba1+ putative interactions with GAD67+ neurons show little difference between sub-regions of IC
We explored whether any of the variables or clusters identified had a relationship to the sub-region of IC in which they were located. MDG did not vary between sub-regions ( Figure 11a) (H(3) = 5.3; p = 0.151). A small difference was found for PSP ( Figure 11b) (H(3) = 9.9; p = 0.019). Post hoc Dunn's tests suggested a potential real difference between DCIC and VCNIC (p = 0.012). However, there was extensive overlap between the distributions so there is a reasonable chance this may not be a true effect. The NCG (Figure 11c We interpret this as being as a result of the greater density of GAD67+ neurons in VCNIC (Figure 4a).

| DISCUSSION
These findings reveal that Iba1+ microglia are present throughout the parenchyma of the IC in the healthy, mature auditory system and form abutments with GAD67+ neuronal somata. Furthermore, we have demonstrated that F I G U R E 4 GAD67+ soma counts and Iba1+ labeling show opposing trends. (a) Box plot showing GAD67+ somata counts in each sub-region across cases. The VCNIC consistently had the highest number of GAD67+ somata while numbers in DCIC and LCIC were lower than mid-CNIC. (b) Iba1+ soma counts were similar between panoramas but (c) Iba1+ percentage field of view showed much more labeling was present in DCIC | 5823 the length of peri-somatic abutments made offers a means to classify between two newly identified GAD67+ cell sub-types.
Taking advantage of the specialized and functionally diverse sub-regions of IC, we found the first evidence, to our knowledge, of differences in microglial morphology between sub-regions. Interestingly, Iba1+ microglia in DCIC, a sub-region known to receive a greater proportion of glutamatergic corticofugal contacts, are more ramified than those in other sub-regions of IC. We also developed a new analysis approach, which was used to investigate the number and extent of abutting microglial processes onto GAD67+ somata. Multivariate cluster analyses were applied to variables and revealed two new types of GAD67+ neuron in IC, distinguished by the extent of Iba1+ microglial abutments onto their soma. Contrastingly, and in agreement with the findings of Hafidi and Galifianakis (2003), we found that GFAP in IC only labels astrocytes lining the glia-limitans externa, and some sparse parenchymal astrocytes in the outermost layers of the DCIC and LCIC.

| Significance of IC sub-regional differences in microglial morphology
Sensory processing can be interpreted through triad models of organization, based on observations that central sub-regions of sensory pathways are dominated by ascending innervation, producing brisk responses at short latencies to simple stimuli, such as in the CNIC. Adjacent are at least two sub-regions that are driven by more complex stimuli at longer response latencies. One of these sub-regions typically receives a diversity of polymodal inputs (LCIC), while the other is primarily driven by descending, top-down afferents (DCIC). These sub-regions exhibit differing connectional as well as cytoarchitectonic and chemoarchitectonic organization (Sweet F I G U R E 5 Process for volumetric Sholl and skeleton analyses to quantify morphological characteristics of Iba1+ microglia. Examples shown from animal 8,210. (a) Example of maximum intensity projection of panorama tiled confocal micrograph showing Iba1+ microglia in IC. The 6-column × 5-row field of view tiled image was taken using 1 μm z-slices through the full extent of the section. Iba1+ cellular labeling was subject to Sholl and skeleton analyses, one at a time (yellow box inset). Scale bar 50 μm. (b-i) Examples from eight extracted slices (of 24 total) from (a) showing high-resolution imaging of soma and processes surrounded by non-cellular labeling, which was cropped from each individual slice. (j) Maximum intensity projection of cropped images. (k) Each slice was thresholded to binarize processes. (l) 3D Sholl analyses were performed at 1μm resolution throughout x, y and z planes from binarized images. (m) Three-dimensional skeletonized framework projected in 2D, to be analyzed for number of branches and max branch length. Scale bars 20 μm (b-m) | 5824 et al., 2005). In IC, there is longstanding evidence for a triad model of organization using a diversity of methodologies, including electrophysiology (Syka et al., 2000), two photon imaging (Wong & Borst, 2019), fMRI (Baumann et al., 2011;De Martino et al., 2013), histology (Faye-Lund & Osen, 1985), immunohistochemistry (Coote & Rees, 2008) and tract tracing of projections (see Malmierca & Hackett, 2010).
The present study provides the first evidence that Iba1+ microglia also exhibit differences between IC sub-regions. There is longstanding evidence that microglia exhibit morphological differences throughout the brain (Lawson et al., 1990) and there are neuroanatomical differences between CNIC, DCIC and LCIC, including increased cytoarchitectural and myeloarchitectural density in CNIC (Faye-Lund & Osen, 1985). There is also a more defined laminar organization of the CNIC than the DCIC, in part owing to the neuronal morphologies present, as well as coursing fibers (Oliver & Morest, 1984). As the DCIC is known to have larger neurons than CNIC (Faye-Lund & Osen, 1985), one might speculate that the increased number of branching ramifications we observed in DCIC ( Figure 6) correlates with larger somata. However, the LCIC is also known to have larger neurons than CNIC, but we found similar number of microglial branching ramifications to mid-CNIC and VCNIC. This suggests that differing microglial morphologies may relate to other aspects of local processing, such as at synapses (Trapp et al., 2007;Tremblay et al., 2010;Wake et al., 2009). The greater number of branching ramifications in DCIC may relate to the differing nature of local synaptic processing therein compared to CNIC and LCIC. Whether this reflects differences in auditory processing or more general maintenance of ongoing signaling at synapses on GABAergic neurons, or indeed neurons of other neurochemical types, remains to be explored.
Another explanation may be the higher density of neurons in CNIC. There have been previous suggestions that microglial density and neuronal density are inversely related throughout the brain (Lawson et al., 1990). However, we found no strong evidence for differences in microglial soma density between IC sub-regions using two different levels of approach (Figures 1b and 4b). This was in spite of much greater numbers of GAD67+ somata in VCNIC than other sub-regions (Figures 3 and 4a). A third possibility is that there are relationships related to myeloarchitectural differences between sub-regions. However, DCIC and VCNIC exhibit similar, high levels of axons of passage compared to mid-CNIC, arguing against such an explanation.

GABAergic neurons
As the auditory pathway contains a large proportion of inhibitory neurons, with around a quarter of neurons being GABAergic (Merchán et al., 2005;Oliver et al., 1994), understanding their structure, function, and organization is a F I G U R E 7 GAD67+ somata were selected that were contained entirely within the x, y, and z coordinates of the imaged tissue section. (a) Example z-slice in VCNIC from animal 8,224, labeled for GAD67 (magenta) and Iba1 (green). Example GAD67+ soma (yellow arrowhead) was selected for analysis. Scale bar 50 µm. Gray, dashed box shows zoomed in regions in panels (b-n). (b) Maximum intensity projection of Iba1+ cell through 13 z-planes in zoomed region around area selected in (a), allowing visualization of microglial processes and soma with a GAD67+ soma. Scale bar 10 µm. (c-n) Example, 1 µm slices though the rostro-caudal axis in the z-plane. Soma perimeter length was measured in each slice (gray outlines). Yellow arrowheads in (i and j) point to magnified region in (o-r). Scale bar 10 µm (c-n). (o) Magnified region from (i) demonstrating an Iba1+ process abutting GAD67+ soma. Scale bar 1 µm (p-r). (p) Thresholding using the Huang algorithm in ImageJ for each label reveals abutment of Iba1+ process onto GAD67+ soma, with no pixels separation and a 1-2 pixel thick plexus of double labeling (white), allowing quantification of the location and length of abutment. The length of abutment was measured along the GAD67+ perimeter. Each distinct process abutting the GAD67+ soma was counted and the total number of Iba1+ somata abutting the GAD67+ soma were measured. (q) Magnified region from (j), showing Iba1+ process in close proximity to but not abutting GAD67+ soma. (r) Thresholding each label reveals clear separation with unlabeled pixels between. This demonstrates that the measures taken facilitated positive and negative identification of abutments through the z-plane | 5826 question of fundamental importance. Previous approaches to classifying GAD67+ neurons in IC have focused on soma size (Ono et al., 2005;Roberts & Ribak, 1987a, 1987b coupled with axo-somatic inputs (Ito et al., 2009), peri-neuronal nets (Beebe et al., 2016) or cytoplasmic calcium binding protein expression (Engle et al., 2014;. The present study shows that while there is merit to these approaches, other features of GABAergic neuronal sub-types exist. Indeed, we have discovered that GAD67+ neurons can be classified into two distinct clusters based on the Iba1+ abutments onto their soma (Figures 8-10).
It may be claimed that cluster analyses, as well as other classification approaches, can produce artificial discriminations between data that reside along continua, such as has been reported for frequency response areas in IC (Palmer et al., 2013).
In many cases, putting data through cluster analyses will produce clustered data, irrespective of whether these clusters represent meaningful differences. However, we permuted our cluster analysis through various iterations and found that not only could the data be best clustered by two clusters, but that those clusters had strong explanatory power, with good silhouette measures of cohesion and separation and good cluster quality (Figure 8e,f). Furthermore, visualization of representative examples of these clusters showed that those in cluster 1 clearly received a greater number and length of abutting Iba1+ processes onto their somata than those in cluster 2 (Figure 9). ROC analyses revealed that the two clusters could be distinguished by three variables that quantified aspects of Iba1+ processes, but not by the number of Iba1+ microglial somata abutting each GAD67+ soma. One of the variables, NLG, had an ROC area under the curve F I G U R E 8 (a) Representative GAD67+ soma (magenta; outlined in gray) receiving abutting somatic processes from two Iba1+ microglia (green) in case 8,212. Scale bar 20 µm. (b) Correlation matrix showing Spearman's rank correlation coefficient for each combination of measures derived from putative interactions between Iba1+ microglia and GAD67+ labelling. Note lower values in first column. (c) Principal component analysis of the five variables in (b). GAD67+ soma diameter clustered separately from the other four variables, which grouped together. (d) As (c) but for each sub-region in IC, showing similar findings in all sub-regions, suggesting robust clustering throughout IC. (e) The largest difference between Schwarz's Bayesian criterion and the ratio of Euclidean distance measures is observed for two clusters, demonstrating best explanatory power. (f) Silhouette measures of cohesion and separation were compared from 2 to 15 clusters, with two clusters demonstrating the highest cluster quality, with 0.5 indicating good quality clustering. These data demonstrate that not only did two clusters exhibit the best available explanatory power, but also good clustering quality, reflecting real underlying differences in Iba1+ abutments onto GAD67+ somata between clusters | 5827 of 0.99, indicating near perfect classification using this variable alone.
Interestingly, the two identified clusters of GAD67+ neurons did not differ in their relative proportion between the sub-regions in IC. Future work may investigate the differing afferent neural inputs to and efferent targets of these cell types, to identify likely physiological and connectional differences between clusters.

| Technical considerations
The use of primary antibodies in less studied species such as guinea pig can be challenging because of potential differences in epitopes and when not adequately controlled for, may lead to spurious observations (Schonbrunn, 2014). This is important when using exploratory approaches as in the present study, to ensure all analyses are predicated on specific and selective labeling (Voskuil, 2017). We therefore conducted extensive control experiments, excluding primary antibody only, secondary antibody only and both antibodies, to ensure analyses were based on true labeling.
The lack of GFAP+ astrocytes in IC parenchyma replicates a similar report in gerbil (Hafidi & Galifianakis, 2003). As in that report, we found extensive GFAP labeling in peri-vascular regions and the glia limitans externa, but the lack of GFAP+ astrocytes in IC parenchyma does not exclude the possibility that astrocytes reside throughout IC. Indeed, there are reports employing S100β and SR101 that have revealed a network of putative astrocytes throughout IC (Ghirardini et al., 2018), that are distinct from GFAP+ astrocytes. However, S100β and SR101 are not selective markers of astrocytes, which hampers interpretability in studies trying to exclusively label astrocytes (Hill & Grutzendler, 2014;Steiner et al., 2007). There are a diversity of other non-GFAP markers that may reveal the distribution of distinct astrocyte sub-types throughout the IC.
There are inherent challenges in measuring Iba1+ abutments onto GAD67+ somata, such as potential loss of information when using 2D measurements from confocal z-slices and quantifying these via 3D measures of entire somata. To minimize the impact of this, we sampled at fine z-slice steps (1 μm). However, our 3D measures did not allow us to distinguish between abutments from one or multiple Iba1+ somata in the z-plane (rostro-caudal axis) which may lead to conflicting interpretations of our measures of Iba1+ length of abutments onto somata. However, this was minimized by implementing two additional measures that are not subject to this same potential bias-namely PSP and NDP. Both variables demonstrated similar ROC area under curve values of 0.91 and 0.97 respectively (Figure 10e). Importantly, NCG did not differ between clusters (Figure 10c), demonstrating that it is the amount of Iba1+ abutment GAD67+ somata receives, and not the origin of these processes, that is important in distinguishing the two clusters.
Functional differences between astrocytes in CNIC and the outer layers of DCIC and LCIC have been suggested previously via 3-chloropropanediol-induced lesions, which selectively destroyed the former but not the latter (Willis F I G U R E 9 Example GAD67+ somata taken from each subregion of IC (maximum intensity projection of five z-slices through the center of each soma in the z-plane, along the rostro-caudal axis), with one example from each designated cluster from the two-step analysis and examples presented from all four animals in the study. Scale bar 10 µm applies to all images. Note the greater length of Iba1+ abutments surrounding GAD67+ somata in cluster 1 (left column), while adjacent Iba1+ processes in cluster 2 make few abutments onto the somata (right column) | 5828 et al., 2003, 2004). The present study leads to the speculation of fundamental gliochemical and physiological differences that may relate to the sub-region-specific roles microglia and astrocytes play in their local milieu (Lawson et al., 1990;Olah et al., 2011). Recently, RT-PCR of single IC astrocytes revealed expression of functional inhibitory neurotransmitter transporters GlyT1, GAT-1, and GAT-3 (Ghirardini et al., 2018). Sub-regional differences in GAD67+ neurons in the present study suggest that GABAergic and glycinergic signaling released from and received by glial cells may also exhibit such variations throughout IC and perhaps in other structures.

| CONCLUSION
We have described, for the first time, that Iba1+ microglia tile the adult IC parenchyma and have discovered subregional differences in the morphology of microglia in IC. Furthermore, multivariate statistical approaches revealed two new clusters of GAD67+ neurons which can be distinguished based on the total amount of Iba1+ abutments they receive onto their soma. Our findings demonstrate morphological diversity amongst IC microglia and differential interactions with GAD67+ somata.