1. Introduction
2D imaging methods, such as optical microscopy, are still preferable to 3D imaging methods due to their high level of specificity and high resolution properties. Histological sections (slices) obtained through 2D imaging methods provide useful information for the diagnosis or the study of pathology. 3D volume reconstruction from these 2D slices is required in order to fully appreciate anatomical structures [1].
Typically, a 3D volume is reconstructed by registering (aligning) the 2D sections with respect to a chosen reference and stacking successive aligned sections [2]. As the acquisition processes of different 2D histological images are performed independently, slice misalignment and deformation is often unavoidable. Figure 1 shows examples of histological slices with noncohorent distortions, tears, hole and missing parts. The deformation varies from section to section and noncohorent distortions may exist in consecutive sections. Choosing an arbitrary slice as a reference slice leads to errors in 3D volume reconstruction, hence, the reference slice should be chosen properly not to contain distortions in order to achive high quality volume reconstruction [3].
Automatic registration of histological slices are necessary because manual registration using interactive alignment is nonreproducible and user dependent, therefore, it cannot be used if the number of sections is large [2, 1]. Among various automated registration methods proposed in the literature, rigid and affine registration methods can only handle global deformations and the transformation recovered from rigid registration has no clinical significance. Since histological slices change smoothly from slice to slice and the section distortions induced by the preparation process are local in nature [4, 3], accurate alignment of these slices can be achived by using elastic registration methods [5, 6, 7, 8, 9].
Although intensity differences in consecutive slices are generally assumed to be small and are often ignored in most of the existing nonrigid registration methods, it is shown that intensity variations for the same tissues can lead to large registration errors and the reconstructed 3D volume will not be smooth [3, 10, 6, 7, 4, 2, 11]. Figure 2 demonstrates that consecutive slices may have different brightness/contrast characteristics even for the same tissue regions. Similar to our previous study on MRI [10], to overcome intensity variations, we perform registration of histological images in standard intensity scale where similar intensities represent similar tissues.
To ensure that reconstructed 3D volume represents the full anatomy, one way is to superimpose the reconstructed volume onto an unsectioned reference volume, if it exists [2]. However, the problem with this approach is that a reference volume is not always available [4]. In the absence of a reference volume, qualitative evaluation of reconstructed volume is error prone and quantitative evaluation is often needed. The smoothness of the reconstructed surface/volume can be used as a quantitative evaluation method in that case [4]. Therefore, we use smoothness based measurement method, Correlation Alignment Measure (CAM), to quantify reconstructed volume.
A novel 3D reconstruction method for histological images is described in this paper. The method tackles three difficult problems in registration of histological images. First, in order to capture intensity variations between slices, we standardize histological images (Section 2
). Second, classical motion estimation based affine and locally affine registration methods are used to register images in feature space and image space respectively (Section
3). Third, we propose an automatic best reference slice (BRS) selection algorithm based on iterative assessment of image entropy and mean squre error (MSE) of the registration process in order to improve the quality of the reconstructed volume (Section 4). Quantitative evaluation of reconstructed volume is given in Section 5. The paper is concluded with future research directions and conclusion (Section 6).2. Standardization of Image Intensity Scale
Image intensity variations are not only influenced by the distribution of light sources, but also the content (different tissues) of the images as different tissues show different intensity levels. It is important to transform the image scale into the standard intensity scale so that for the same body region, intensities will be similar [3, 10, 12]. This is called the intensity standardization. As can be seen from Figure 2, some edges are clearly visible in one image, but are not visible in the other. This shows the need for intensity standardization.
Standardization is a nonlinear preprocessing technique which maps image intensity histogram (scale) into a standard intensity histogram through a training and a transformation step. In the training step, a set of images of the same body region are given as input to ”learn” histogramspecific parameters, or the landmarks. In the transformation step, any given image is standardized with the estimated histogramspecific landmarks obtained from the training step.
It has been proven that image histograms of the same body region are always of the same type, and most of the histograms of biomedical images are bimodal [13]. Since the histological rat brain images we study also produce bimodal intensity histograms, we extract histogram specific landmarks according to bimodal distribution as suggested in [12, 13].
In bimodal histograms, one of the histogram specific landmarks is the mode () representing the main foreground object in the image, as depicted in Figure 3(a).
Other histogram specific landmarks denoted by and are extracted according to the range of intensity of interest (IOI) by setting minimum and maximum percentiles values, namely and . For any image , we will consider histogram specific landmarks as:
(1) 
where we can increase the number of histogram specific landmarks by setting more percentiles values for the main foreground [12, 13].
In the training step, for each image, the landmarks obtained from the histogram are mapped into the standard scale by mapping intensities from to where and are minimum and maximum intensities on the standard scale respectively. The formula for mapping to is the following [12].
(2) 
Figure 3(b) shows the two linear mappings. The first from to and the second from to . Overall mapping, , from to can be summarized as follows:
(3) 
where converts any number y to the closest integer Y such that Y or . Further details and proofs can be found in [13].
Figure 4 shows slices before and after standardization. The first row shows the original data displayed using the default window setting. The second row shows the same slices after standardization displayed using the ”standard” window settings with the parameters defined in Section 5.
3. Registration in Feature Space
Registration of histological slices requires serial registration procedure which is a combination of transformations. Let be the transformation that warps the source image to the target image . The transformation is computed serially as follows:
(4) 
where represents composition.
3.1. Locally Affine Nonlinear Transformation
Since consecutive slices are not exactly the same, rather slices vary smoothly, locally affine globally smooth (LAGS), registration algorithm fits well to the problem. LAGS^{1}^{1}1also known as elastic registration registration algorithm uses 8affine parameters to fully represents changes between 2D images. Two of these affine parameters are needed to capture local brightness and contrast patterns [8] and six affine parameters are used to capture local deformations for 2D images. As described in our previous study [10], there is no need to use these 2 affine parameters because the standardization procedure is used to remove intensity variations among the same tissue types. Readers are encouraged to read [8, 10] to understand theory of LAGS and the modified algorithm which takes into account the standardization procedure.
Choice of the feature space plays a vital role in image registration especially if the similarity metric bases on the optimization function independent of spatial information (i.e mutual information). Since these kind of registration methods do not take into consideration the spatial information of pixel/voxel intensity distribution/variation, the optimization algorithm may get stuck in local maximum resulting in misalignment. Defining a feature space capturing variations of graylevel characteristics will overcome the drawbacks of intensity based approaches. To align the images globally, we used a particular feature space which represents an image by continuous variables, called edgeness, and describes the intensity variance of a predefined region over the image
[14, 3].3.2. Notation and Formulation
We represent an image (section/slice) by a pair where is a twodimensional (2D) array of pixels and is intensity function defined on , assigning an integer intensity value for each pixel . A particular feature space that allows the representation of an image by continuous variables is called the edgeness space [14], which describes the intensity variance of local regions in an image. At image coordinate , the edgeness is represented by
(5) 
where represents a fixed radius. Note that this is not just to determine whether a specific voxel/pixel is edge or not [14]. Instead, within a specified radius , the image feature content is forced to stay beyond a variation level which prevents the registration process from getting stuck in local maxima. Edgeness maps for various slices are shown in Figure 5.

