Development of a High-Throughput Pipeline to Characterize Microglia Morphological States at a Single-Cell Resolution

Visual Abstract

We thank the Ciernia Lab and Pavlidis Lab members for their thoughtful feedback and suggestions during lab meetings throughout the progression of this project.We also thank Wai Hang (Tom) Cheng, whose help was instrumental with learning how to image on the Axioscan slidescanner and getting started with microglia morphology analysis; Nicholas Michelson, whose help was invaluable when troubleshooting code in ImageJ for various features of MicrogliaMorphology; and Dylan Terstege, who generously provided the materials for FASTMAP alignment to the Allen Brain Atlas before they were published.We also thank Dr. Brian MacVicar for sharing his lab's Cx3cr1-GFP mice with us, which we used for the 2xLPS in vivo experiments.We are grateful for the resources made available through the Dynamic Brain

Introduction
Clinical and postmortem studies support a role for altered microglial function as a critical component of multiple brain disorders, including schizophrenia, autism spectrum disorder, and Alzheimer's disease (Tetreault et al., 2012;Suzuki et al., 2013;Hansen et al., 2018;Zhuo et al., 2023).In addition to their resident immune functions, microglia play critical roles to establish and maintain normal brain function including the regulation of neuronal cell number (Cunningham et al., 2013), shaping of brain circuitry (Schwarz et al., 2012;Jobling et al., 2018;Anderson et al., 2019;Nguyen et al., 2020;Mehl et al., 2022), and fine-tuning of neuronal connections (Schafer et al., 2012;Zhan et al., 2014;Wang et al., 2020), processes that have all been shown to be fundamentally disrupted in many brain disorders (Salter and Stevens, 2017;Bachiller et al., 2018;Dubbelaar et al., 2018;Hammond et al., 2018;Lenz and Nelson, 2018;Zhu et al., 2023).Microglia directly communicate with other cell types and modulate their function by releasing and responding to various molecular substrates in the brain environment including cytokines, chemokines, and neurotransmitters (Mead et al., 2012;Liu et al., 2016;Vainchtein et al., 2018;Ferro et al., 2021;Logiacco et al., 2021;Salvador et al., 2021).Thus, dysregulated microglial responses can disrupt the homeostatic brain environment and normal communication across cell types, ultimately altering brain function and behavior.
As rapid responders to their local environments, microglia exhibit a dynamic range of phenotypes defined by multiple parameters including transcriptomic signatures (Hammond et al., 2019;Li et al., 2019), epigenomic regulation (Ayata et al., 2018;Ciernia et al., 2018;Meleady et al., 2023), and changes in cellular morphology.Microglia are classically thought to engage in immune functions that are mirrored in morphology, where homeostatic microglia exhibit a ramified morphology that switches to an ameboid, unramified form in response to inflammatory signals in their environment.This morphological switch is thought to allow for increased mobility to sites of infection or Circuits cluster and the NeuroImaging and NeuroComputation Core at the UBC Djavad Mowafaghian Centre for Brain Health (RRID:SCR_019086) that supported this work.injury, efficient phagocytosis, and release of cytokines into the microenvironment, all functions that are characteristic of a "pro-inflammatory" state.However, a mechanistic link between reduced branching and increased inflammatory function has only been shown recently, where increased P2RY12 potentiation of the THIK-1 channel caused decreased microglial ramifications and increased activity of the Il1b inflammasome in response to tissue damage (Madry et al., 2018).Cdk1-mediated microtubule remodeling has also recently been shown to be required for efficient cytokine trafficking and release and transformation of microglia from ramified to ameboid forms after LPS exposure in vitro and in situ (Adrian et al., 2023).In contrast, findings from other studies display a reverse relationship, where ramified microglia have been shown to phagocytose synapses during adult neurogenesis (Sierra et al., 2010;Kamei and Okabe, 2023), and ameboid microglia show reduced phagocytic capabilities in epilepsy (Abiega et al., 2016;Paolicelli et al., 2022).The relationship between form and function is clearly not as dichotomous as historically thought, and there has been an increasing effort in the field to move away from dualistic classification of microglial function and toward a clearer understanding and appreciation of the heterogeneous states of microglia (Dubbelaar et al., 2018;Paolicelli et al., 2022).
Microglia morphology has been shown to be highly context and signal dependent, displaying various forms ranging from ameboid-like under inflammatory states to hyper-ramified in mouse models of stress-induced depression, accelerated aging, and Alzheimer's disease (Beynon and Walker, 2012;Raj et al., 2014;Hellwig et al., 2016;Madry et al., 2018).Subsets of microglial populations have also been shown to display rod-like morphologies characterized by long, thin processes protruding from oval-shaped somas that retract neuron-adjacent planar processes in response to diffuse brain injury (Taylor et al., 2014).Furthermore, microglia are highly motile and never truly quiescent even in healthy conditions, constantly extending out protrusions to scan their environments for pathogens and other harmful signals, as revealed by in vivo two-photon imaging studies of the mouse brain (Davalos et al., 2005;Nimmerjahn et al., 2005;Bernier et al., 2020).Different morphological forms have been observed to be conserved across species (Geirsdottir et al., 2019) and spatially distributed across brain regions (Lawson et al., 1990;Wierzba-Bobrowicz et al., 1998;Mittelbronn et al., 2001;Tan et al., 2020;Brandi et al., 2022), neurodevelopmental timepoints (Paolicelli et al., 2011;Parakalan et al., 2012;Cengiz et al., 2019), and various pathological contexts (Davies et al., 2017;Savage et al., 2019;Meng et al., 2024).Microglial morphological forms including ramified, rod-like, hypertrophic, and ameboid have been commonly observed in both humans and mice, and microglia have been shown to display similar dendritic morphology across species (Geirsdottir et al., 2019;Savage et al., 2019;Paolicelli et al., 2022).
Microglia morphology can be explained using various measures of branch length, branching complexity, territory span, and cell circularity, which together comprise distinct sets of morphology measures that are changing in conjunction with each other.Nevertheless, studies often selectively report changes in individual features such as the number of branches or cell area alone, ultimately depicting an incomplete or biased representation of changes in microglia morphology that fail to capture true morphological states.Therefore, an analysis approach that considers contributions of all feature measures together to explain various axes of morphology is necessary to gain a better understanding of a microglia's actual morphological state and relationship to cellular function.A plethora of tools exist to analyze microglia morphology (Fernández-Arjona et al., 2017;York et al., 2018;Young and Morrison, 2018;Salamanca et al., 2019;Clarke et al., 2021) but often involve significant effort to extract meaningful information at a larger scale, as many approaches that allow for the analysis of morphological information at a cellular resolution involve manually choosing and segmenting individual cells within an image.This introduces a potential bias for which cells are selected and vastly limits the feasible sample size for analysis.While there has been considerable progress in the development of toolsets which automate these timeconsuming steps (Clarke et al., 2021;Melo et al., 2023), there has been less development, transparency about, and availability of methods to analyze the resulting morphological measures (Reddaway et al., 2023).The underlying code and datasets for the majority of published toolsets are not openly available nor well documented, further limiting the uptake and progression of the most up-to-date toolsets by the larger research community (Reddaway et al., 2023).
Here, we describe an accessible and open-source morphology analysis toolset: MicrogliaMorphology (ImageJ tool) and MicrogliaMorphologyR (R package), which supply the field with new tools to systematically evaluate the heterogeneity of microglia morphological states by considering 27 different measures of microglia morphology.To demonstrate use cases for our toolset, we characterized and analyzed microglia morphology in an experimental model of repeated immune stimulation by peripheral lipopolysaccharide (LPS) administration, a commonly used model (Wendeln et al., 2018;H. Jung et al., 2023) which induced population shifts in the four major classes of microglia morphology in our dataset: ramified, hypertrophic, ameboid, and rod-like.Application of MicrogliaMorphology and MicrogliaMorphologyR by the scientific community will yield novel insights into microglia morphology differences in the brain at a single-cell resolution and in a spatially resolved manner across various experimental models and research questions.Additionally, the toolset is not limited to microglia morphology alone but can also be applied in the same way to characterize morphology of other cell types.

