QUANTIFICATION OF THE 3D MORPHOLOGY OF THE BONE CELL NETWORK FROM SYNCHROTRON MICRO-CT IMAGES PEI DONG ,1,2, ALEXANDRA PACUREANU4, MARIA A. ZULUAGA1,2,3, CECILE OLIVIER1,2, QUENTIN GRIMAL5, FRANÇOISE PEYRIN1,2 1.CREATIS, Université de Lyon ; CNRS UMR 5220; INSERM U1044; Université Lyon 1; INSA Lyon, F69621 Villeurbanne Cedex, France; 2European Synchrotron Radiation Facility, X-Ray Imaging Group, 38043 Grenoble, France; 3Centre for Medical Image Computing, University College London, WC1E BT, London, UK; 4Centre for Image Analysis and Science for Life Laboratory, Uppsala, University Box, 337SE-751 05 Uppsala, Sweden; 5UMPC Univ Paris 6, UMR 7623, Laboratoire d’Imagerie Paramétrique, 75005 Paris, France e-mail: Pei.Dong@creatis.insa-lyon.fr; alexandra.pacureanu@cb.uu.se; M.ZuluagaValencia@cs.ucl.ac.uk; cecile.olivier@creatis.insa-lyon.fr; quentin.grimal@upmc.fr; peyrin@esrf.fr (Received October 31, 2013; revised April 13, 2014; revised May 29, 2014; accepted May 29, 2014) ABSTRACT In the context of bone diseases research, recent works have highlighted the crucial role of the osteocyte system. This system, hosted in the lacuno-canalicular network (LCN), plays a key role in the bone remodeling process. However, few data are available on the LCN due to the limitations of current microscopy techniques, and have mainly only been obtained from 2D histology sections. Here we present, for the first time, an automatic method to quantify the LCN in 3D from synchrotron radiation micro-tomography images. After segmentation of the LCN, two binary images are generated, one of lacunae (hosting the cell body) and one of canaliculi (small channels linking the lacunae). The binary image of lacunae is labeled, and for each object, lacunar descriptors are extracted after calculating the second order moments and the intrinsic volumes. Furthermore, we propose a specific method to quantify the ramification of canaliculi around each lacuna. To this aim, a signature of the numbers of canaliculi at different distances from the lacunar surface is estimated through the calculation of topological parameters. The proposed method was applied to the 3D SR micro-CT image of a human femoral mid-diaphysis bone sample. Statistical results are reported on 399 lacunae and their surrounding canaliculi. Keywords: 3D image analysis, cortical bone, morphology of lacuno-canalicular network, ramification of canaliculi, synchrotron micro-CT The osteocyte system is hosted in the lacuno- INTRODUCTION canalicular network (LCN). Although the LCN is Bone diseases severely affect the quality of life raising increasing interest, it is a structure difficult to around the world. According to a survey of Inter-assess since it is deeply embedded in calcified bone national Osteoporosis Foundation, 1 of 3 women and matrix and it is composed of a huge number of 1 of 5 men with age over 50 suffer from osteoporosis, nanometric structures. Lacunae are flat ellipsoidal which is associated to bone fractures. Bone fragility voids housing the osteocyte bodies and with a density remains only partially understood despite decades of between 26000 to 90000/mm3 (Cardoso et al., 2013). research in this area. Recently, the crucial role of the Canaliculi are narrow tunnels, enclosing the cell osteocyte system at the cellular scale was highlighted processes of the osteocytes with a reported diameter (Bonewald, 2011). The osteocyte system, which is around 100–700 nm (You et al., 2004). The canaliculi composed of osteocytes communicating through connect lacunae and allow the circulation of interstitial dendrites, is essential in bone mechanosensation and fluid. Osteocytes are hypothesized to be stimulated bone mechanotransduction and is supposed to play an by load-induced interstitial fluid displacement circu­important role in orchestrating bone adaptation lating in the LCN. Computational transport models have (Burger and Klein-Nulend, 1999; Han et al., 2004; been proposed to estimate shear stresses (Weinbaum Gu et al., 2007). et al., 1994; Anderson et al., 2005; Anderson and Knothe Tate, 2008), but they are generally based on idealized lacuno-canalicular geometries. Therefore, the assessment of the LCN geometry is important to help understanding the strain-sensing mechanisms of the osteocyte network. However, the 3D assessment of the LCN is still limited. Conventionally two-dimensional imaging tech­niques, such as light microscopy, atomic force micro­scopy (AFM), scanning electron microscopy (SEM) or transmission electron microscopy (TEM), have been proposed to visualize and quantify the LCN (Marotti et al., 1995; Kamioka et al., 2009; Lin and Xu, 2011; Sharma et al., 2012). Confocal microscopy can provide three-dimensional images, however, the limited penetration of light restricts the thickness of the visualized specimen and the spatial resolution is anisotropic (Sugawara et al., 2005). Recently, new 3D imaging techniques such as ptychography and serial FIB/SEM have been reported to image the LCN respectively at spatial resolutions of 65 nm (Dierolf et al., 2010) and 30 nm (Schneider et al., 2011). Due to the very high spatial resolution of these techniques, the field of view remains limited to a few cells, which restricts the possibility to obtain statistically significant results on the osteocyte network. 3D synchrotron radiation microtomography (SR micro-CT) has been postulated to be an ideal technique for visualizing the LCN (Müller, 2009). In previous works, we have demonstrated the feasibility of SR micro and nano-CT to image the LCN at the European Synchrotron Radiation Facility (ESRF) (Langer et al., 2012; Pacureanu et al., 2012). SR CT presents a number of advantages over conventional cone-beam CT. In particular, since the photons flux is several orders of magnitude higher than the one of conventional X-ray tubes, sub-micrometric spatial resolution can be reached while keeping good signal to noise ratio in relatively short time. With a 3D iso­tropic voxel size of 280 nm, this technique can provide a relatively large field of view (FOV) of about 6003 µm3. Compared to the other 3D approaches, a few thousands of cells may be imaged in a single scan. However, due to the novelty of these images and the complexity of the LCN, the current quantification tools cannot provide an automatic analysis of the LCN. Previous efforts have been made for the segmentation of LCN from the reconstructed images (Pacureanu et al., 2009, 2010, 2011; Zuluaga et al., 2011). In order to quantify this network, techniques for extraction of dedicated parameters need to be designed. Our aim is to propose new automated methods to quantify the LCN from such 3D SR micro-CT images. The next section briefly describes image acquisition and segmentation and presents the methods used to quantify the lacunae and canaliculi. In particular we present a new scheme to assess the distribution of canaliculi around each lacuna. The last section reports the validation of the method on a simple geometrical phantom, an isolated lacuna and finally results on a large volume encompassing 399 lacunae. MATERIALS AND METHODS SAMPLE PREPARATION The specimen was extracted from a woman (78 years old) femoral mid-diaphysis, obtained from a multi-organ collection. The sample was prepared fol­lowing a series of treatments (Biobank, Presles en Brie, France), which underwent delipidation, elimi­nation of medullary protein and sterilization. Then, the sample was stored at a room temperature until the imaging acquisition. The collection of sample followed the procedure of the Human Ethics Committee of the “Centre du don des Corps” at the University Rene Descartes (Paris, France) and is in accordance with legal clauses stated in the French Code of Public Health. SYNCHROTRON RADIATION MICRO­ CT IMAGING Image acquisition was performed on the 3D parallel beam SR micro-CT setup developed at beam-line ID19 at ESRF (Salomé et al., 1999). The pixel size was fixed to 280 nm providing a FOV of about (600 µm)3. The X-ray beam energy was set to 18 keV. To reduce the effect of radiation dose on the sample, it was necessary to optimize the detector setup (Pacureanu et al., 2012). A 4.5 µm thick LSO scintillator was used to transform the X-ray into the visible light which was detected with an E2V CCD camera. During the scan, the samples were rotated around a vertical axis with a total angle of view of 180° and 2000 projections were recorded. The total time for each scan was around 30 minutes. After scanning, the 3D image was reconstructed using a standard filtered backprojection algorithm. The recon­structed image has an isotropic voxel size of 280 nm. Fig. 1a illustrates a region of interest from a trans­versal slice of a reconstructed image (904 × 649 pixels). The organization of the lacuno-canalicular network concentrically distributed around the Haver­sian canal can be observed. SEGMENTATION The segmentation of the LCN from the micro-CT images is challenging due to the small size of cana­liculi compared to the voxel size in the reconstructed image. Here, we followed a method introduced in a previous work to perform the segmentation (Pacureanu et al., 2010). First, a mask including the relevant parts of the bone tissue and excluding Haversian canals is computed. This was generated by using a large median filter (radius = 25) followed by a simple thres­holding (see Fig. 1b). The segmentation of the LCN was performed through a variational region-growing method (Pacureanu et al. 2010; Revol-Muller et al., 2012) based on an energy functional combining grey level information from the original image and shape information extracted via a 3D tube enhancement filter based on the eigenvalues of the Hessian (Sato et al., 1998; Pacureanu et al., 2013). Subsequently, the lacunae and canaliculi are separated by a 3D opening operation, with a 3 × 3 × 3 ball-shaped structuring element. After this step, two binarized images, of the lacunae and the canaliculi are available. Fig. 1c illustrates a composite image of the segmented LCN. We can formalize the segmented LCN as the union of each lacuna Ln with its set of canaliculi Cn: M LCN = U(Ln U Cn ), (1) n=1 where M is the total number of lacunae in the image. QUANTIFICATION METHOD Labeling In order to quantify each lacuna Ln with its set of canaliculi Cn, the segmented lacunae are labeled by using a connected component analysis. Here, a fast labeling method proposed by Hoshen and Kopelmanm (1976) was implemented in 3D. This method allows labeling of all objects by only scanning the image twice in a raster pattern, i.e. from top to bottom and bottom to top. Quantification of lacunae After labeling, the total number of the lacunae (N.Lc) and the lacunar density with respect to the bone volume (N.Lc/BV) and the tissue volume (N.Lc/TV) can be easily calculated. The bone volume (BV), the tissue volume (TV) and the canal volume fraction (Ca.V/TV) were quantified from the mask image. Then, a number of descriptors were extracted from each labeled lacuna Ln as described below. (a) (b) (c) Fig. 1. (a) A slice of reconstructed image of SR micro-CT (voxel size 280 nm, FOV: 904 × 649 pixels), including a Haversian canal, lacunae and canaliculi. (b) A slice of mask image which excludes the region of Haversian canal. (bone tissue in white and porosity in black) (c) The corresponding seg­mented lacunae (in yellow) and canaliculi (in green). Fig. 2. (a) Segmented image with 3 flat elliptic-like lacunae and many tube-like canaliculi. (b) Labeled lacunae after the connected component analysis. Different colors show the different labels. The 3D renderings were generated with the software Avizo®. Moment based descriptors Considering the ellipsoidal shape of the osteocyte lacunae, we used a 3D moment-based approach to fit each lacuna to an ellipsoid (Denis et al., 2008; Flusser et al., 2009). Let µpqr(Ln) be the second-order central moments of Ln, with p+q+r=2 for (p,q,r).[0,2] and M(Ln) be the corresponding 3 × 3 second order central moment matrix. We denote .k (k = 1,2,3) the eigenvalues of M(Ln) sorted by decreasing value. The half axes of the fitting ellipsoid having the same moment matrix are given by: ak = 5.k Lc.V (L ), k =1,2,3 , (2) n where Lc.V(Ln) is the lacunar volume, which equals to the total mass of the object (zero-th order moment). From the fitting ellipsoid, the osteocyte lacunar length (Lc.L1(Ln) = 2a1), width (Lc.L2(Ln) = 2a2), and depth (Lc.L3(Ln) = 2a3) were calculated. In addition, the lacunar anisotropy was quantified as Lc.L1(Ln)/Lc.L2(Ln) and Lc.L2(Ln)/Lc.L3(Ln). The ratio between the volume of the fitting ellipsoid and the actual one is given by: .= 20. 5... 3Lc.V (Ln )5/2. (3) 123 Besides, we also calculated the lacunar volume ratio with respect to the bone volume (Lc.V/BV) and the tissue volume (Lc.V/TV). Finally, the distance between nearest lacunae (Lc-Lc) was estimated by computing for each lacunae the minimum distance between its gravity center to that of the closest osteo­cyte lacuna. Intrinsic volumes In addition, we calculated the intrinsic volumes to further quantify the surface area Lc.S(Ln), the Euler number Lc..(Ln) and the structure model index (SMI) Lc.SMI(Ln) of the lacuna Ln (Ohser and Schladitz, 2009). The intrinsic volumes are invariant geometric functions serving as a basis of object features. In our three-dimensional image with the orthorhombic cubic lattice of spacing ., the four intrinsic volumes Vj(Ln), j .[0,3], respectively representing the Euler number, integral of mean curvature, surface area and volume of the lacuna, can be efficiently calculated from the approximation V)3-k of V3-k , expressed as: .^3-k (k ) V3-k (Ln ) = 2..vl hl (Ln ) , (4) l=0 where l is an index spanning local configurations of voxels, hl (Ln) is the number of configuration l in the set Ln and vl (k ) is the weight coefficient of the local configuration. The surface area Lc.S(Ln), the integral of the mean curvature Lc.M(Ln), and the Euler number Lc..(Ln) of the lacuna Ln can be esti­mated as: ) 2 (1) Lc.S(L ) = 2V (L ) = 4. vh , n 2 n l ) Lc.M (L ) =.V (L ) = 2..v(2)h , (5) n 1 n l ) Lc..(L ) = V (L ). n 0 n 4 (1) The surface area weight vl , integral of the mean curvature weight 2.vl (2) , and the weights vl (3) for the Euler number for each configuration are given in (Ohser et al., 2009). In addition, we also calculated the structure model index of each lacuna (Lc.SMI(Ln)) from the intrinsic volumes (Ohser and Schladitz, 2009), expressed as: Lc.SMI (L ) n . (6) = 12Lc.V (Ln )Lc.M (Ln ) Lc.S (Ln )2 It characterizes the lacunar shape with values of 0, for a pure plate, 3 for a rod and 4 for a sphere (Hildebrand and Rüegsegger, 1997). Quantification of ramification of canaliculi We propose a new method to characterize the 3D ramification structure of canaliculi issued from a lacu- na. Considering that the canaliculi are branching after radiating from the surface of lacunae, our goal is to give a signature of the number of canaliculi per la­cuna at different distances from the lacuna boundary. This method is based on the use of the 3D Euler number. Let us recall that the Euler number of an object in a three-dimensional space can be expressed as a function of the Betti numbers ß0, ß1, and ß2 as: . = ß0 -ß1 + ß2, (7) where ß0, ß1, ß2 represent respectively the number of connected components, the number of tunnels and the number of cavities in the object (Odgaard, 1997). In practice, when the object is represented as a discrete set of voxels, the 3D Euler number can be calculated as (Serra, 1982; Toriwaki et al., 2002): . = n0 -n1 + n2 - n3, (8) where n0, n1, n2, and n3 are respectively the number of vertices, edges, faces, and voxels. To calculate the number of canaliculi per lacuna at different distances from the lacunar surface, first we extract the bounding surface of the lacuna using morphologic derivative operations. The bounding surface .rLn can be formulated as the subtraction between the r + 1 and r times dilated lacuna: . L =. L \ . L rn (r+1) n rn , (9) where . denotes dilation. This bounding surface will be used to count the number of canaliculi at the distance r, from the surface of the lacuna. Each boun­ding surface .rLn has a single connected component with one cavity inside. Second, a special bounding surface with holes is generated as: H L =. L \(. L I C ) . (10) rn rn rnn After that, the Euler number of the bounding sur­ face with holes H (L ) is considered as: rn .(H L ) =ß (H L ) -ß (H L ) +ß (H L ) , (11) rn 0 rn 1 rn 2 rn where ß0(H L ) =1 (one connected component), rn ß1(H L ) = the number of holes on the bounding sur­ rn face, ß2(H L ) = 1 (one cavity inside). rn Finally, the number of canaliculi per lacuna Lc.NCa(L)(r) at a given distance r can be calculated n as: Lc.NCa(L )(r) =ß (H L ) = 2 -.(H L ) . (12) n 1 rn rn The dilation parameter reflects the distance at which the count is performed. Thus for each lacuna, the set {Lc.NCa(Ln )(r) r = r0,rmax } gives a signature of the distribution of canaliculi around the lacuna. The rmax was chosen smaller than the mean inter lacunar dis-tance Lc-Lc. The evolution of the number of cana-liculi with r reflects the ramification of the canaliculi. To further characterize the branching of canaliculi, we also calculated the ratio between the maximum and minimum number of canaliculi per lacuna at different distances from the surface of lacunae. Besides, we estimated an average length of the primary canaliculi for a given lacuna rp(Ln). The calculation was done under the assumption that the probability of any canaliculi to bifurcate is uniform in a certain range of distance. When considering only the first bifurcation, the number of canaliculi with respect to the distance (r) to the surface of lacunae increases linearly. Therefore, the relation between the number of canaliculi versus the distance was fit by a straight line Lc.NCa(Ln)(r) ~ a(Ln)r + b(Ln). Then, rp(Ln) is estimated as the length for which 1.5 times the initial number of canaliculi Lc.NCa(Ln)(r0) and by cor­recting the distance step .r/2: 1.5Lc.NCa(L )(r ) - b(L ) .r n 0 n rp (Ln ) =+ . (13) a(Ln )2 Finally, the descriptive statistics of each parameter were calculated on a whole population of lacunae in a given bone sample. RESULTS VALIDATION OF THE CANALICULI COUNTING METHOD The quantification of numbers of canaliculi was validated on two test volumes. The first one is a simple geometrical phantom constituted of three ellipses each with six cylinders, simulating three lacunae interconnected by canaliculi (cf. Fig. 3a). The volume size is 643 voxels. The second is the seg­mentation of one isolated lacuna with all canaliculi interconnected, shown in Fig. 3b. The volume size is 149 × 149 × 85 voxels. The application of the method to the first simple phantom provided 6 canaliculi for each ellipse as expected by construction. Concerning the second experimental image of an isolated lacuna, the validation was carried out by manual check and by a test program. The test program consisted in applying a connected component analysis after removing the lacuna. The method was applied with two choices of the dilation parameter r. For r equals 1 (respectively 15), the number of canaliculi found was 22 (respectively 32). The comparison with the test program in these two cases (Fig. 3c,d) shows that the correct number of canaliculi was obtained. The difference between the two results reveals that canaliculi are branched as can be observed in the 3D renderings. Thus, the dilation parameter is useful to put in evidence the ramification of canaliculi. (a) (b) (c) (d) Fig. 3. 3D rendering of test volumes used for the validation. (a) A phantom made of three ellipses each with six rod-like branches, volume size: 64^3 voxels. (b) An isolated lacuna with fully interconnected canaliculi, volume size: 149 × 149 × 85 voxels. (c)­ (d) Canaliculi counted by test program and manual check, showed in different colors. (c) 22 canaliculi for dilation parameter r = 1 (d) 32 canaliculi for dilation parameter r = 15 (4.2 µm away from the lacunar surface). APPLICATION TO SR MICRO-CT IMAGE The proposed method was applied to a large SR­micro-CT image acquired at the ESRF at a spatial resolution of 280 nm. A Volume of Interest (VOI) made of 904×649×998 voxels was extracted from a (2048)3 reconstructed volume. The VOI includes an osteon with hundreds of lacunae Fig. 1a. In the seg­mented osteocyte lacunae image, unacceptable objects such as micro-cracks, ring artifacts, and some irregu­larly shaped fragments can be extracted as well. From the literature, it is known that the size and shape of lacunae resides in a certain range (Mc-Creadie et al., 2004; Dong et al., 2014). Thus, the segmented lacu­nae were filtered by thresholding some descriptors as follows: 50 µm3 < Lc.V < 1000 µm3, Lc.L1/Lc.L2 < 5, and 0 . Lc... 2. In addition, we also removed the outliers, by excluding the lacunae for which the canalicular regression slope was smaller than 0.5 or higher than 50. (a) (a) Fig. 4. 3D rendering of the labeled lacunae (displayed in different colors) concentric distributed around the osteon. (a) top view (b) side view. The segmented image included 399 lacunae and their associated canaliculi. Fig. 4 illustrates the la­beled lacunae in different colors. We report the bone tissue histomorphometric parameters in Table 1. Besides, the statistical results on the lacunar des­criptors with the mean, standard deviation, minimum and maximum are reported in Table 2. The average distance between nearest lacunae was calculated and was found to be 23.2 µm. The mean and standard deviation of lacunar volume and surface was 216.4 ± 84.7 µm3 and 238.9 ± 66.4 µm3. The average lengths of the three axes of lacuna are 15.2 µm, 7.8 µm and 4.0 µm, and length ratio between the three axes is about 4:2:1. The mean of lacunar SMI is 3.2. To quantify the numbers of canaliculi, we used 8 different dilation parameters (r), from 5 to 40 with step of 5 corresponding to distances from 1.4 µm to 11.2 µm with step of 1.4 µm. Boundary lacunae were excluded from analysis, due to incomplete LCN. The statistical results are reported in Fig. 5. The average number of canaliculi per lacuna was found to increase between 41.5 and 139.1 with the distance r. The linear regression explains 98.0% of the variability of the number of canaliculi around their mean (see Fig. 5). The average length of the primary canaliculi is about 5.6 µm. The median of ratio between the maximum and minimum number of canaliculi per lacuna is 3.4. Fig. 6 illustrates 3 dif­ferent lacunae with different number of canaliculi. Fig. 6a shows a lacuna with large number of 98 canaliculi, Fig. 6b, the one with most frequent number of 65 canaliculi, and Fig. 6c, a lacuna with a small number of 49 canaliculi. The three ramification patterns calcu­lated for each of these lacunae are also shown in Fig. 6d. Fig. 5. The mean and standard deviation of the number of canaliculi counted at the eight different distances from the surface of lacunae (n = 240; boundary lacunae were excluded from analysis). DISCUSSION In this paper, we proposed a technique to auto­matically calculate LCN descriptors from a segmen­ted 3D SR micro-CT image. Both lacunae and cana­liculi were quantified. Considering the lacunae, we calculated both conventional shape parameters and parameters based on ellipsoidal fit. An original method to assess the ramification pattern of the canaliculi was proposed based on mathematical morphology operators and the Euler number. It consists in counting the numbers of canaliculi per lacunae at different distances. Conversely to a method that would try to identify each node of the canaliculi network, this method has the advantage to be more flexible to imperfections in segmentation which are unavoidable with the limited spatial resolution compared to the canaliculi size. The estimation of the numbers of canaliculi per lacunae was checked on a simple phantom and a complete lacuna with fully interconnected canaliculi. Then, this method was applied to a larger experimen­tal image of human bone tissue. This is the first time that the numbers of cana­liculi at different distances are measured in 3D. We have extracted these measurements on several hundreds of cells in human femoral bone. Traditional investiga­tions are based on 2D histological sections analyzed manually. Compared to the extrapolation of 2D mea­surements to 3D, the present method does not require any assumption on idealized lacuno-canalicular models. Thus this method is expected to yield unbiased para­meters. Since the method is fully automated, it can be applied to bone samples within large field of view including several osteons. The canaliculi were counted in 240 lacunae, after excluding boundary lacunae. Keeping all lacunae would provide biased results due to truncation of LCN. The average number of canaliculi near the surface of lacunae is in agreement with the results published by Beno (Beno et al., 2006) who reported an estimated 3D number of canaliculi of 41 in human bone. The low and high estimates are also consistent with Beno’s work. In a more recent work, the number of canaliculi per lacuna was calculated from confocal microscopy images on a limited depth of 25 µm (Sharma et al., 2012). The authors reported about 85 primary and 387 secondary canaliculi per lacuna on rats’ tibial metaphysis. These higher numbers can be explained by differences in method, species and bone site. In our method, the various dilation parameters permitted to put in evidence the increased branching of canaliculi with the radial distance to the lacuna. Our results clearly indicate an increase of the number of canaliculi with the distance to the lacunae. In addition, for most lacunae, the increase could be explained by a linear regression between the number of canaliculi and the distance suggesting the bifur­cation of the canaliculi. (a) (b) (c) (d) Fig. 6. 3D renderings of three isolated lacunae cropped from the segmented image. (a) A lacuna with 98 canaliculi counted at 5.6 µm to surface of lacuna. (b) A lacuna with 65 canaliculi counted at 5.6 µm to surface of lacuna. (c) A lacuna on the boundary of the osteon with 49 canaliculi counted at 5.6 µm to surface of lacuna. (d) Plots of the number of canaliculi calculated at eight different distances from the surface of the three selected lacunae. Table 1. Histomorphometric parameters of the bone tissue from the human femoral sample. Lc.TV BV TV BV/TV Ca.V/TV Lc.TV/BV Lc.TV/TV N.Lc/BV N.Lc/TV N.Lc (mm3) (mm3) (mm3) (%) (%) (%) (%) (mm-3) (mm-3) 399 8.63 × 10-5 12.39 × 10-3 12.85 × 10-3 96.4% 3.6% 0.70% 0.67% 32200 31042 N.Lc – number of osteocyte lacunae; Lc.TV – total lacunae volume (mm3); BV – bone volume (mm3); TV – tissue volume (mm3); BV/TV – bone volume fraction (%); Ca.V/TV – canal volume fraction or bone porosity (%); Lc.TV/BV and Lc.TV/TV– lacunar volume density (mm-3); N.Lc/BV and N.Lc/TV– lacunar number density (mm-3) Table 2. Osteocyte lacunar descriptors from the human femoral sample. Lc-Lc Lc.V Lc.S Lc.L1 Lc.L2 Lc.L3 Lc.L1/Lc.L2 Lc.L1/Lc.L3 Lc.SMI (µm) (µm3) (µm2) (µm) (µm) (µm) mean±std 23.2 ± 4.9 216.4 ± 84.7 238.9 ± 66.4 15.2 ± 3.7 7.8 ± 1.8 4.0 ± 1.0 2.1 ± 0.8 4.1 ± 1.5 3.2 ± 0.2 [min, max] [10.4, 38.8] [52.8, 611.9] [79.2, 883.3] [6.4, 49.7] [3.6, 14.5] [1.9, 7.0] [1.1, 4.9] [1.4, 11] [2.5, 3.8] Lc-Lc: nearest distance between two lacunae (µm); Lc.V – lacuna volume (µm3); Lc.S – lacuna surface area (µm2); Lc.L1, Lc.L2 and Lc.L3 – length, width and depth of lacuna (µm); Lc.L1/Lc.L2 and Lc.L1/Lc.L3 – anisotropy of lacuna; Lc. SMI – structural model index of lacuna The results strongly rely on image segmentation which itself is dependent on the quality of the acquired images. Further advances in optimizing the imaging technique may improve image quality, while progress in the segmentation methods would reduce noise and improve the connectivity of the lacuno­canalicular network. Although the numbers of cana­liculi that we found are realistic, it is difficult to evaluate the method due to the lack of comparison data and due to the complexity of the network. These preliminary observations need to be confirmed by further works. The presented method has the advantages of producing 3D measurements in an automated and model independent manner. We believe that this can contribute to improving our knowledge on the LCN. Furthermore, future work will also be pursued to extract additional characteristics of the LCN, such as the canaliculi length and the canaliculi density per each lacuna. ACKNOWLEDGEMENTS The authors acknowledge the support of the ESRF for providing beamtime on beamline ID19 for syn­chrotron micro-CT imaging in the context of the Long term Project (LTP) MD431 and the help of the ID19 group during data acquisition. This work was partly supported by the project "Microtomos" (DVO­20081013492) from the FRM (Fondation de la Recherche Médicale). It was also done in the frame­work of the French GdR Stic-Santé (CNRS 2643) and the LabEx PRIMES (ANR-11-LABX-0063) of Uni­versité de Lyon, within the program (ANR-11-IDEX­0007) operated by the French National Research Agency. REFERENCES Anderson EJ, Kaliyamoorthy S, Alexander JID, Tate MLK (2005). Nano-microscale models of periosteocytic flow show differences in stresses imparted to cell body and processes. Ann Biomed Eng 33:52–62. Anderson EJ, Knothe Tate, ML (2008). Idealization of pericellular fluid space geometry and dimension results in a profound underprediction of nano-microscale stres­ses imparted by fluid drag on osteocytes. J Biomech 41:1736–46. Beno T, Yoon Y-J, Cowin SC, Fritton SP (2006). Estimation of bone permeability using accurate microstructural measurements. J Biomech 39:2378–87. Bonewald LF (2011). The amazing osteocyte. J Bone Miner Res 26:229–38. Burger EH, Klein-Nulend J (1999). Mechanotransduction in bone–role of the lacuno-canalicular network. FASEB J 13(Suppl):S101–12. Cardoso L, Fritton SP, Gailani G, Benalla M, Cowin SC (2013). Advances in assessment of bone porosity, per­meability and interstitial fluid flow. J Biomech 46: 253–65. Revol-Muller C, Grenier T, Rose J-L, Pacureanu A, Peyrin F, Odet F (2012). Region Growing: Adolescence and adulthood ;Two visions of region growing: in feature space and variational framework. VISAPP 1, Rome, Italy, 286–97. Denis EP, Barat C, Jeulin D, Ducottet C (2008). 3D complex shape characterization by statistical analysis: Application to aluminium alloys. Mater Charact 59: 338–43. Dierolf M, Menzel A, Thibault P, Schneider P, Kewish CM, Wepf R, Bunk O, Pfeiffer F (2010). Ptychographic X-ray computed tomography at the nanoscale. Nature, 467:436–9. Dong P, Haupert S, Hesse B, Langer M, Gouttenoire P-J, Bousson V, Peyrin F (2014). 3D osteocyte lacunar morphometric properties and distributions in human femoral cortical bone using synchrotron radiation micro-CT images. Bone 60:172–85. Flusser J, Suk T, Zitov B (2009). Moments and moment invariants in pattern recognition. J. Wiley, Chichester, West Sussex, U.K.; Hoboken, N.J. Gu G, Kurata K, Chen Z, Väänänen KH (2007). Osteocyte: a Cellular Basis for Mechanotransduction in Bone. J Biomech Sci Eng 2:150–65. Han Y, Cowin SC, Schaffler MB, Weinbaum S. (2004). Mechanotransduction and strain amplification in osteo­cyte cell processes. PNAS, 101:16689–94. Hildebrand T, Rüegsegger P (1997). Quantification of bone microarchitecture with the structure model index. Comput Methods Biomech Biomed Engin 1:15–23. Hoshen J, and Kopelman R (1976). Percolation and cluster distibution. I. Cluster multiple labeling technique and critical concentration algorithm. Phys Rev B 14:3438-45. Kamioka H, Murshid SA, Ishihara Y, Kajimura N, Hasegawa T, Ando R, et al. (2009). A method for observing silver-stained osteocytes in situ in 3-microm sections using ultra-high voltage electron microscopy tomography. Microsc Microanal 15:377–83. Langer M., Pacureanu A, Suhonen H, Grimal Q, Cloetens P, Peyrin F (2012). X-ray phase nanotomography re­solves the 3D human bone ultrastructure. PLoS ONE, 7: e35691. Lin Y, Xu S (2011). AFM analysis of the lacunar­canalicular network in demineralized compact bone. J Microsc 241:291–302. Marotti G, Ferretti M, Remaggi F, Palumbo C (1995). Quantitative evaluation on osteocyte canalicular density in human secondary osteons. Bone 16:125–28. McCreadie BR, Hollister SJ, Schaffler MB, Goldstein SA (2004). Osteocyte lacuna size and shape in women with and without osteoporotic fracture. J Biomech 37: 563–72. Müller R (2009). Hierarchical microimaging of bone struc­ture and function. Nat Rev Rheumatol 5:373–81. Odgaard A (1997). Three-dimensional methods for quanti­fication of cancellous bone architecture. Bone 20:315– 28. Ohser J, Nagel W, Schladitz K (2009). Miles formulae for boolean models observed on lattices. Image Anal Ste­reol 28:77–92. Ohser J, Schladitz K (2009). 3D images of materials structures: processing and analysis. Wiley. Pacureanu A, Langer M, Boller E, Tafforeau P, Peyrin F (2012). Nanoscale imaging of the bone cell network with synchrotron X-ray tomography: optimization of acquisition setup. Med phys 39:2229–38. Pacureanu A, Larrue A, Peter Z, Peyrin F. (2009). 3D non­linear enhancement of tubular microscopic bone poro­sities. IEEE ISBI, Boston, USA, 602–5. Pacureanu A, Revol-Muller C, Rose J, Ruiz MS, Peyrin F (2010). Vesselness-guided variational segmentation of cellular networks from 3D micro-CT. IEEE ISBI, Rotterdam, The Netherlands, 912–5. Pacureanu A, Rollet J, Revol-Muller C, Buzuloiu V, Langer M, Peyrin F (2011). Segmentation of 3D celluar networks from SR-micro-CT images. IEEE ISBI, Chicago, USA, 1970–73. Salomé M, Peyrin F, Cloetens P, Odet C, Laval-Jeantet AM, Baruchel J, et al. (1999). A synchrotron radiation microtomography system for the analysis of trabecular bone samples. Med Phys 26:2194–204. Sato Y, Nakajima S, Shiraga N, Atsumi H, Yoshida S, Koller T, et al. (1998). Three-dimensional multi-scale line filter for segmentation and visualization of curvi­linear structures in medical images. Med Image Anal 2:143–68. Schneider P, Meier M, Wepf R, Müller R (2011). Serial FIB/SEM imaging for quantitative 3D assessment of the osteocyte lacuno-canalicular network. Bone 49: 304–11. Serra J (1982). Image analysis and mathematical morphology. Academic Press. Sharma D, Ciani C, Marin PAR, Levy JD, Doty SB, Fritton SP (2012). Alterations in the osteocyte lacunar-canali­cular microenvironment due to estrogen deficiency. Bone 51:488–97. Sugawara Y, Kamioka H, Honjo T, Tezuka K, Takano-Yamamoto T (2005). Three-dimensional reconstruction of chick calvarial osteocytes and their cell processes using confocal microscopy. Bone 36:877–83. Toriwaki J, Yonekura T (2002). Euler number and connec­tivity indexes of a three dimensional digital picture. Forma 17:183–209. Weinbaum S, Cowin SC, Zeng Y (1994). A model for the excitation of osteocytes by mechanical loading-indu­ced bone fluid shear stresses. J Biomech 27:339–60. You L, Weinbaum S, Cowin SC, Schaffler MB (2004). Ultrastructure of the osteocyte process and its pericel­lular matrix. Anat Rec A Discov Mol Cell Evol Biol 278:505–13. Zuluaga MA, Dong P, Pacureanu A, Orkisz M, Peyrin F (2011). Minimum cost path approach for the seg­mentation of bone canalicular network from nano-CT images. IEEE NSS/MIC Valencia, Spain.