Deformably Registering and Annotating Whole CLARITY Brains to an Atlas via Masked LDDMM
The CLARITY method renders brains optically transparent to enable high-resolution imaging in the structurally intact brain. Anatomically annotating CLARITY brains is necessary for discovering which regions contain signals of interest. Manually annotating whole-brain, terabyte CLARITY images is difficult, time-consuming, subjective, and error-prone. Automatically registering CLARITY images to a pre-annotated brain atlas offers a solution, but is difficult for several reasons. Removal of the brain from the skull and subsequent storage and processing cause variable non-rigid deformations, thus compounding inter-subject anatomical variability. Additionally, the signal in CLARITY images arises from various biochemical contrast agents which only sparsely label brain structures. This sparse labeling challenges the most commonly used registration algorithms that need to match image histogram statistics to the more densely labeled histological brain atlases. The standard method is a multiscale Mutual Information B-spline algorithm that dynamically generates an average template as an intermediate registration target. We determined that this method performs poorly when registering CLARITY brains to the Allen Institute’s Mouse Reference Atlas (ARA), because the image histogram statistics are poorly matched. Therefore, we developed a method (Mask-LDDMM) for registering CLARITY images, that automatically finds the brain boundary and learns the optimal deformation between the brain and atlas masks. Using Mask-LDDMM without an average template provided better results than the standard approach when registering CLARITY brains to the ARA. The LDDMM pipelines developed here provide a fast automated way to anatomically annotate CLARITY images; our code is available as open source software at http://NeuroData.io.
Further author information: (Send correspondence to K.S.K.)
K.S.K.: E-mail: email@example.com Kwame S. Kutten, Joshua T. Vogelstein, Nicolas Charon, Li Ye, Karl Deisseroth, Michael I. Miller, “Deformably registering and annotating whole CLARITY brains to an atlas via masked LDDMM,” Optics, Photonics and Digital Technologies for Imaging Applications IV, Peter Schelkens, Touradj Ebrahimi, Gabriel Cristóbal, Frédéric Truchetet, Pasi Saarikko, Editors, Proc. SPIE 9896, 989616 (2016).
Copyright 2016 Society of Photo Optical Instrumentation Engineers. One print or electronic copy may be made for personal use only. Systematic electronic or print reproduction and distribution, duplication of any material in this paper for a fee or for commercial purposes, or modification of the content of the paper are prohibited.
Imaging whole brains at the cellular level without disturbing their underlying structure has always been challenging. All cells are surrounded by a phospholipid bilayer which scatters light, rendering most biological tissues opaque to the naked eye. Thus it is often necessary to physically slice brains in order to use light microscopy. Sectioning tissue has two major drawbacks for researchers interested in building a whole brain connectome. First, slicing can dislocate synapses and axons necessary for tracing neuronal circuitry. Second, the inter-sectional resolution will always be much lower than the intra-sectional resolution, making neurite tracing difficult .
CLARITY avoids these problems by converting the brain into a translucent hydrogel-tissue hybrid. In the procedure, the brain is first perfused with hydrogel monomers and formaldehyde. When heated, the monomers and formaldehyde polymerize to form a molecular mesh which crosslinks amine groups of biological molecules. Since phospholipids lack amine groups, they do not crosslink with the mesh and can be eluted away with a strong detergent. The remaining hydrogel-brain hybrid is relatively translucent and permeable to fluorescent antibodies, making it amenable to labeling and interrogation by light-sheet microscopy. [2, 3] An axial slice through a CLARITY volume and magnified cutout are shown in Fig. 1.
1.2 NeuroData Cluster
CLARITY image volumes are often over 1 terabyte in size, far too large to be visualized and analyzed on a personal computer. Furthermore, extracting meaningful information from “big” images can be both time-consuming and difficult. Thus, the NeuroData cluster was created to address these challenges. As part of the project, open source software for access, visualization, and analysis of terabyte-scale images was developed. Images are stored in a multiresolution hierarchy with level 0 being the native resolution, and each subsequent level at half the previous level’s resolution. Thus the infrastructure is optimized for applying computer vision algorithms in parallel across multiple scales.
Data can be accessed through a RESTful API, a stateless interface which allows end users to download image cutouts or upload data using specific URLs. In the API, images are identified by a unique token, each of which can have one or more channels. The Connectome Annotation for Joint Analysis of Large data (CAJAL) package provides access to this API through MATLAB  while NeuroData Input/Output (ndio) provides access through Python. Images can also be visualized in a web browser using NeuroDataViz.
By annotating an entire image volume, one can draw conclusions on the texture and shape of a given brain structure. Since manual labeling is time-consuming, the most efficient annotation method is registration to a standard atlas. Spatially normalizing several subject brains into an atlas space makes it easier to determine how one brain differs from another in any given structure or location. Furthermore, this can also aide in visualization. Raw CLARITY images are often acquired in an oblique plane, making it difficult for observers to identify structures on a 2D display. Aligning the brains with an atlas solves this problem by allowing brain visualization in one of the three standard planes (axial, sagittal, and coronal). Transforming atlas labels to a subject’s space facilitates analysis of image features within brain structures. (Fig. 2).
1.4 Previous Work
Registering CLARITY images to each other or to an atlas has recently become a topic of interest. A preceding study described the development of a pipeline for registering the Allen Institute’s Mouse Reference Atlas (ARA) to images of transsynaptic viral traced brains.  The ARA is a widely used mouse brain atlas which includes Nissl-stained reference images and over 700 manually defined brain structures. A test image was acquired using serial-two photon (STP) tomography, a technique which pairs a two-photon microscope with a vibratome for automatic tissue slicing.  Since the intensity profiles of the Nissl-stained ARA and the test image differed greatly, the images were registered using their corresponding brain masks. Masks were aligned first using affine registration, followed by deformable registration using Large Deformation Diffeomorphic Metric Mapping (LDDMM), an algorithm which computes smooth invertible transforms between image volumes. Qualitative observation showed that the registration results were acceptable in most parts of the brain, although alignment of deeper structures were less accurate. In a different study, 25 CLARITY images were registered and averaged to create a single reference template using a mutual information metric with B-spline transforms. This template was then used to construct an atlas for experiments combining CLARITY with transsynaptic viral tracing. 
The Insight Segmentation and Registration Toolkit (ITK) is an open source library funded by the National Library of Medicine. It was developed by Kitware Inc. and has been used widely within the medical imaging community.  ITK’s registration framework is designed to be modular. To register a moving image to a fixed image, the user selects a metric (e.g. mean squared error) to compare the images, and the transform (e.g. affine) applied to the moving image. An optimizer (e.g. gradient descent) is used iteratively to improve the transform parameters. Additionally, the user can select the type of interpolation (e.g. nearest neighbor) for resampling the moving image.  ITK version 4 includes support for time-varying velocity field transforms often used in diffeomorphic registration algorithms such as LDDMM. This greatly facilitates the implementation of these types of algorithms in ITK.  This functionality can be used in Python through SimpleITK, an easy-to-use interface for ITK’s algorithms. 
In addition to their large size, CLARITY images present several unique challenges for image registration. In the CLARITY images of this study, a neuron’s brightness is proportional to its activity, which means that CLARITY images have a functional component. Regions which appear bright in one CLARITY brain may appear dark in another (Fig. (a)a). Registration is further complicated by brain deformation introduced in the clarifying process (Fig. (b)b) and missing data (Fig. (c)c).
The registration pipeline from the preceding study  was reimplemented using SimpleITK. This pipeline, known henceforth as the Mask-LDDMM pipeline, registered images using their masks. Additionally, an Image-LDDMM pipeline which directly registered images using intensity was developed.
2.1 Image Acquisition
12 CLARITY mouse brains (5 wild type controls and 7 behaviorally challenged) were imaged using CLARITY-Optimized Light-sheet Microscopy (COLM) (whole brain COLM imaging and data stitching performed by R. Tomer, in preparation). In brief, raw volumes were acquired in 0.585 m x 0.585 m resolution slices with a slice spacing of 5 to 8 m. Images were stored in the NeuroData cluster at 6 resolution levels, with level 0 being the full resolution and level 5 being the lowest resolution. To avoid registration complications, four CLARITY brains which were not missing any data (Control239, Challenged178, Challenged199, and Challenged188) were selected to test the pipelines.
2.2 Mask Generation
Let be the background space. Let be the template (moving) image which will be deformed to match target (fixed) image . Thus for CLARITY-ARA registrations is the the Nissl-stained ARA, are the corresponding ARA annotations, and is a CLARITY image. Since the resolution of ARA version 2 is 25 m isotropic, CLARITY images were downloaded from the NeuroData cluster and resampled to the same resolution. ARA mask was generated from by taking the union of all foreground labels. For each CLARITY image , mask was generated using the following procedure. First, was binary-thresholded to remove the background. Next, this rough mask was opened using a 50 m radius ball-shaped kernel to remove foreground grains. Finally, the mask was closed by a 125 m radius kernel to remove large foreground holes.
In the Image-LDDMM Pipeline, preprocessing consisted of two steps. In the first step, the brain masks and were applied to the images and respectively. Next, the masked was histogram-matched to the masked . In the matching procedure, 32-bin histograms were calculated for both the template and target images. The histograms were matched exactly at 8 quantile points, and by interpolation at all other intensities between these points. In the Mask-LDDMM pipeline, no preprocessing was done. Instead and were replaced by their corresponding masks, and , during registration.
Registration was then done in a three step process with rigid alignment, followed by affine alignment, and finally deformable registration. In rigid registration, parameters for the rigid transformation augmented matrix were optimized using gradient descent on the mean squared error between the transformed histogram-matched and . For affine alignment, the same optimization scheme and image metric were used to find affine matrix between and . The results of these intermediate steps, the Mask-Rigid/Mask-Affine in the Mask-LDDMM pipeline and Image-Rigid/Image-Affine in the Image-LDDMM pipeline, were stored for quantitative evaluation.
Let and . Deformable registration was done by the LDDMM algorithm which used gradient descent to minimize the objective function
where is the velocity of the flow from to and is a kernel which ensures that is sufficiently smooth. The greater is, the smoother the transform. . By convention and were used. The algorithm was implemented in ITK by building upon the library’s time-varying velocity field registration method. To minimize computational cost, velocity field was discretized into only 4 time steps. Displacement was found by integrating using a 4th order Runge-Kutta algorithm. After LDDMM the final transform and its inverse were found.
Deformed labels were resampled to a level 5 resolution and fed into the NeuroData cluster. The infrastructure automatically propagated the annotations to higher resolution levels. Thus the full resolution images with ARA labels overlaid could be visualized from a web browser in NeuroDataViz (Fig. 4).
2.6 NeuroData Registration in Python
Image processing and registration code from the pipelines were turned into a python module called NeuroData Registration (ndreg). Functions from ndreg store results as NIfTI image files to leave a record of each step for easy debugging. Internally, these functions called SimpleITK or custom binaries. The module also includes convenience functions for downloading and uploading data to the NeuroData cluster through ndio. The ndreg module is available as open source software at http://NeuroData.io.
2.7 Quantitative Evaluation
The registrations were quantitatively evaluated using aligned template-to-target mutual information, surface error, and landmark error. Let and be random variables representing the intensities of deformed template image and target image respectively. The Mutual Information (MI) between these images is
where , , and are the deformed template histogram, target histogram, and joint histogram respectively. Computing MI directly may yield unstable results. Therefore it was estimated using the Viola-Wells method, as implemented in ITK. As recomended in the doccumentation, a standard deviation of 0.4 was used to smooth the histograms after normalizing the template and target image intensities to a zero mean and unity standard deviation. Densities , , and were then estimated from 1000 foreground samples using a Gaussian distribution-based Parzen window.
Landmark-based methods were also used to evaluate the results. Specifically, landmarks were chosen and placed on the ARA and CLARITY images using the MRI Studio software suite’s DiffeoMap program (https://www.mristudio.org/). Of the landmarks, 28 were placed on the surface of the brain while the remaining 27 were placed on internal structures. The error in position of the landmark after registration is
where and are the positions of the template and target landmarks for .
The Hausdorff Distance (HD) between the deformed template and target surfaces is defined as
where and denote the surfaces of the deformed template and target respectively. Since the Hausdorff Distance is a maximum distance between surfaces, it was usually too sensitive to outlier surface points. Therefore it was more convenient to compute the median of the surface distances between all points in and
Three experiments were performed. In the first experiment, CLARITY images were registered directly to the ARA. In the next experiment, CLARITY images were registered to other CLARITY images of different conditions (Challenged to Control). In the final experiment, CLARITY images were registered to other CLARITY images of the same condition (Challenged to Challenged). Registration results from the Image-LDDMM and Mask-LDDMM pipelines were compared to a MI-BSpline pipeline described in an upcomming work by Tomer et al. The pipeline registers images using their Mutual Information under a B-Spline spatial transformation. The results of these experiments are summarized in Fig. 5-6.
Fig. (a)a and the green lines in Fig. 6 show the results of CLARITY to ARA registrations. The Mask-LDDMM pipeline consistently outperformed the MI-BSpline and Image-LDDMM pipelines in surface and landmark distances. Interestingly, the Image-LDDMM pipeline consistently yielded higher MI values than the Mask-LDDMM and MI-BSpline pipelines. Despite this, visual inspection of Image-LDDMM results in Fig. (a)a revealed that the alignments were relatively poor. This can likely be attributed to the great difference in appearance between the Nissl-stained ARA and CLARITY images.
3.2 CLARITY-CLARITY between conditions
3.3 CLARITY-CLARITY within a condition
Fig. (c)c and the magenta lines in Fig. 6 show the registration results of CLARITY to CLARITY registrations within a condition. As in the inter-condition registrations, Image-LDDMM yielded better MI values than Mask-LDDMM. But unlike the inter-condition registrations, Image-LDDMM and Mask-LDDMM had lower median landmark errors than the MI-BSpline method. Once again Mask-LDDMM, Image-LDDMM, and MI-BSpline gave comparable surface distance results.
In most cases, the Mask-LDDMM pipeline outperformed both Image-LDDMM and MI-BSpline in aligning CLARITY brains with the ARA. The MI-BSpline pipeline gave better results than Image-LDDMM in CLARITY-CLARITY transforms between conditions. Image-LDDMM outperformed MI-BSpline in CLARITY-CLARITY transforms within a condition. However, there are some limitations to these findings. Human error in landmark placement may have been a factor in these results. Furthermore, the pipelines were tested on brains without missing data. When only partial data is available, complete brain masks cannot be constructed and Mask-LDDMM should not be used.
Acknowledgements.This work was graciously supported by the Defense Advanced Research Projects Agency (DARPA) SIMPLEX program through SPAWAR contract N66001-15-C-4041 and DARPA GRAPHS N66001-14-1-4028.
- Kim, S.-Y., Chung, K., and Deisseroth, K., “Light microscopy mapping of connections in the intact brain,” Trends in Cognitive Sciences 17(12), 596–599 (2013).
- Chung, K., Wallace, J., Kim, S.-Y., Kalyanasundaram, S., Andalman, A. S., Davidson, T. J., Mirzabekov, J. J., Zalocusky, K. A., Mattis, J., Denisin, A. K., Pak, S., Bernstein, H., Charu Ramakrishnan, L. G., Gradinaru, V., and Deisseroth, K., “Structural and molecular interrogation of intact biological systems,” Nature 497, 332–337 (2013).
- Tomer, R., Ye, L., Hsueh, B., and Deisseroth, K., “Advanced CLARITY for rapid and high-resolution imaging of intact tissues,” Nature Protocols 9(7), 1682–1697 (2014).
- Burns, R., Roncal, W. G., Kleissas, D., Lillaney, K., Manavalan, P., Perlman, E., Berger, D. R., Bock, D. D., Chung, K., Grosenick, L., Kasthuri, N., Weiler, N. C., Deisseroth, K., Kazhdan, M., Lichtman, J., Reid, R. C., Smith, S. J., Szalay, A. S., Vogelstein, J. T., and Vogelstein, R. J., “The Open Connectome Project data cluster: Scalable analysis and vision for high-throughput neuroscience,” Proceedings of the 25th International Conference on Scientific and Statistical Database Management (2013).
- Roncal, W. R. G., Kleissas, D. M., Vogelstein, J. T., Manavalan, P., Lillaney, K., Pekala, M., Burns, R., Vogelstein, R. J., Priebe, C. E., Chevillet, M. A., and Hager, G. D., “An automated images-to-graphs framework for high resolution connectomics,” Frontiers in Neuroinformatics 9 (2014).
- Kutten, K. S., Eacker, S. M., Dawson, V. L., Dawson, T. M., Ratnanather, T., and Miller, M. I., “An image registration pipeline for analysis of transsynaptic tracing in mice,” in [SPIE Medical Imaging ], Proc. SPIE 9788 (2016).
- “Technical white paper: Allen reference atlas - version 2,” tech. rep., Allen Institute for Brain Science (2011).
- Ragan, T., Kadiri, L. R., Venkataraju, K. U., Bahlmann, K., Sutin, J., Taranda, J., Arganda-Carreras, I., Kim, Y., Seung, H. S., and Osten, P., “Serial two-photon tomography for automated ex vivo mouse brain imaging,” Nature Methods 9(3), 255–258 (2012).
- Menegas, W., Bergan, J. F., Ogawa, S. K., Isogai, Y., Venkataraju, K. U., Osten, P., Uchida, N., and Watabe-Uchida, M., “Dopamine neurons projecting to the posterior striatum form an anatomically distinct subclass,” Elife 4(e10032) (2015).
- Johnson, H. J., McCormick, M. M., and Ibáñez, L., The ITK Software Guide Book 1: Introduction and Development Guidelines. The Insight Software Consortium, 4 ed. (January 2016).
- Johnson, H. J., McCormick, M. M., and Ibáñez, L., The ITK Software Guide Book 2: Design and Functionality. The Insight Software Consortium, 4 ed. (July 2016).
- Avants, B. B., Tustison, N. J., Song, G., Wu, B., Stauffer, M., McCormick, M. M., Johnson, H. J., and Gee, J. C., “A unified image registration framework for ITK,” in [Biomedical Image Registration ], 7359, 266–275, Springer-Verlag (2012).
- Lowekamp, B. C., Chen, D. T., Ibáñez, L., and Blezek, D., “The design of SimpleITK,” Frontiers in Neuroinformatics 7 (2013).
- Beg, M. F., Miller, M. I., Trouvé, A., and Younes, L., “Computing large deformation metric mappings via geodesic flows of diffeomorphisms,” International Journal of Computer Vision 61(2), 139–157 (2005).