In vivo experiments
All experiments were conducted in accordance with the University of British Columbia guidelines, with approval from the Canadian Council on Animal Care Committee.Mice were housed in groups of two to four on a regular 12 h light/dark cycle, and all experiments were performed in regular light during the mouse's regular light cycle.CX3C motif chemokine receptor 1 (Cx3cr1) is a commonly used microglial marker that is expressed on microglia and other immune cells, and cells under control of the endogenous Cx3cr1 locus express GFP (IMSR_JAX:005582; S. Jung et al., 2000).Eight-week-old male and female Cx3cr1-GFP mice bred on a C57BL/6J background (n = 3/condition, 2 females and 1 male) were intraperitoneally injected with 0.5 mg/kg lipopolysaccharide (LPS; lipopolysaccharides from E. coli O55:B5, Sigma-Aldrich L5418) or vehicle solution (PBS; phosphate-buffered saline, Fisher BioReagents BP3991) once every 24 h for 2 d.
1xLPS experiments (Extended Data Fig. 3-1).Mice were housed in groups of two to four on a reversed 12 h light/dark cycle.Because mice are nocturnal animals, all experiments were performed in red light during the mouse's dark cycle when they are most active.Eight-week-old male and female C57BL/6J mice (n = 2/sex/condition) were intraperitoneally injected once with 1.0 mg/kg lipopolysaccharide (LPS; lipopolysaccharides from E. coli O55:B5, Sigma-Aldrich L5418) or vehicle solution 1xPBS (phosphate-buffered saline, Fisher BioReagents BP3991).

Tissue collection
Three hours (2xLPS experiments) or 24 h (1xLPS experiments) after the final injection, mice were quickly anesthetized with isoflurane and transcardially perfused with 15 ml of 1xPBS before brains were extracted for downstream immunohistochemistry experiments.Extracted brains were immersion fixed in 4% paraformaldehyde for 48 h before cryoprotecting in 30% sucrose for 48 h prior to cryosectioning.Cryoprotected brains were then sectioned at 30 µm on the cryostat, collected in 1xPBS for long-term storage, and processed for immunohistochemistry.

Imaging
All mounted brain sections were imaged on the ZEISS Axioscan 7 microscope slide scanner at 20× magnification with a step size of 0.5 µm using the z-stack acquisition parameters within the imaging software (ZEISS ZEN 3.7).During image acquisition, Extended Depth of Focus (EDF) images were created using maximum projection settings and saved as the outputs in the final .czifiles.Maximum projection EDF images only compile the pixels of highest intensity at any given position in a z-stack to construct a new 2D image which retains the 3D information.Using ImageJ, we created .tiffimages of each fluorescent channel from .czifiles and selected and saved .tiffs of brain regions of interest (ROIs) to use as input for downstream morphological analysis in MicrogliaMorphology and MicrogliaMorphologyR.We focused our analyses on multiple brain regions including the hippocampus, frontal cortex, and striatum, as well as subregions within them (Extended Data Fig. 4-1D).Images of coronal brain sections containing these regions were aligned to the Allen Brain Atlas (mouse.brain-map.org;Pinskiy et al., 2015) using the ImageJ macro FASTMAP (Terstege et al., 2022).

