Quantitative 3-D morphometric analysis of individual dendritic spines

The observation and analysis of dendritic spines morphological changes poses a major challenge in neuroscience studies. The alterations of their density and/or morphology are indicators of the cellular processes involved in neural plasticity underlying learning and memory, and are symptomatic in neuropsychiatric disorders. Despite ongoing intense investigations in imaging approaches, the relationship between changes in spine morphology and synaptic function is still unknown. The existing quantitative analyses are difficult to perform and require extensive user intervention. Here, we propose a new method for (1) the three-dimensional (3-D) segmentation of dendritic spines using a multi-scale opening approach and (2) define 3-D morphological attributes of individual spines for the effective assessment of their structural plasticity. The method was validated using confocal light microscopy images of dendritic spines from dissociated hippocampal cultures and brain slices (1) to evaluate accuracy relative to manually labeled ground-truth annotations and relative to the state-of-the-art Imaris tool, (2) to analyze reproducibility of user-independence of the segmentation method, and (3) to quantitatively analyze morphological changes in individual spines before and after chemically induced long-term potentiation. The method was monitored and used to precisely describe the morphology of individual spines in real-time using consecutive images of the same dendritic fragment.

microscopy technique, difficulties in dendritic spine segmentation from dendrites, and the identification of true spine boundaries [22][23][24] pose challenges to the accurate quantitative analysis of spines in contemporary neurobiology. Progress in imaging technologies allowed the researchers to acquire information about dendritic spines from images obtained both in 2-D and 3-D. Thus, it is possible to examine the complete biological context including interactions with neighboring cells, connectivity, and chemical and protein composition of dendritic spines. The final piece needed to complete the puzzle is a robust and unbiased image analysis tool for quantitative assessment of the information stored in dendritic spines morphology 25,26 .
There are predominantly two kinds of approaches when dealing with the analysis of dendritic spines changes. Some methods utilize 2-D MIP (Maximum Intensity Projection) images of dendritic spines providing structural details of individual spines. However, accurate analysis of dendritic spines based on 2-D MIP images is nearly impossible. Therefore, in the present work, we focused on analyzing dendritic spine images based on 3-D volume that were generated from the confocal image stack (see Fig. 1).
A few previous studies addressed the issue of individual spine morphometry. Imaris software 27 is a commercially available tool for the four-dimensional (4-D) analysis of dendritic spines. Although Imaris software is good for analyzing the overall spine population, it fails to accurately model the three-dimensional (3-D) morphology of individual spines. Swanger et al. 28 also reported an automated method for the 4-D analysis of dendritic spine morphology. However, they used the same Imaris pipeline to develop the automated tool. Both methods generally fail to assess individual spine plasticity.
Among the existing 3-D approaches dealing with morphological analysis of the structure of dendritic spines, Janoos et al. proposed a 3-D reconstruction based on skeletonization. In this approach, first the center of the neuron is extracted, then the longest line is considered as the backbone of the dendrite and the shorter lines are considered as the centerlines of dendritic spines 29 . Such a method is time consuming to get the centerlines from the 3-D volume with complex shapes. Additionally, even when images are at the limits of confocal laser scanning microscopy resolution the precision of segmentation is limited due to skeletonization method and by the quantization errors 30 .
A 3-D neuron analysis approach was also proposed by Rodriguez et al. 31 . A 3-D reconstruction algorithm uses the Rayburst diameter 31,32 , where Rayburst is defined as casting a multidirectional core of rays from an interior point to the surface of a solid, allowing quantification of anisotropic and irregularly shaped 3D structures. The Rayburst diameter in each layer of a spine is calculated and the head and neck of each spine is defined according to the distribution of the diameters. The aforementioned method detects and classifies the spines in an efficient way, however it may be not accurate in segmentation of very complex structures. The spine shape is defined only by the head to neck ratio and a global threshold of this parameter is not adaptive when performing an analysis of different kind of images 30 .
The latest morphological analysis algorithm for dendritic spines was introduced by Shi et al. 30 and it is based on a semi-supervised learning (SSL) approach. In this framework, first the dendrite backbone is tracked on a 2-D plane then all the dendrites' surface with meshes is reconstructed in a 3-D. Next, dendritic spines segmentation based on wavelet transform is performed. The segmenting positions are located where the wavelet response on a spine section quickly changes. Features of spines such as the head and neck diameter, spine length, volume, etc. are extracted after the spines segmentation. In the last part of analysis, a small portion of detected spines are selected by a neurobiology expert, and chosen as the training set for classification. The labels of the rest of the spines can be calculated after the training process and all the detected spines can be classified by the learning framework 30 . However the accuracy and performance of the above methods strongly depend on the size of the training dataset and the features included in training vectors. Here, we propose the fast and accurate methodology which allows studying the morphology of individual spines. Since our method is not based on machine learning approach it does not require the training data set selected by a neurobiologist. We used basic mathematical notations to define different key spine compartments (e.g., spine head and spine neck) and experimentally verified that the quantitative analysis of the newly defined spine attributes accurately modeled spine plasticity. The approach that we developed allows the user to mark specific dendritic spines, segment the spines as 3-D volumes, and extract relevant morphometric features with high accuracy and minimal user intervention.