4. Automatic Best Reference Slice Selection
The quality of the reconstructed 3D volume mostly depends on the choice of the reference slice. The reference slice is used as a target image and all the remaining slices are being considered as source images to be registered onto the target image. If the reference slice is distorted or noisy, reconstructed 3D volume will not be optimal. Once the reference slice is identified as target image, registration based fusion methodology can be applied for reconstruction [2].
4.1. Metrics
Selecting best reference slice can be based on high confidence image features such as MSE, entropy, edge, texture, color, intensity histograms, etc.

MSE:
In the case of distortions, structural discontinuity is not minimum even for the consecutive slices. When affine registration is performed for global alignment of images, the optimization procedure tries to minimize MSE between images but due to distortion, it will not reach low MSE values. Furthermore, it is also known that with small SNR values, alignment is difficult, leading to high registration errors. Therefore, MSE can be used as a tool for checking whether the slices are distorted or not. While high MSE values indicate most probably distorted and noisy slices, low MSE values indicate strong similarity between consecutive images.

Edge: In feature space, we emphasise edgeness features of an image by mapping image space into the feature space where edgeness parameters hold both edge information and spatial variations of pixel intensities over all regions in the image. Therefore, we assume that MSE between any image pair already includes high confidence information related to edges.

Contrast/Brightness: Contrast/brightness patterns also play important role in image contents. Since standardization method has been used to correct intensity variations, intensity for the same tissue is the same for all images.