MicrogliaMorphology
MicrogliaMorphology is designed to be a user-friendly ImageJ macro that wraps around existing ImageJ plugins AnalyzeParticles, Skeletonize (2D/3D), and AnalyzeSkeleton and is written using the ImageJ macro (IJM) language (Schneider et al., 2012).All supporting code for MicrogliaMorphology is available on GitHub at https://github.com/ciernialab/MicrogliaMorphology, and a detailed video tutorial which includes relevant troubleshooting steps is available on YouTube at https://www.youtube.com/watch?v=YhLCdlFLzk8.After MicrogliaMorphology is installed into the user's ImageJ plugins folder as described in the GitHub repository, it will appear in the Plugins dropdown menu from the ImageJ toolbar, where it can be clicked on to begin the user prompts.In Step 1, users are prompted to measure datasetspecific parameters, which are critical because every imaging dataset is prepared and acquired differently and thus requires user input to determine what parameters most appropriately and accurately capture microglia morphology within individual datasets.Image thresholding parameters including the method and radius considered for auto local thresholding are important to determine the most appropriate thresholding method which captures full, single microglial cells without losing branching connectivity.Area ranges that accurately describe single microglial cells are important to exclude any artifacts of 2D representation in the EDFs such as cell particles and incomplete cells or multiple cells that are overlapping.Dataset-specific image thresholding and cell area parameters determined by these initial steps are then called to within the macro to inform downstream steps of MicrogliaMorphology.In the subsequent steps, the only user input involves following prompts to select input folders to call from and output folders to write to, with the option of batch processing.All protocols, computation, and analysis described below have been written to be automated within MicrogliaMorphology, unless otherwise specified.
Step 2 after determining dataset-specific parameters is thresholding input .tifffluorescent images.MicrogliaMorphology cleans up and thresholds input images according to standard protocol (Young and Morrison, 2018): images are binarized and converted to grayscale before the brightness and contrast are enhanced, unsharp mask filter applied to clarify existing detail, despeckle function applied to remove noise, and auto local or auto thresholding applied to the resulting images.Then, a second despeckle step is applied before dilation, and connection steps are applied to connect branches, after which outliers are removed.The final images are then used as input for Step 3, which uses AnalyzeParticles and ROI manager functions to create and save new images of single cells which pass the area criteria specified in Step 1. Generation of single-cell images in this step allows for the measurement of morphology measures at an unprecedented cellular resolution from hundreds to thousands of different microglia cells, each with unique identifiers.
The single-cell images are used as input for Step 4, which uses Skeletonize (2D/3D; Lee et al., 1994) and AnalyzeSkeleton (Arganda-Carreras et al., 2010) to generate measures of different morphology features including maximum branch length, average branch length, and numbers of end point voxels, junction voxels, triple points, branches, junctions, slab voxels, and quadruple points for every cell.MicrogliaMorphology saves these outputs as individual .csvfiles, which contain all skeleton measures for every individual cell in the dataset, marked by unique identifiers.Step 5 involves FracLac (Karperien, 1999(Karperien, , 2013)), a plugin separate from MicrogliaMorphology which uses fractal analysis to measure additional morphology features including the width of bounding rectangle, maximum radius from hull's center of mass, maximum span across hull, diameter of bounding circle, maximum radius from circle's center of mass, perimeter, mean radius, mean radius from circle's center of mass, area, foreground pixels, height of bounding rectangle, max/min radii from circle's center of mass, relative variation (CV) in radii from circle's center of mass, span ratio of hull (major/minor axis), max/min radii from hull's center of mass, relative variation (CV) in radii from hull's center of mass, density of foreground pixels in hull area, and circularity.Because FracLac is incompatible with the IJM language, it was not integrated into our MicrogliaMorphology macro, and Step 5 must be completed using FracLac-specific user prompts and common parameters outlined in Young and Morrison (2018).Steps to batch process the single-cell images generated from MicrogliaMorphology using FracLac are written out in detail on our GitHub page for MicrogliaMorphology. Importantly, the unique identifiers for each cell are retained in the FracLac output, allowing for integration with the MicrogliaMorphology measures.The final AnalyzeSkeleton output from MicrogliaMorphology (Step 4) and FracLac (Step 5) are merged using a custom function in MicrogliaMorphologyR to generate a final, master .csvfile containing measures for 27 different morphology features for every individual cell.
An additional feature within MicrogliaMorphology is the ColorByCluster feature, which allows the user to color the microglia cells in the original .tiffinput images by their k-means cluster identification (see below, MicrogliaMorphologyR).This is a unique feature of MicrogliaMorphology that allows the user to visually validate their morphological clusters and gain insight into their spatial distribution in the brain (Fig. 3E).All original input immunofluorescent .tiffimages containing FASTMAP subregion ROIs, thresholded images generated by MicrogliaMorphology, thresholded single-cell images generated by MicrogliaMorphology, and all analysis code for the data presented in this paper are available on the Open Science Framework (OSF) at https://osf.io/8t5c2/.Finally, tidied up datasets that were used for analysis in this paper are included as part of MicrogliaMorphologyR and can be loaded with the package (data_1xLPS_mouse, data_2xLPS_mouse, data_2xLPS_mouse_fuzzykmeans, and data_ImageTypeComparison).