Theory and Definitions
Basic Definitions and Notations. A 3-D cubic grid, or simply a cubic grid, is represented by | Z Z 3 , the set of integers. A grid point, often referred to as a point or voxel, is an element of Z 3 and represented by a triplet of integer coordinates. Standard 26-adjacency 33 is used here, and two adjacent voxels are often referred to as neighbors of each other. The set of 26-neighbors of a voxel p, excluding itself, is denoted by N(p).
An object  is a fuzzy subset is the membership function. The support O of an object  is the set of all voxels with non-zero membership (i.e., is the background). We use a calibrated capital letter to denote a fuzzy subset, whereas a regular capital letter is used to denote its support, which is a binary set. Images are always acquired with a finite field of view. Thus, we will assume that an object always has a bounded support.
Let S denote a set of voxels. A path π in S from ∈ p S to ∈ q S is a sequence 〈 = =  p p p p q , , , l 1 2 〉 of voxels in S such that every two successive voxels are adjacent. A link is a path 〈p, q〉 of exactly two adjacent voxels. The length of a path π = 〈 〉  p p p , , , l 1 2 in a fuzzy object , denoted as π П ( )  , is the sum of lengths of all links along the path: where − p q is the Euclidean distance between p, q. The fuzzy distance 34,35 between two voxels ∈ p q Z , 3 in an object , denoted by  ω p q ( , ), is the length of one of the shortest paths from p to q: is the set of all paths from p to q. The fuzzy distance transform (FDT) of an object  is an image is the set of positive real numbers, including zero, that is the fuzzy distance from the background:

q O
With regard to spine morphology, the challenges are (1) separating the fuzzy objects  Spine from  Dendrite , which are fused at various unknown locations and scales, and (2) morphologically defining spine compartments in the segmented fuzzy objects  Spine .
The first challenge is addressed using the MSO algorithm 36 in two sequential steps.
Step 1: segmentation of the combined region   ∪ Spine D endrite from the background: Step 2: separation of  Spine and  Dendrite . The first step may be trivially achieved using simple thresholding 37 and connectivity analysis 38 . Let  be the fuzzy segmentation of the combined region that is obtained in Step 1. All subsequent analyses will be confined to the support O of , which will be the "effective image space. " Let → I O I I : [ , ] min max be the image intensity function over O. In the second step, segmentation is modeled as the opening of two fuzzy objects that are mutually fused at different unknown regions and scales in the shared intensity space, I. Here, the main challenge is to determine the local size of the suitable morphological operator that can separate multiple small mutually disconnected structures (e.g., spines) from a large connected structure (e.g., dendrite). The developed MSO algorithm combines FDT 39 and fuzzy connectivity 40 to iteratively open the two objects in I. 36 is to gradually erode the assembly of two fused objects until those two objects become mutually disconnected, thus creating two separate objects. The first iteration starts with two sets of seed voxels, S Spine and S Dendrite , and a set of common separators, S S . The initial FDT map Ω Spine,0 for the first object is computed from O, except that the voxels in ∪ S S S Dendrite are added to the background. The FDT map Ω Dendrite,0 for the other object is computed similarly. It is reasonable to assume that the sets S Spine , S Dendrite , and S S are mutually exclusive.