Entropy: Entropy is another measure often used to characterise the information content of a data source. It has been used as a metric for image registration in the form of mutual information. Large mutual information between images implies high similarity and vice versa.
A common way to define entropy for the image is:
(6) 
where
is the probability density function for intensity value
and is set of all possible gray values in . It is certain that if the amount of information in the image is large, the entropy value will be high. From this point, we propose the combination of entropy information with MSE as a high confidence feature measure to automatically identify BRS for 3D volume reconstruction.4.2. Best Reference Slice (BRS)
To find the BRS according to the feature sets defined in Section 4.1, an iterative framework is developed to compare slices that are not consecutive and to avoid errors arising from noisy and distorted images. The framework divides slices into subvolumes. Let be total number of slices, the subvolumes where and the number of slices in subvolume such that: .
Since MSE is inversely related to the similarity of images and entropy information is directly related to the content of the image, we define:
(7) 
where is feature space based mean square error after registering the image with image and
(8) 
for subvolume .
5. Experiments and Results
5.1. Rat Brain Data
We have registered a stack of 350 Nisslstained slices acquired by cyrosectioning coronally an adult mouse brain with a resolution of 590x520 pixels at a resolution of 15m and 24bit color format [15].
5.2. Standardization Parameters
Based on the experiments in [12, 13], minimum and maximum percentile values are set to and respectively. In the standard scale, and are set to and . Examples of intensity mapping for a couple of images, source and target, can be found in [10] where it can be seen that after standardization, histograms are more similar in shape and location. It means that intensities have tissuespecific meaning after standardization [3, 10, 12, 13].
5.3. Implementation
Briefly, registration is performed initially for slices in each subvolume separately (See Section 4.2). Three kinds of registration are performed in the reconstruction process: rigid, affine and LAGS. MSEs are calculated according to affine registration in edgeness space and is used to select BRSs for each subvolume. Affine registration is performed in a serial manner combining transformation functions. Then, LAGS registration is performed to capture local deformations in each subvolume with respect to the chosen reference. Once LAGS registration processes are completed, subvolumes are registered to each other in a rigid manner.
In order to achieve low computational cost and accelerate the registration process, coarsetofine multiresolution framework is used. Registration in finer level is performed with the result of the previous level as initial condition. This process continues until the finest level is reached. To ensure a more accurate solution, we perform standardization after each warping/interpolation. Either small or large, intensity changes caused by the interpolation are captured by standardization
[10].5.4. Evaluations
Quantitative evaluation of the results of the reconstruction process is often difficult. It has been shown in [16] that an ideal measure of the quality of the reconstruction is the smoothness of the reconstructed surfaces. In that work, they propose a new measure based on evaluation of smoothness of the reconstructed volume called Correspondence Alignment Measure (CAM).
The CAM measure relies on the assumption that if a point is perfectly aligned, it lies midway between its corresponding points on neighbors’ sections. To compute the CAM measure for a given image, first of all, corresponding points for specifed control points in the image are identified. The associated confidence values in two adjacent images are then calculated. If the confidence is greater than a predefined threshold
, square root of the summation of the deformation vectors are added to the cumulative sum. Finally, the cumulative sum is normalized by the number of pixels which have contributed. Note that CAM gives one value for each image, therefore, mean or standard deviation of CAM values of serial images are needed to compare reconstructions. Reconstucted volume is smooth if the mean or the standard deviation of CAM measures are low and vice versa.
Summary of the changes in mean and standard deviation in CAM values is given by Table 1. The values in Table 1 are obtained by considering the worst case which uses all the slices instead of just a few slices from the middle of the stack as defined in [16], and is set to . Even for the worst case, CAM values indicate that a smooth volume is constructed with the proposed framework. While mean values dropped by % and %, the standard deviation values dropped by % and % for affine and locally affine registered stacks respectively with respect to rigid registered stacks. Figure 6 shows CAM values for each section. Registered stacks has lower CAM values which means that smoother reconstructed surface/volume is obtained.
–  Rigid Reg.  Affine Reg.  LAGS Reg. 