MicrogliaMorphologyR
MicrogliaMorphologyR is an R package that wraps several existing packages including tidyverse, Hmisc, pheatmap, factoextra, lmerTest, lme4, Matrix, SciViews, ggpubr, glmmTMB, DHARMa, rstatix, and gridExtra (Bates et al., 2015(Bates et al., , 2023;;Auguie, 2017;Brooks et al., 2017;Kuznetsova et al., 2017;Kolde, 2019;Wickham et al., 2019;Kassambara, 2020;Grosjean, 2022;Hartig, 2022;Harrell, 2023;Kassambara, 2023a,b).While our ImageJ macro, MicrogliaMorphology, facilitates the semiautomated measurement of 27 individual morphology features at a single-cell level, our complementary R package, MicrogliaMorphologyR, allows for analysis and visualization of this data to characterize microglia morphological states and gain insight into their relevance in experimental models.Functions within MicrogliaMorphologyR integrate correlation analyses and statistical modeling approaches and are used in conjunction with dimensionality reduction by principal component analysis and k-means clustering to characterize morphological states and quantify population shifts in the experimental model of choice.MicrogliaMorphologyR also includes exploratory data analysis functions to generate heatmap and boxplot visualizations of data in flexible ways including at the single-cell level, animal level, and experimental condition level.Furthermore, MicrogliaMorphologyR includes functions for generating quality control metrics on input data such as identifying values that dominate and disproportionately skew feature distributions, data normalization options, and performing linear mixed effects modeling, ANOVA, and other statistical analyses on the input dataset.All source code for MicrogliaMorphologyR and descriptions of functions can be found on GitHub at https://github.com/ciernialab/MicrogliaMorphologyR.

Statistical analysis of microglia markers
To test for LPS-induced morphological population shifts in cluster membership, we first calculated the percentage of cells in each morphology cluster for every brain region and antibody for every mouse using the "clusterpercentage" function within MicrogliaMorphologyR. To assess how cluster membership changes with LPS treatment across microglia markers, we filtered for each individual brain region before fitting a generalized linear mixed model using a beta distribution to model the percentage of cluster membership as a factor of cluster identity, treatment, and antibody interactions with MouseID as a repeated measure ["percentage ∼ Cluster*Treatment*Antibody + (1|MouseID)"] using the "stats_cluster.animal" function from MicrogliaMorphologyR, which is wrapped around the glmmTMB R package (Brooks et al., 2017).The beta distribution is suitable for values like percentages or probabilities that are constrained to a range of 0-1.We fit separate models for each brain region individually because the model would not converge when considering Cluster* Treatment*Antibody*BrainRegion interactions together, and we were more interested in Cluster*Treatment*Antibody interactions for our analysis.Tests between treatments (PBS vs 2xLPS) were corrected for multiple comparisons across clusters and antibodies (∼Treatment|Cluster|Antibody) using the Bonferroni's method and q values <0.05 were considered statistically significant (Extended Data Table 4-2).
To highlight the utility of MicrogliaMorphologyR for assessing changes in individual morphology measures, we used the "stats_morphologymeasures.animal" function, which is wrapped around the "lmer" function from the lme4 package in R, to fit linear mixed effects models for each brain region to the measure values as a factor of treatment and antibody interactions with MouseID as a repeated measure ["Value ∼ Treatment*Antibody + (1|MouseID)"].We focused on cell area as an example.Tests between treatments (PBS vs 2xLPS) were corrected for multiple comparisons across antibodies (∼Treatment|Antibody) using Bonferroni's method and q values <0.05 were considered statistically significant (Extended Data Table 4-4).

Code accessibility
The datasets generated and/or analyzed during the current study are freely available to download, and all relevant materials including input images, data, and supporting code used in this study are available on the OSF website at https://osf.io/8t5c2/.The ImageJ and R code underlying both MicrogliaMorphology and MicrogliaMorphologyR, as well as the R code for the analysis demonstrated in this paper are all available through GitHub repositories at https://github.com/ciernialab/MicrogliaMorphology and https://github.com/ciernialab/MicrogliaMorphologyR.

Morphology analysis toolset: ImageJ macro MicrogliaMorphology
Using ImageJ plugins (Karperien, 1999(Karperien, , 2013;;Young and Morrison, 2018), we have developed an accessible and userfriendly ImageJ tool, MicrogliaMorphology, that automates the characterization of a vast range of morphology features from hundreds to thousands of individual microglial cells (Fig. 1A).We compiled the ImageJ code into an ImageJ macro format such that users can simply click through options for their morphology analysis and specify where to read and write files.The following steps are automated within the MicrogliaMorphology ImageJ macro such that users can easily perform morphology analysis by following user prompts within ImageJ.Briefly, immunofluorescence images are binarized, cleaned up to remove any background noise, and thresholded according to standard protocol (Young and Morrison, 2018) before individual cell images are created based on area measurements to exclude any artifacts that arise as a product of 2D image representation.The newly generated single-cell images are then used downstream as input for ImageJ plugins AnalyzeSkeleton (Arganda-Carreras et al., 2010) and FracLac (Karperien, 1999(Karperien, , 2013;;Young and Morrison, 2018; Fig. 1A) to measure 27 unique morphology features from individual cells in a high-throughput and semiautomated manner (Fig. 1A).Importantly, MicrogliaMorphology saves the ROI information of each individual cell in the output file so that the user has the option of linking the individual cells back to their spatial locations in the original input images.We have enabled this option through the complementary ColorByCluster functions in MicrogliaMorphology and MicrogliaMorphologyR.
As MicrogliaMorphology is wrapped around the ImageJ plugin FracLac (Karperien, 1999(Karperien, , 2013)), which can only handle 2D input, it is limited to the analysis of 2D images.While 3D reconstructions of microglia offer the benefit of finer-grained detail for morphology analysis, the size of the generated datasets and the time and resources necessary to construct and process such data limits the applicability of 3D measurements to smaller areas, inevitably making the analysis comparatively low throughput.Two-dimensional image options include maximum projection EDF images, which only compile the pixels of highest intensity at any position in a z-stack to construct a new 2D image which preserves some 3D information, versus individual 2D images sampled that only capture single-plane information.To assess how accurately microglial morphology is represented in 3D image forms as compared with 2D image forms (EDF, single-plane 2D), we analyzed just the skeletal measures (Fig. 2E) from 20 individual microglial cells manually isolated from z-stack images of the hippocampus of female mice.AnalyzeSkeleton, one of the ImageJ plugins that MicrogliaMorphology is wrapped around, is able to measure skeletal features in both 2D and 3D, enabling direct comparisons across different image forms.Five cells of each morphological type (ameboid, ramified, hypertrophic, rod-like) were manually classified and selected out from the original 3D z-stacks, from which EDF images were generated and 2D single-plane images from the center of the stack were saved to create the test dataset.As expected, the raw values of all skeletal measures quantified decreased considerably when measured from 3D to either of the 2D image forms (Extended Data Fig. 2-1A).While cells with more ramification and cell branching complexity (hypertrophic, ramified classes) were better captured in 3D (Fig. 2A), relative differences in morphology across the four different forms were well conserved across image types (Fig. 2B).After dimensionality reduction, PC2, which mostly captured variability described by the maximum branch length (Extended Data  correlated between 3D and EDF forms but not 3D and 2D forms (Fig. 2C, Extended Data Fig. 2-1B, Extended Data Table 2-1).Cell branching complexity, as described by numbers of junctions, end point voxels, branches, slab voxels, and triple points, was highly correlated across image types (Fig. 2D, Extended Data Fig. 2-1B).Although the number of junctions was also highly correlated across image types, junction voxels, or the numbers of actual pixels which make up the junctions, had low correlation scores when comparing 3D images with either of the 2D image types (Fig. 2D), which is unsurprising, as a 3D image would retain much more of this kind of information.Numbers of quadruple points, which would describe the most complex type of junction, was the least well captured in 2D representations (Extended Data Fig. 2-1B).As expected, these results together indicate that EDF images better retain 3D skeletal information than 2D forms do but that relative group differences between ameboid, hypertrophic, ramified, and rod-like morphologies are still maintained across image forms.

Morphology analysis toolset: R package MicrogliaMorphologyR
Once the 27 morphological features are measured from individual microglia using the MicrogliaMorphology ImageJ macro and FracLac, they are concatenated into a final output data file which can be analyzed further to gain insight into microglia morphology changes in any given experimental model.We have provided an R package, MicrogliaMorphologyR, which contains a set of functions that implement one set of approaches for such an analysis in R. Using MicrogliaMorphologyR, the user can conduct exploratory data analysis to generate visualizations of their own data in flexible ways including heatmaps of how morphological features vary across morphological clusters and boxplots of how morphological populations shift at the subject-level across treatment conditions.Functions within MicrogliaMorphologyR are used in conjunction with principal component analysis and k-means clustering to gain further insight into microglia morphology features, classify individual cells by their morphological states, and allow for the quantification of morphological population shifts in experimental contexts.MicrogliaMorphologyR also includes functions for generating quality control metrics on input data such as identifying values that dominate and disproportionately skew feature distributions, data normalization options, visually exploring different sources of variability in the dataset, and performing ANOVA analysis, linear mixed effects modeling, and other statistical analyses on the input dataset.

Application to 2xLPS mouse dataset
To demonstrate the utility of MicrogliaMorphology and MicrogliaMorphologyR, we describe an experimental dataset collected from the brains of 8-week-old Cx3cr1-GFP mice.Male and female mice were injected peripherally with two daily injections of vehicle or 0.5 mg/kg lipopolysaccharide (LPS), a major structural component of gram-negative bacteria that is commonly used to induce and study microglial responses in the brain (Fig. 1B).To capture a diverse range of microglia morphologies across multiple brain regions, we focused on profiling the frontal cortex, striatum, and hippocampus from six individual mice (n = 3/treatment, 2 females and 1 male).Using MicrogliaMorphology, we were able to quantify 27 morphological features from a total of 43,332 individual microglial cells, which made up our input dataset for analysis.
Within the pool of morphology measures, we defined distinct sets of highly correlated features that describe different aspects of morphology, including branching complexity, area and territory span, branch length, and cell shape (Fig. 3A, Extended Data Table 3-1).Spearman correlation analysis across all 27 features revealed that those which describe branching complexity (number of end point voxels, junction voxels, triple points, branches, junctions) were highly correlated to each other compared with the other features [abs(R) ≥ 0.8; p < 0.05; Extended Data Table 3-1].Similar correlations were observed for features that describe area and territory span (width of bounding rectangle, maximum radius from hull's center of mass, maximum span across hull, diameter of bounding circle, maximum radius from circle's center of mass, perimeter, mean radius, mean radius from circle's center of mass, area, number of slab voxels, foreground pixels, height of bounding rectangle), branch length (maximum branch length, average branch length), and cell shape (max/min radii from bounding circle's and hull's center of mass, relative variation in radii from bounding circle's and hull's center of mass).As expected, cell circularity was highly negatively correlated with span ratio of the bounding hull (major/minor axis), a measure whose higher value indicates greater cell oblongness (Fig. 3A, Extended Data Table 3-1).The relationships observed among the 27 morphology features were consistently captured in another LPS dataset collected under entirely different conditions (reversed light cycle, single 1.0 mg/kg LPS exposure, 24 h collection time; Extended Data Fig. 3-1A, Extended Data Table 3-2), demonstrating that MicrogliaMorphology is able to consistently and robustly capture different aspects of microglia morphology across experimental models.

Dimensionality reduction and soft clustering
To define morphological states from our 27-feature dataset, we performed dimensionality reduction using principal component analysis followed by fuzzy k-means clustering on the first three principal components (PCs), which together explained 84.6% of the variability in the dataset (Extended Data Fig. 4-1B).Spearman's correlation of the first three PCs to the 27 features showed that each PC was differentially correlated to and described by different sets of morphology features [abs(R) ≥ 0.75; p < 0.05; Fig. 3B, Extended Data Table 3-4].PC1 was highly positively correlated to features describing branching complexity and territory span, meaning that individual cells with greater branching complexity or area had higher PC1 scores (Fig. 3B, Extended Data Table 3-4).PC1 was also highly positively correlated to density of foreground pixels in hull area, which describes a cell's occupancy within its territory and can be a proxy for soma and/or branch thickness.Taking these correlations together, PC1 captured the variability in the dataset driven by branching complexity, territory span, and territory occupancy.In a similar manner, PC2 captured variability driven by cell circularity and cell shape and PC3 captured variability driven by average branch length (Fig. 3B, Extended Data Table 3-4).In line with our feature analysis (Fig. 3A, Extended Data Fig. 3-1A), we also observed that the PCs were similarly described by the same distinct sets of features in the 1xLPS dataset (Extended Data Fig. 3-1B, Extended Data Table 3-3).
The first three PCs were used as input downstream for fuzzy k-means clustering (Cebeci, 2019), a soft clustering method that is similar in concept and algorithm to k-means clustering, which partitions data points within a given dataset into defined numbers of clusters based on their proximity to the nearest cluster's centroid.In fuzzy k-means, data points are not exclusively assigned to just one cluster, but rather given scores of membership to all clusters, allowing for "fuzziness" or overlap between two or more clusters.This allows for additional characterization of high-scoring cells within each cluster, cells with more ambiguous identities, and other cases that the user might be interested in, which might be informative to their specific dataset.Fuzzy k-means also assigns a final "hard" cluster assignment based on the class with the highest membership score, which can be used as input for downstream analysis.These final cluster assignments were then used for the analysis of the 2xLPS mouse dataset in this paper, unless otherwise specified.Using exploratory data analysis methods including the within sum of squares and silhouette methods (Extended Data Fig. 4-1C), we found that a clustering parameter of four yields the highest degree of within-cluster similarity and was thus the most optimal parameter to use for our example 2xLPS dataset.

Cluster characterization and analysis
Once cluster membership was defined using k-means clustering, we further explored what features describe the different clusters (Fig. 3D) and how cells belonging to each cluster visually look using the ColorByCluster feature in MicrogliaMorphology (Fig. 3E).We recommend that users always perform these steps in addition to the initial clustering optimization steps (Extended Data Fig. 4-1C) to verify that the clusters defined within their datasets are morphologically distinct and in line with expected differences in microglia morphology.We computed the average values for all 27 morphology measures, scaled across clusters, to characterize how each morphological cluster was differentially defined by the various morphology measures relative to the other clusters (Fig. 3D).Cluster 1 had the lowest branching complexity and territory span, resembling the classic ameboid shape in the original images upon visual confirmation using the ColorByCluster feature in MicrogliaMorphology. Cluster 2 had the greatest oblongness and branching inhomogeneity and resembled rod-like shapes; Cluster 3 had the highest branch lengths and density of foreground pixels in the hull with average territory span values relative to the other clusters and appeared hypertrophic; and Cluster 4 had the greatest branching complexity, territory span, and circularity and appeared ramified (Fig. 3D,E).Cluster 1 (ameboid), 2 (rod-like), and 4 (ramified) cells had relatively lower overlap in PC space with each other compared with Cluster 3 (hypertrophic) cells, which highly overlapped with Cluster 1 and Cluster 2 cells (Fig. 3C).This was expected, as hypertrophic cells represent a state between ameboid and rod-like forms on the morphological spectrum (Fig. 3C,E).Clusters were similarly described by the different morphology measures in the independent 1xLPS dataset (Extended Data Fig. 3-1C,D), further pointing to MicrogliaMorphology and MicrogliaMorphologyR as a robust means to characterize and analyze microglia morphologies.

LPS-induced changes in cluster membership and individual morphology measures
To highlight the utility of MicrogliaMorphologyR for assessing changes in cluster membership, we fit statistical models (see Materials and Methods) to test for LPS-induced morphological population shifts using the "stats_cluster.animal"function.Analysis of deviance (Type II Wald chi-square tests) on the models revealed significant contributions of cluster, treatment, and antibody interactions for all three brain regions [FC: X 2 (6, n = 6) = 17.676,Pr(>Chisq) = 0.007; HC: X 2 (6, n = 6) = 50.684,Pr(>Chisq) = 3.428 × 10 −9 ; STR: X 2 (6, n = 6) = 77.402,Pr(>Chisq) = 1.228 × 10 −14 ; Extended Data Table 4-1].Using our toolset, we were able to characterize morphological population shifts within brain regions using different microglial markers in our experimental mouse model.Within the hippocampus and striatum, LPS-induced decreases in ameboid cluster membership and increases in hypertrophic cluster membership were maintained across all three microglia antibodies, while within the frontal cortex, changes in ameboid cluster membership were only seen in Iba1 and P2ry12-stained datasets (Fig. 4A, Extended Data Table 4-2).LPS-induced changes in ramified and rod-like cluster membership were mostly maintained between Cx3cr1-and Iba1-stained datasets, compared with changes in the P2ry12-stained dataset (Fig. 4A, Extended Data Table 4-2).Antibody-specific differences were also apparent upon examination of the immunofluorescent images for each of the antibodies, where in the baseline PBS condition, P2ry12 distribution was less concentrated in the cell bodies and more spread throughout the cell branches (Fig. 4B).
To highlight the utility of MicrogliaMorphologyR for assessing changes in individual morphology measures, we fit statistical models (see Materials and Methods) to test for LPS-induced shifts in cell area as an example using the "stats_morphologymeasures.animal" function.Analysis of deviance (Type II Wald chi-square tests) on the models revealed significant contributions of treatment and antibody interactions for all three brain regions [FC: X 2 (2, n = 6) = 46.489,Pr(>Chisq) = 8.035 × 10 −11 ; HC: X 2 (2, n = 6) = 189.851,Pr(>Chisq) = 5.949 × 10 −42 ; STR: X 2 (2, n = 6) = 255.691,Pr(>Chisq) = 3.002 × 10 −56 ; Extended Data Table Within the hippocampus and striatum, we observed LPS-induced increases in area for P2ry12-stained datasets that were not present in .Taken together, our findings from the cluster-level and morphology measure-level analyses show that different antibodies produce different patterns of morphology changes but all work effectively with our toolset.

Discussion
MicrogliaMorphology and MicrogliaMorphologyR, a high-throughput pipeline to characterize microglia morphological states at a single-cell resolution Microglia exhibit a dynamic range of morphologies including, but not limited to, ramified, ameboid, rod-like, and hypertrophic forms that are highly context specific and often rapidly changing in response to local environmental cues (Savage et al., 2019;Paolicelli et al., 2022;Reddaway et al., 2023).There has been a concerted effort as a field to move away from dualistic characterization of all microglia as "resting" or "activated," which is often described in terms of morphological differences, and toward a clearer understanding and appreciation for the heterogeneous "states" of microglia that coexist in the brain in any given context (Dubbelaar et al., 2018;Paolicelli et al., 2022).In line with these efforts, there have been many recently published tools that classify and analyze microglia morphological subpopulations in an automated and high-throughput manner (York et al., 2018;Salamanca et al., 2019;Clarke et al., 2021;Leyh et al., 2021;Colombo et al., 2022;Melo et al., 2023;Reddaway et al., 2023;van Weering et al., 2023).Using our toolset, MicrogliaMorphology and MicrogliaMorphologyR, we take a data-informed approach to characterize different populations of microglia morphologies and to statistically model how membership across all morphological states dynamically changes in experimental contexts and across brain regions in an automated and high-throughput manner, which offers a great advantage over more labor-intensive morphological approaches which employ manual categorizations of cells or assessment of individual measures of morphology rather than morphological states (Torres-Platas et al., 2014;Young and Morrison, 2018).Furthermore, the ColorByCluster feature within MicrogliaMorphology and functions within MicrogliaMorphologyR facilitate comparisons of morphological measures across clusters and together provide a thorough validation for verifying cluster identities both visually and analytically compared with existing tools.We demonstrate that MicrogliaMorphology and MicrogliaMorphologyR are able to reproducibly detect both subtle and pronounced changes in microglia morphology and together provide a robust method to characterize morphological states across a wide range of experimental and disease models.While our dataset was too underpowered to quantify sex differences in morphology, our tool could also be used to explore known sex differences in microglia in various contexts such as early brain development (Sullivan and Ciernia, 2022).
Importantly, we made both MicrogliaMorphology and MicrogliaMorphologyR free and open source resources.Our toolset only relies on software that is open source and freely available to download, and all relevant materials including input images, data, and supporting code used in this study are available on the OSF website at https://osf.io/8t5c2/.We have also included all of the single-cell images that were generated for this study at the OSF link, which with further processing can provide a unique, benchmarking dataset for researchers interested in applying other approaches such as machine learning methods to classify microglia morphology.The ImageJ and R code underlying both MicrogliaMorphology and MicrogliaMorphologyR, as well as the R code for the analysis demonstrated in this paper, are all available through GitHub repositories at https://github.com/ciernialab/MicrogliaMorphologyand https://github.com/ciernialab/MicrogliaMorphologyR in the hopes that the larger research community can openly share troubleshooting tips, benefit from discussion, and continue to expand upon our work and develop our toolsets for broader use.
MicrogliaMorphology and MicrogliaMorphologyR can be applied to quantify morphology from datasets stained for 3 different microglia markers Cx3cr1, Iba1, and P2ry12 are all antibody markers that are commonly used to visualize and study various aspects of microglia including morphology.Of these three markers, P2ry12 is the most microglia specific, as both Cx3cr1 and Iba1 also label other macrophages (Paolicelli et al., 2022).While all three of these markers are considered "homeostatic" markers, their expression can change in disease-associated and inflammatory states (Kenkhuis et al., 2022) and infiltration of peripheral macrophages into the brain in diseases like Parkinson's disease and multiple sclerosis (Prinz and Priller, 2017) could complicate the analysis of morphology in datasets stained with nonspecific markers such as Cx3cr1 and Iba1.Thus, careful consideration of morphological markers should be taken depending on the experimental context in which microglia are being studied (Paolicelli et al., 2022).In our analyses of LPS-induced shifts in morphological populations and changes in individual morphology measures from sample-matched Cx3cr1, Iba1, and P2ry12-stained images, we found that P2ry12 showed unique differences in the percentage of morphological populations present across brain regions, the directionality of shifts across morphological populations, and the specific morphological features such as cell area that change with LPS administration (Fig. 4A,C).However, the objective of generating the triple-labeled datasets in this study was to provide a proof-of-concept for the application of MicrogliaMorphology and MicrogliaMorphologyR to different microglia antibodies.Thus, the preliminary findings presented here should be further explored in a dataset with sufficient sample size to be able to more comprehensively characterize and compare marker-specific effects.

Future directions
To classify microglia by their morphological characteristics in any approach, hard cutoffs are used to define where one class starts and the next begins.This type of binning of morphologies is limited in that microglia are inherently dynamic and more realistically exist along a continuum of morphological forms.To represent the dynamic nature of microglia morphology, we also demonstrate that MicrogliaMorphology and MicrogliaMorphologyR can be optionally integrated with a soft clustering approach using fuzzy k-means clustering, which is similar in concept and algorithm to k-means clustering.Soft clustering approaches such as fuzzy k-means clustering not only yield final cluster (or class) identities for every cell but also membership scores of belonging to any given cluster.This allows for additional characterization of high-scoring cells within each cluster (i.e., quintessential "rod-like", "ameboid", "hypertrophic", or "ramified" cells), cells with more ambiguous identities (e.g., a cell that is 5% rod-like, 5% ameboid, 45% hypertrophic, and 45% ramified), and other cases that the user might be interested in which might be informative for their specific dataset.Fuzzy k-means also assigns a final hard cluster assignment for every cell based on the class with the highest membership score, so the user can also use these final assignments as input for downstream analysis.While we used the hard cluster assignments for the analysis in this paper, we provide an example of using the soft clustering assignments from fuzzy k-means to analyze just the highscoring cells for each morphological class in our GitHub page for MicrogliaMorphologyR.While we used k-means clustering approaches in this study, our toolset is highly flexible and can also be integrated with other clustering approaches such as hierarchical clustering or Gaussian mixture models.
Microglia have long been known as a highly heterogeneous cell type as defined by their morphology, electrophysiological properties, transcriptomic profiles, and surface expression of immune markers (Hammond et al., 2019;Li et al., 2019;Masuda et al., 2020;Paolicelli et al., 2022).Context-specific regulation of morphology further emphasizes the need to probe microglial phenotypes from multiple angles in conjunction with morphology to gain more clarity on the relationship between microglial form and function (Dubbelaar et al., 2018;Paolicelli et al., 2022).While the majority of studies of microglia morphology have yielded observational insights into the range of forms present in various contexts, only a few (Parakalan et al., 2012;Madry et al., 2018, p. 1;Adrian et al., 2023) have actually explored how different morphological states directly contribute to microglial function in the brain.The rise of single-cell sequencing technologies has provided vast new insight into the molecular mechanisms that shape heterogeneous microglial responses and has granted us a better understanding of microglial "states" in homeostatic, developmental, and disease-relevant contexts (Hammond et al., 2019;Li et al., 2019;Masuda et al., 2020).However, transcriptomic characterization alone does not capture the diversity of changes that microglia exhibit, and we still lack a direct understanding of whether morphologically different microglia populations are transcriptomically distinct.Parakalan et al. (2012) is one such study that directly explored these relationships by identifying over 2,000 differentially expressed genes with unique sets of biological functions between ameboid and ramified microglia laser-dissected and pooled from rat brains.However, while it is agreed upon that microglia exhibit a wide range of morphological forms across various biological contexts, it is still unclear whether we can transcriptomically define the heterogeneous morphological states that exist outside of the "resting" versus "activated" morphological dichotomy and in what ways these transcriptomic signatures relate to microglial function.
The advent of spatially resolved transcriptomics and development of methods for integrating multiple data modalities has opened new avenues to explore these relationships more directly.The ability to map morphologically classified microglia back to their spatial locations in their original input images using the complementary ColorByCluster functions in MicrogliaMorphology and MicrogliaMorphologyR allows for not only the visual verification and exploration of morphological cluster identity across tissue sections but also facilitates the direct integration of spatial transcriptomics data to morphological data at a cellular resolution.Our toolset serves as a resource that can complement new tools and approaches such as spatial transcriptomics to answer questions about the relationship between microglia morphology and microglia function more directly.

Requirements of MicrogliaMorphology and MicrogliaMorphologyR
Our toolset requires the use of 2D image forms (EDF or single-plane 2D images) and depends on some user input to determine dataset-specific parameters for MicrogliaMorphology. Introductory skills coding in the R language are also necessary to be able to use MicrogliaMorphologyR, and larger computational resources may be required for the analysis of larger datasets.However, many of these requirements exist for alternative approaches for morphology analysis that have been presented as well.

Figure 1 .
Figure 1.Study overview.A, Outline of steps involved in MicrogliaMorphology and MicrogliaMorphologyR. B, Experimental mouse model used for dataset described throughout paper, made using Biorender.Example images from the dorsal hippocampus with individual microglia insets for each treatment condition.Microglia (Iba1) in yellow and DAPI nuclear stain in blue.Full size images scale bar 200 µm, insets 30 µm.

Figure 2 .
Figure 2. Comparison of 3D versus EDF versus single-plane 2D image types.A, Samples represented in principal components space and colored by image type or morphology class.Each point is either a 2D, 3D, or EDF representation of 1 of 20 different cells.B, Comparison of changes across morphological classes when cells are represented in 2D, 3D, or EDF forms.Values on plots are z-scores (centered and scaled) calculated within image type.C, Spearman's correlation of PCs 1-2 after dimensionality reduction across image types.D, Individual Pearson's correlations between image types for specific morphology features measured using AnalyzeSkeleton.E, Visual description of morphology features measured using AnalyzeSkeleton.See Extended Data Figure 2-1 and Extended Data Tables 2-1 and 2-2 for additional data.

Figure 3 .
Figure 3. Characterization of morphological clusters in 2xLPS dataset.A, Spearman's correlation matrix of 27 features measured by MicrogliaMorphology [*abs(R) ≥ 0.8; p < 0.05].B, Spearman's correlation of morphology measures to first three PCs after dimensionality reduction [*abs(R) ≥ 0.75; p < 0.05].C, Cluster classes displayed in PCs 1-2 space.D, Average values for all 27 morphology features, scaled across clusters.E, Individual cells spatially registered back to original images and visually annotated by morphological class using ColorByCluster feature.See Extended Data Figure 3-1 and Extended Data Tables 3-1-3-4 for additional data.

Figure 4 .
Figure 4. Analysis of morphological clusters and individual morphology measures across brain regions and antibody markers in 2xLPS dataset.A, LPS-induced shifts in morphological populations across antibodies and within brain regions (*p < 0.05; Bonferroni).B, Immunofluorescent images of the same microglial cells stained with Cx3cr1, Iba1, and P2ry12 in PBS and 2xLPS conditions.Scale bars, 50 µm.C, LPS-induced changes in cell area across antibodies and within brain regions (*p < 0.05; Bonferroni) See Extended Data Figure 4-1 and Extended Data Tables 4-1-4-4 for additional data.