Multi-scale Opening Algorithm. The basic idea of the MSO algorithm
Let us now consider the coupling of two objects, where a dendritic segment (green) and the spines (red), with significant intensity overlap (illustrated in Fig. 2), are fused with each other at different unknown locations and scales. Let μ Dendrite and μ Spine denote the dendrite and spine membership functions, defined as the following:  min max is the image intensity function over O. I Spine , and I Dendrite is the representative spine and dendrite intensities that define the respective transition between pure and shared intensity bands (see Fig. 2 Spine be the set of voxels that fall inside the pure intensity band for dendrite and spine respectively. Thus, the set of voxels that fall within the shared intensity band is Spine . A fuzzy representation of the composite object may be obtained by taking the fuzzy union of the two membership functions that are shown in Equations 4 and 5. The iterative approach of the multi-scale opening of two structures takes several iterations to grow the path-continuity of an object, starting from its seed voxels (commonly added in large-scale regions), to a peripheral location with fine-scale details.
Note that after the iterative propagation of the MSO algorithm, the dendrite region is segmented as a single connected component. O Spine represents one or more disjointed spine regions (R i ), such that , where K is the total number of disjointed spine segments in O Spine , and each such segmented spine region R i contains at least one spine seed ∈ p S Spine .
Morphological definitions for the spine regions. Once the spines are segmented from the dendrite, the challenge is to assess the morphological attributes accurately. The morphology of the dendritic spines is complex and difficult to quantify. In most cases, they are described by simple parameters, such as length and head width 41 , in 2-D MIP images of the z-stacks that are acquired from the confocal images of dissociated hippocampal cultures. Previously, we developed a convolution kernel-based approach for the segmentation of spines from 2-D MIP images 42 . However, high-resolution confocal microscopy and effective 3-D rendering allow the visualization of complex structures in greater detail.
In the present study, we defined several key morphological features of 3-D dendritic spines for plasticity analysis. Specifically, we defined four key spine features that are related to the base and head of a spine using standard notations of digital topology and geometry 43,44 .
where ⋅ is the cardinality of a set). The head and tip of a spine are defined using the FDT map 43,44 Here, = l 2 is used to avoid noisy local maxima.
Definition 3. The center of the head CH i of a spine R i is the locally deepest point in the spine. In a situation where multiple locally deepest points satisfy the farthest distance criterion, their centroid is used.
In a situation where multiple points of R i satisfy the farthest distance criterion, their centroid is used.
Note that CBP i , CH i , and T i play key roles in estimating spine attributes, such as length of the spine, neck-length, neck-width, head-width, etc., for each individual spine R i . To estimate these features, we further extended the   Figure 3 shows an illustration of these key spine attributes with respect to a segmented spine.
We now estimate the neck length   points are usually concentrated around the spine-head regions. We utilize this attribute to identify the Filopodia type of spines. Then, to classify the Mushroom type spines we estimate the base-to-head distance BH ( ) i of a spine. We have observed that the Mushroom spines have short base-to-head distance in comparison to the overall length of spine. Finally, the ratio BH L / i i decisively classifies the Mushroom spines from the Spine-head Protrusions. Spine classification results on different dendritic segments are shown in Fig. 4. Apart from the four spine categories we have also observed some branched or bifurcated spines (see Fig. 4E), which could not be decisively classified into any of the four classes. Although they are few in numbers, we have placed them in a separate category for future analysis.

Experimental Results
The methodology that we developed is useful in a variety of applications that involve the accurate volumetric assessment of spine plasticity. Two specific challenges are involved in this process: (1) accurate 3-D segmentation of individual spines from the dendritic segment and (2) quantitative analysis of individual spine morphology for the effective assessment of structural changes in dendritic spines. Specifically, we analyzed spine plasticity using the 3-D morphological features that are presented in the previous section. The confocal light microscopy images of dendritic spines from dissociated hippocampal cultures were used for (1) the analysis of accuracy relative to ground-truth annotations that were generated by experimental biologists and the available state-of-the-art Imaris tool, 2) the analysis of reproducibility of user-independence of the segmentation results, and (3) the quantitative analysis of morphological changes in spines. In the following section, we describe the image acquisition protocol and image sets that were used for experimental validation. We then present the accuracy analysis of the currently developed method relative to ground-truth annotations and comparisons with the state-of-the-art Imaris tool. We also present the reproducibility analysis using three mutually blinded experts. Finally, we discuss structural changes in spines.