Mean  55.911  51.832  45.461 
Std  12.223  9.232  8.833 
6. Conclusion
In this study, we present a novel method for reconstructing 3D rat brain volumes from 2D histological images. The framework bases on three fundamental premises. (1) All histological images must be standardized for accurate registration leading to 3D volume reconstruction. (2) For accurate and succesful registrations in consecutive slices, a reliable feature space must be taken into account. (3) For automatic 3D volume reconstruction, the reference slice must be chosen properly by avoiding slices with high noise, distortions and other factors. To validate the reconstructed volume, the smoothness of the volume is considered. Experimental results indicate that the reconstructed volume is highly accurate.
7. Acknowledgements
This research is funded by the European Commission Fp6 Marie Curie Action Programme (MESTCT2005021170) under the CMIAG (Collaborative Medical Image Analysis on Grid) project.
We wish to thank Prof. Jayaram K. Udupa, Image Processing Group of University of Pennsylvania, for his comments and suggestions.
References
 [1] S. Ourselin, A. Roche, G. Subsol, X. Pennec, and N. Ayache, “Reconstructing a 3d structure from serial histological sections,” Image and Vision Computing 19, pp. 25–31, September 2000.

[2]
G. Malandain, E. Bardinet, K. Nelissen, and W. Vanduffel, “Fusion of autoradiographs with an MR volume using 2D and 3D linear transformations,”
NeuroImage 23, pp. 111–127, September 2004.  [3] U. Bağcı, “Fundamental issues of registration: Applications for change analysis in health and disease,” CMIAG Annual Report , 2007.
 [4] T. Ju, J. Warren, J. Carson, M. Bello, I. Kakadiaris, W. Chiu, C. Thaller, and G. Eichele, “3d volume reconstruction of a mouse brain from histological sections using warp filtering,” J. Neurosci. Methods (156), pp. 84–110, 2006.
 [5] J. C. Gee, D. Haynor, L. Briquer, and R. Bajcsy, “Advances in elastic matching theory and its implementation,” in CVRMed, pp. 63–72, 1997.
 [6] C. Davatzikos, “Spatial transformation and registration of brain images using elastically deformable models,” Computer Vision and Image Understanding: CVIU 66(2), pp. 207–222, 1997.
 [7] H. Lester and S. Arridge, “A survey of hierarchical nonlinear medical image registration,” Pattern Recognition (32), pp. 129–149, 1999.
 [8] S. Periaswamy and H. Farid, “Elastic registration in the presence of intensity variations,” IEEE Transactions on Medical Imaging 22, pp. 865–874, July 2003.
 [9] S. Wirtz, B. Fischer, J. Modersitzki, and O. Schmitt, “Superfast elastic registration of histologic images of a whole rat brain for 3d reconstruction,” in Medical Imaging 2004: Image Processing. Edited by Fitzpatrick, J. Michael; Sonka, Milan. Proceedings of the SPIE, Volume 5370, pp. 328334 (2004)., 5370, pp. 328–334, May 2004.
 [10] U. Bağcı and L. Bai, “Multiresolution elastic medical image registration in standard intenstiy scale,” in 20th Brazilian Symposium on Computer Graphics and Image Processing (SIBGRAPI), October 2007.
 [11] S. Ourselin, C. Sattonnet, A. Roche, and G. Subsol, “Automatic alignement of histological sections for 3D reconstruction and analysis,” Analytical Cellular Pathology 18(3), p. 123, 1999.
 [12] Y. Ge, J. Udupa, L. Nyul, L. Wei, and R. Grossman, “Numerical tissue characterization in ms via standardization of the mr image intensity scale,” Journal of Magnetic Resonance Imaging 12, pp. 715–721, November 2000.
 [13] G. N. Laszlo and J. Udupa, “On standardizing the mr image intensity scale,” Magnetic Resonance in Medicine 42(6), pp. 1072 – 1081, 1999.
 [14] T. Butz and J.P. Thiran, “Affine registration with feature space mutual information,” in MICCAI ’01, pp. 549–556, (London, UK), 2001.
 [15] A. MacKenzieGraham, E. Lee, I. Dinov, M. Bota, D. Shattuck, S. Ruffins, H. Yuan, F. Konstantinidis, A. Pitiot, Y. Ding, G. Hu, R. Jacobs, and A. Toga, “A multimodal, multidimensional atlas of the c57bl/ 6j mouse brain,” J.Anat (204), pp. 93–102, 2004.
 [16] E. Guest and R. Baldock, “Automatic reconstruction of serial sections using the finite element method,” BioImaging 3, pp. 154–167, May 2001.
Comments
There are no comments yet.