Description of the datasets. This study was carried out in accordance with the Ethical Committee on
Animal Research of the Nencki Institute, based on the Polish Act on Animal Welfare and other national laws that are in full agreement with EU directive on animal experimentation. All effort was made to minimize animal suffering. Dissociated hippocampal cultures from postnatal day 0 Wistar rats were prepared as described in Basu et al. 42 . Cultured hippocampal neurons were transfected 14 days in vitro (DIV) with Syn-GFP plasmid to visualize neuronal morphology. Live-cell imaging was performed on 20-22 DIV. Dendritic segments that were decorated with dendritic spines were imaged at time 0, before stimulation, and then cLTP was induced by bath application of a mixture of 50 μM forskolin, 50 μM picrotoxin, and 0.1 μM rolipram (each dissolved in dimethylsulfoxide [DMSO]) in maintenance media. Dendritic segments were imaged 10 min and 40 min after cLTP induction. Images were acquired using a Carl Zeiss LSM780 confocal microscope with a C-Apochromat 40×/1.2 NA water immersion objective using a 488 nm wavelength argon laser at 3% transmission and 70 nm/pixel resolution. A series of z-stacks were acquired at 0.2 μm steps.
In the first dataset, three different neurons from rat dissociated hippocampal cultures were imaged using a confocal light microscope, before and after cLTP induction. All of the images were captured three times: at baseline (before cLTP) and 10 and 40 min after cLTP induction. In the second dataset, three different neurons from rat dissociated hippocampal cultures were similarly imaged at baseline and 10 and 40 min after mock cLTP induction (i.e., only the solvent, DMSO, was used). During image pre-processing, we took the confocal z-stack and performed Gaussian de-noising on the 3-D image stack. The pre-processed images at time 0 are labeled as T0, and the images that were captured at 10 and 40 min are labeled as T10 and T40, respectively.
The preparation of the brain slices was based on Magnowska et al. 45 . For the visualization of changes in the shape of dendritic spines, of 1,1′-Dioctadecyl-3,3,3′,3′-Tetramethylindocarbocyanine Perchlorate (DiI) staining in stressed and control mice were performed. The mice were anesthetized and the transcardinal perfusion with 1,5% paraformaldehyde was performed. Then the brains were dissected and sliced using vibratome. Slices (140 μm thick) that contained the different brain structures recovered for at least 1.5 h at RT. Random dendrite labeling was performed using 1.6 μm tungsten particles (Bio-Rad, Hercules, CA, USA) coated with propelled lipophilic fluorescent dye (DiI; Invitrogen) that were delivered to the cells by gene gun (Bio-Rad) bombardment. Images of dendrites were acquired under 561 nm fluorescent illumination using a confocal microscope (63×objective, 1.4 NA) at a pixel resolution of 1024 × 1024 with a 3.43 zoom, resulting in a 0.07 μm pixel size. A series of z-stacks were acquired at 0.2 μm steps. Accuracy Analysis. For accuracy analysis, we first generated ground-truth spine segmentation results for all of the neuronal images by manually labeling ideal spine regions of a sample spine population with the help of experimental biologists using the open-source image processing software Fiji 46 and ITKSnap 47 . Although both Fiji and ITKSnap are general-purpose image analysis tools, neither of them are capable of extracting morphological attributes that are specific to a dendritic spine.
For the quantitative analysis of spine segmentation accuracy, we considered the most generic features, such as the volume and length of a spine. Comparative assessment was performed with regard to the ground-truth annotations that were performed by experts in this domain. The agreement in the estimated feature values and error difference over a sample spine population are illustrated in Fig. 5A. We calculated Pearson's correlation coefficients to assess the mean agreement of the estimated feature values between the currently developed method and the ground-truth annotations. The respective Pearson's correlation coefficients for volume and length, estimated over a sample spine population with manually annotated ground-truth results, were 0.89 and 0.82, respectively. Bland-Altman plots (difference plots) between the estimated values and ground-truth annotations for the two features are shown in Fig. 5C. This is a method of data plotting that is used to analyze agreement between two different data series. Most of the data samples fell within the µ σ ± .
1 5 range, representing a significant correlation among the segmented results and ground-truth annotations.
For qualitative comparisons of the 3-D segmentation results, we considered the state-of-the-art Imaris tool 27 that was applied over the same set of images that were considered in this study. Figure 6 shows the qualitative comparisons of our segmentation results with Imaris. Note that Imaris is a model-based segmentation and feature-extraction tool (outside in) that can lead to the corruption of data. Our proposed method utilizes segmentation-based modeling and feature extraction methodology (inside out). Therefore, the segmentation results and features that are derived are more accurate and robust when considering user variations. Moreover, the spines that are segmented using the currently developed software can be automatically segmented into one of four categories: stubby, filopodia, mushroom, spine head protrusions (Fig. 4) which represent morphological categorization scheme used in EM literature [48][49][50] . However, the proposed algorithm do not classify the spines as thin because this is an intermediate form of spine 51 and the mathematical definition of such spine topology does not exist.
Notably, unlike other approaches, our method distinctively identifies spine-head protrusions (Fig. 4D). Moreover, there is another type of spines that we have noted during the segmentation process, namely branched spines (see Fig. 4E). This type of spines are rare and therefore we do not classify them automatically using our present approach, but propose to include in our future work.

Reproducibility Analysis.
In this section, we present multi-user reproducibility of the currently developed method to assess the reliability and robustness of the segmentation methodology. The reproducibility analysis was performed on a sample spine population with three mutually blinded, independent experimental biologists.  Figure 7 shows the quantitative and qualitative analysis of volume, length and head width estimations over the sample spine population using three mutually blinded users. The percent standard deviations relative to the mean feature estimation for the three independent users were ± ± 11%, 17% and ±6% for volume, length and head-width respectively.
Analysis of Spine Plasticity. For detailed qualitative and quantitative analyses of dendritic spines, we assessed the segmentation results of the currently developed methodology over the same dendritic segments that were obtained from different sets of the T0, T10, and T40 images. Figure 8 shows the qualitative assessment of spine plasticity over a sample dendritic segment at the three time points. We present detailed 3-D morphological changes in dendritic spines over time (i.e., before and after cLTP). To investigate the relative changes in morphology of the segmented spines one important feature was considered from the overall experiment. Figure 9A shows the average quantitative changes in dendritic spine volume 10 and 40 min after cLTP induction.
We also show a three-class (stubby, filopodia and mushroom) spine plasticity analysis over three time points (see Fig. 9B). It shows that the 4% the of stubby spines changes to filopodia and 2% changes to mushroom, 8% of the filopodia changes to stubby and 4% changes to mushroom, and 12% of the mushroom changes to stubby from T0 to T10 after cLTP induction. Within the time span of T10 to T40, 3% of the stubby spine changes to filopodia and 9% changes to mushroom, 4% of the filopodia changes to stubby, 10% of the mushroom changes to stubby and 2% changes to filopodia.
We also investigated the changes in dendritic spines volume from brain slices (Fig. 9C) Sample segmentation results on the images acquired from brain slices are shown in Fig. 10.
We do not provide statistics since our aim is not to assess any particular biological problem, but we aim to show that our method is able to detect differences in spine structure between certain experimental conditions.

Conclusion and Discussion
There are two critical features related to neuronal signal transmission -connectivity that is related to dendritic spine number and the synaptic strength expressed by dendritic spines structural changes. In the present study we focused on the accurate 3-D segmentation and quantitative assessment of dendritic spine structure. We do not provide the method for dendritic spine density assessment since achieving accurate spine density requires several thousands of spines analyzed, thus it is extremely time consuming to perform with the 3D approach. Analysis of dendritic spine number can be done manually or with the use of 2-D software (e.g. Cheng et al. 52 ), and such approach is impossible when user wants to get both "real" 3-D spine shape and density.
The segmentation method that we developed was validated for individual spines using real-time experiments and consecutive images of the same dendritic fragment. Our results are consistent with other studies that reported spine head growth upon cLTP induction 45,53 . The plasticity analysis was performed 10 and 40 min after cLTP induction relative to baseline images. The currently developed method was also able to quantitatively assess changes in volume upon stimulation that were consistent with previous studies that reported spine head growth upon cLTP.
Accuracy was evaluated relative to manually labeled ground-truth annotations and relative to the stateof-the-art Imaris tool. We have shown, that the accuracy of the analysis done using Imaris is lower and does not reproduce the exact spine structure and it also fails to capture detailed morphology of individual spines. Imaris program fits only few standard shapes whereas their real multiplicity is much greater. In contrast our method allows to determine the real dendritic spine shape.
To assess the reproducibility of the segmentation results, three blinded experts separately assessed the efficacy of the methodology. The existence of experimental limitations (staining efficiency/quality, image resolution, and detection of fluorescence signal) impedes identification of true spine boundaries whereas our software provides the possibility to minimize user bias. Although this experiment exclusively used confocal light microscopy images of dendritic spines, the method may be extended in the future for use with other super-resolution imaging techniques, such as photo-activated localization microscopy. The present 3-D segmentation method may be used for different experimental protocols that study the structural dynamics of dendritic spines in vitro and in vivo under various physiological and pathological conditions. Moreover, we believe that this will be a research tool that enables the detection of even subtle changes in the 3-D dendritic spine structure. Such methodological advances in spine morphological studies will allow more precise analyses and better interpretations of biological data regarding structural plasticity. The present method facilitates accurate, unbiased spine segmentation results and can significantly improve the way we study dendritic spine plasticity. However, it strictly depends on the user knowledge and experimental design which spines are going to be incorporated into the analysis. Data Availability. The datasets generated during and/or analyzed during the current study are available from the corresponding authors on reasonable request.