INTRODUCTION
Osteoporosis represents one major health condition for our growing elderly population[1,2]. It accounts for severe morbidity and increased mortality in postmenopausal women and it is becoming an emerging health concern even in aging men[3-5].
Despite a lack of well-known critical hormonal changes, responsible for a fast degenerative action on bone tissue[2,6], the male structure undergoes a slow progressive loss of bone mass and the risk of being affected by osteoporosis in late age is increasing even in men[7-9]. Nevertheless, perimenopausal female population still now shows the highest risk of early osteoporosis onset. The increased lifespan in the industrialized world accounts for the increasing incidence of osteoporosis and bone fractures with a perspective of additional years (at least 20) of disability osteoporotic women have to face in their later life.
Osteoporosis, as defined by World Health Organization, is a systemic disease of the skeleton characterized by low bone mineral density (BMD), and microarchitectural deterioration of bone tissue with consequent increased bone fragility that predisposes to fracture risk. Due to the silent progression of bone structure degeneration, osteoporosis diagnosis often follows a painful fracture event. To date, only a small percentage of individuals have been known to be osteoporotic while the condition of most pathologic people had remained undiagnosed until a fracture occurs[10-12]. It is worth noting that the first fracture event increases the risk and accelerates the onset of new ones[13].
Screening of the population at risk for bone degenerative processes and treatment of osteoporotic patients to prevent or reduce bone fragility fractures would improve quality of life in aging people. Early diagnosis of bone deterioration by tools sensitive to osteoporosis onset and progression would represent a promising approach to prevent disability and reduce mortality in osteoporotic patients.
OSTEOPOROSIS AND MEDICAL IMAGING
BMD estimate by means of dual-energy X-ray absorptiometry (DXA) is normally used in clinical practice for osteoporosis diagnosis. BMD alone, however, has been found to be not a good predictor for fracture risk[14]. Recently, detection of BMD by means of ultrasound techniques has been obtaining more attention as a possible less invasive diagnostic tool for bone degenerative disorders. At present, echographic investigation in clinical practice remains the “alternative” way to assess bone mass only in those subjects who cannot undergo X-ray irradiations. Quantitative ultrasound (QUS) detection of bone degeneration is becoming ever more promising and potentially useful mainly in screening studies[15]. This method, in fact, offers advantages over X-ray-based techniques, including low cost, portability, and non-ionizing radiation exposure[16-18]. In addition, it is known that ultrasound transmission depends not only on mineral content, but also on biomechanical properties of the bone, thus providing extremely useful information about bone structure properties. As a matter of fact, several studies have highlighted that bone parameters measured by QUS provide information about fracture risk independent of BMD[19-22].
The incidence of bone fractures is not always correlated to diminished BMD[23]. In fact, too high of a mineral content could predispose to a higher risk of bone fractures, probably due to an increased rigidity of the structure. Moreover, pharmacological treatment with agents acting by preventing and/or reducing bone loss is not always able to reduce bone fracture risk[24]. Bone microarchitecture has been found to be a determinant of bone fragility independent of bone density, thus representing an interesting aspect of bone strength to give insight into patterns of bone degenerative processes in both aging and osteoporosis. From a clinical point of view, understanding the role microarchitecture plays in the mechanisms of bone fragility as well as in the action of drugs to prevent fractures would improve clinical management of osteoporotic pathology[25].
For a long time, the unique font of information about bone microarchitecture in humans had been coming from studies on biopsies of bone tissue from autoptic samples.
In the last few decades, the impressive diffusion into clinical practice of medical imaging techniques has opened up new perspectives to characterize bone structure in vivo. The rapid development of information technology together with advances in power calculation capacity has solicited proliferation of different new methods to assess bone microarchitecture changes with aging and pathology.
Several processing methods have been reported as suitable tools to quantify bone tissue in both classical sites at high risk of fracture, such as spine and femur head, and mirror sites, such as calcaneus, tibia, and radius[26,27]. Most of them uses axial radiographic images[28,29] or computer tomography, CT, projections[30] to in vivo characterize trabecular bone architecture (TBA).
Magnetic resonance imaging
Magnetic resonance imaging (MRI) has recently emerged as a suitable tool to measure trabecular bone structure in vivo[31]. In particular, MRI-based diagnosis of TBA deterioration could be used to complement standard BMD-DXA measurement for assessing osteoporosis onset and progression.
Several aspects of this technology, in fact, candidate MRI as a non-invasive, non-ionizing tool for in vivo study of bone tissue in human beings[32]. MRI is based on the interaction between a high-gradient magnetic field, radiofrequency pulse transmission, and protons in the tissues under investigation. Most clinical MRI systems have field strengths of 1.5 T, 3.0 T being quite rare. Bone tissue contains few mobile water and fat protons and appears as low signal intensity (hypointense) on MRI. Bone marrow, on the other hand, has abundant fat and water protons and appears as high signal intensity (hyperintense). Hypointense cortical and trabecular bone appear sharply outlined against hyperintense bone marrow and juxatacortical fat. High field strengths and improved coil technology have made it possible to achieve an in-plane resolution of about 150 μm[33]. However, because of current signal-to-noise constraints, minimal slice thickness is usually about 500 μm[34]. Although this resolution is larger than trabecular thickness, trabecular spacing (800-1000 μm) is still much larger than the size of a single MR voxel, so that trabecular parameters can be reasonably measured or estimated with MRI.
MR-derived structural dimensions are similar, but not identical, to histological or micro-computed tomography (CT) dimensions and a two-dimensional rather than three-dimensional approach is applied[34,35]. A magnetic susceptibility difference of about 3 ppm existing between bone and marrow leads to signal dampening at the bone-marrow interface and, as a result, an artificial overestimation of trabecular dimension occurs. This artefact, known as trabecular broadening, is responsible for an apparent increase in trabecular thickness that vary with pulse sequence applied and field strength[36]. MR pulse sequence strongly affects MR susceptibility artefact but it is not clear yet which pulse sequence is the best for trabecular bone imaging[37,38]. Interestingly, trabecular broadening is potentially advantageous to MRI since it enables visualization of small trabeculae that would not normally be seen due to partial volume effect. Most researchers currently use either gradient-echo or spin-echo imaging. Spin-echo technique has the advantage to reduce trabecular broadening artefact[38,39].
Trabecular structural parameters measured with 3.0 T are better defined than those at 1.5 T using micro-CT as a reference standard. Although correlations between 3.0 T and micro-CT are better than 1.5 T, increased trabecular broadening has been observed at 3.0 T due to an increase in susceptibility effect[36].
High correlation has been found also for trabecular number and spacing by comparing 3.0 T MR and high resolution pQCT. Histomorphometry and other different structural analysis techniques applied to both modalities provide different absolute values but high correlation (r > 0.8) for all bone structure parameters[34].
Micro-MR units with special high-field-strength (7 T, 11 T or higher) are used in experimental settings to obtain high-resolution MR images of small animals in vivo or bone specimens in vitro. Images obtained at about 90 μm yield TBA scalar parameters similar to those obtained with histomorphometry methods[40,41].
MRI AND CLINICAL ASSESSMENT OF BONE STRUCTURE
A recent review of osteoporosis imaging literature highlights that high-resolution MRI quantification of TBA deals with two different aspects of osteoporosis assessment: prevalence and incidence of osteoporotic fractures and TBA response to drug therapies. Early studies suggest that trabecular parameters as measured on MR images are able to separate patients with and without osteoporotic fractures better than BMD[42,43]. These two groups are better characterized by indices that quantify trabecular shape transformation from plate-like to rod-like by using structural model index method or plate-to-rod ratio method. Also finite element model produces good results. However, classic morphometric parameters such as bone volume fraction (BV/TV), trabecular number (Tb.N), and trabecular thickness (Tb.Th) show best results when compared to BMD[32,44].
Measuring the effect of drug therapies for osteoporosis, incidence of bone fragility fractures represents a crucial point as its assessment requires much time. BMD has then become the best single surrogate marker of bone strength. Nevertheless, it is widely demonstrated that BMD presents severe limitations. As a matter of fact, in the MORE study, including 7700 women treated with raloxifene, a 40% reduction of fracture risk has been highlighted with only 4% associated to increased BMD[45].
TBA parameters from MRI would monitor changes induced by antiresorptive therapies better than BMD. TBA analysis by using MRI has been demonstrated to be more sensitive than BMD even in monitoring several therapeutic effects: salmon calcitonin in different skeleton sites[46] testosterone at distal tibia[47]; alendronate in peripheral districts[48]. By comparing CT and MRI[32,49], CT has the advantage to visualize bone tissue with a higher space resolution while has the disadvantage to high doses of radiation when applied to central skeleton districts such as vertebras and femur neck, the two best sites to estimate fracture risk in primitive (postmenopausal and senile) osteoporosis. MRI with high field modalities (3, 7 T, and higher) have been introduced to significantly improve signal-to-noise resolution of the image; nevertheless, higher field strengths also increase magnetic susceptibility induced effects responsible for alterations in bone structure parameters[34,50]. However, bone structural parameters measured by means of MRI techniques have been widely correlated to in vitro histological ones. Therefore, MRI emerges as a promising non-invasive non-ionizing tool for in vivo characterization of TBA potentially useful in early diagnosis of fracture risk in osteoporosis. Its adoption and diffusion in clinical practice are limited because of texture analysis, mainly based on methods producing a set of numerous calculated parameters that make difficult their interpretation[51]. Most methods, in fact, are based on classical histomorphometric techniques that provide for a large number of measurements (up to 25) to be analyzed[51,52]. A promising MRI-based approach to define bone structure uses wavelet techniques[53]. Wavelet analysis is generally applied in functional medical imaging and finds successful application with MRI[54]. Several applications of wavelet techniques for texture analysis of different biological structures have been already described[55,56]. Attempts to reduce computational charge and speed up the methods are still under study; however, the time-cost remains too high for a routine application in clinical practice[55].
Recently, a method of MR image analysis has been proposed to provide a limited number of parameters sensitive to bone microarchitecture changes in aging and pathology[57-59]. The method has been developed taking into account biocomplexity of human beings and fractal properties of many physio-anatomic structures[60-66] and it is the result of a proper combination of image acquisition, texture analysis, and mathematical solution to the study of TBA.
A NEW APPROACH TO THE STUDY OF TRABECULAR BONE ARCHITECTURE
The introduction in biomedical field of paradigms, such as theory of complexity, chaos and fractals, suggests new approaches and provides innovative tools to the study of morpho-functional degenerative processes with aging and age-related pathologies, among is osteoporosis[61,65,67-72]. Combining biocomplexity analysis tools together with advanced techniques of image processing and image analysis, it is possible to develop computerized methods able to provide e limited number of parameters sensitive to age-related changes as well as to pathology onset and progression[65,66].
Recently, an original method of TBA analysis has been developed able to draw out one numerical index sensitive to physio-pathologic changes of the structure. The basic idea to build such a method stems from the complexity of living being and fractal nature of many physio-anatomic structures[31]. Several structures and functions of human body have got fractal properties[62,72,73]. Trabecular bone, in addition to be a good model for fractal analysis in biological structures, offers the opportunity of analyzing fractal lacunarity[74]. In this particular context, fractal lacunarity analysis appears the most suitable approach to define trabecular bone network. Fractal lacunarity, in fact, by measuring space-filling capacity of a complex object, has the potential to describe both bone network discontinuity and sizes of trabecular spaces (bone marrow), whose changes represent an index of increased bone fracture risk.
The method has been developed on MR images of lumbar vertebras. Spin-echo technique has been used to visualize the trabecular structure in the inner portion of vertebral body. The procedure adopted to set up the standard version of the method is described and discussed hereinafter and is schematically represented in Figure 1. It provides for fractal lacunarity analysis in a region of interest (ROI) after a pre-processing of the image to optimize trabecular network visualization. The parametric characterization of curvilinear graph, result of fractal analysis, is calculated by using our original bio-mathematical model based on hyperbola function model. In fact, hyperbola formula contains three coefficients (α, β, γ) that represent our suitable numerical indices, where α correlate with the fractal dimension, and β, related to the concavity, characterizes the lacunarity. The result is a triplet of parameters (α* β* γ*) that univocally characterizes any single TBA analyzed[57,58].
Figure 1 Schematic representation of fractal lacunarity analysis.
A: 1.5 T MRI spin-echo image (512 × 512 pixel, pixel size equal to 0.469 mm) of fifth axial section of lumbar vertebra; B: Rectangular ROI within the inner perimeter of vertebral body in a intermediate step of binarization procedure on the reverse gray-scale image. The plot (right bottom) represents the result of gliding box algorithm application (dotted line) as fitted by hyperbola function (solid line) used to calculate the triplet of parameters α, β, γ. MRI: Magnetic resonance imaging; ROI: Region of interest.
Application of the method to several TBA from 25 female subjects with different age and physio-pathologic status (4 young, 5 pre-menopause, 6 post-menopause, 10 osteoporotic with an age-range of 31-37, 42-52, 51-81, and 59-75 respectively) has highlighted that, among the three coefficients, parameter β is particularly sensitive to both age and physio-pathologic changes. In particular, parameter β correlates with physio-pathologic status and assumes decreasing values from healthy young to perimenopausal to osteoporotic patients. Results also show that parameter β is statistically different (probability significance value P≤ 0.5 for Student-Newman-Keuls multiple range test) in the three classes of TBA considered (healthy young, perimenopausal, and osteoporotic)[57,58]. A correlation between parameter β and age can also be observed with a decreasing trend of β values from young to old subjects[58]. It has to point out that the healthy old subject (in this context without clinical signs of osteoporosis) shows a β value higher than the younger osteoporotic patients. This evidence candidates parameter β as a standard for TBA characterization as well as an index of structure integrity, potentially predictive of bone fracture risk. Low values of β correspond to a decreased structural integrity linked to increased fracture risk[57,58,65].
Results from the application of three software prototype versions[57-59,66], differing in the pre-processing step, on the same set of TBA images confirm the potentiality of the proposed bio-mathematic model. Note that in Zaia et al[57,58] this method has been presented for binary images arising from a pre-processing step of image I by using image J program. The simple extension considered[58,59] can be used also for gray-scale images. The efficiency of such an extension is usually improved by a different pre-processing step through a sigmoid function (Figure 2). In particular, sigmoid function operates a sort of rescaling of gray tones by weighting the gray level of each single pixel: this procedure allows limiting information loss due to binarization procedure[59,66].
Figure 2 Representation of sigmoid function for different choices of parameters k, σ.
A: Sigmoid graph obtained by varying values of parameter k; B: Graph obtained by varying values of parameter σ. This new pre-processing step prescribes to consider image J, in place of image I, where the pixels are defined as follows: J(i, j) = 1/1 + esp(-k (I(I,j) - σ), i = 1,2,…M, j = 1,2,…N, and k, σ > 0 are two given parameters. Note that the procedure goes toward a complete binarization by increasing parameter k, related to sigmoid rectangularization.
Figure 3 shows comparison of results from the first version of the algorithm, (described hereinafter) developed to work on binary images (black/white) obtained by pre-processing original gray-scale MR images in a different computational environment[57,58], and version 2 where the algorithm has been improved to directly work on original gray-scale images in a unique computational environment[66]. In both cases correlation of parameter β with both age and physio-pathologic changes is confirmed. It is worth noting that a fixed pair of values for parameters k and σ has been chosen such that for the triplet of coefficients α, β, γ, from lacunarity analysis we obtain results comparable with those of the previous studies[57,58]. In particular, k equal to 7 and σ equal to 0.7 have been used to obtain results reported in Figure 3. This procedure increases the sensitivity of the method as suggested by the large range of values parameter β can assume when compared with those from black/white version.
Figure 3 Representation of trabecular bone architecture changes with aging and osteoporosis from lacunarity texture analysis on binary and gray-scale images.
Results are expressed as parameter β calculated on the same set of original images with two different software prototype versions: A: On binary images after a pre-processing step in a different computational environment; B: Directly applied to original gray-scale images. Black circles represent TBA lacunarity in non-osteoporotic subjects of different age. Gray circles represent TBA in osteoporosis (BMD-based diagnosis). Decreasing values of parameter β correspond to increasing lacunarity related to a higher microarchitecture deterioration. TBA: Trabecular bone architecture; BMD: Bone mineral density.
When version 2 of the method, more sensitive than the first one, is applied to a larger sample size (59 female subjects: 13 young, 17 pre-menopause, 14 post-menopause, 15 osteoporotic with an age-range of 23-39, 40-52, 50-81, and 57-78 respectively), a certain degree of overlapping of β values in osteoporotic and non-osteoporotic subjects can be observed (Figure 4). These results confirm the scientific evidence that TBA can represent a predictor of fracture risk independent of BMD.
Figure 4 Effect of sample size on trabecular bone architecture changes with aging and osteoporosis.
Results are expressed as parameter β calculated from lacunarity texture analysis with gray-scale version 2 of the method. Larger sample size (B) of 59 subjects comprises the smaller one (A). Black circles represent TBA lacunarity in non-osteoporotic subjects of different age. Gray circles represent TBA in osteoporosis (BMD-based diagnosis). Decreasing values of parameter β correspond to increasing lacunarity related to a higher microarchitecture deterioration. Overlapping of data from osteoporotic and “healthy” subjects (BMD-based diagnosis) stresses that TBA represents a determinant of bone fragility independent on BMD. TBA: Trabecular bone architecture; BMD: Bone mineral density.
Version 3 of the software prototype always works on original gray-scale images but, by omitting gray-scale inversion and introducing a different analysis of pixel gray levels by sigmoid function, allows correlating parameter β directly to fracture risk: that is, increasing β values correspond to increased TBA deterioration, index of increased fracture risk[59]. Lacunarity β values, shown in Figure 5, are such that a clear age-related healthy status trend can be drawn, osteoporotic-related values being much higher. In this case, parameter β ranges within a larger set of values and makes the method much more sensitive, such that it identifies three different clusters of values that can be related to three different levels of structural integrity: a lower cluster corresponding to healthy status (low fracture risk); a higher one corresponding to osteoporotic status (high fracture risk); an intermediate one where pathologic and non-pathologic subjects are mixed and can be associated to a similar fracture risk (Figure 5). Once again, it is stressed that TBA could predict fracture risk better than BMD alone. Note that four out of ten menopausal subjects show β values similar to the lowest osteoporotic ones. This intriguing observation suggests that these menopausal subjects could have a high risk of pathology onset while the osteoporotic patients could have a lower fracture risk than the other osteoporotic ones. Further studies enrolling a large number of subjects differing for both age and healthy status need to support this hypothesis.
Figure 5 Results from trabecular bone architecture texture lacunarity analysis of magnetic resonance images with version 3 of the method.
Lacunarity of TBA from 30 subjects differing in age and physio-pathologic status is expressed by parameter β. Circles: age-related alteration of TBA in “healthy” subjects; Triangles: pathologic-related changes of TBA in osteoporotic BMD-based diagnosis subjects. Parameter β values directly correlate with fracture risk. TBA: Trabecular bone architecture; BMD: Bone mineral density.
From these preliminary studies, the proposed method emerges as a potential diagnostic tool for an improved characterization of osteoporotic pathology useful in early diagnosis of bone fragility fractures as well as in the assessment of therapy efficacy in preventing or decreasing fracture risk.
The methodological approach adopted to characterize TBA texture is highly innovative at both clinical and technological levels. The result from the proposed original mathematical solution is a triplet of parameters (α, β, γ); parameter β, particularly sensitive to TBA physio-pathologic alterations, emerges as a potential standard for TBA characterization as well as a parametric index for a best clinical management of osteoporosis[65].
The improved versions of the method, based on gray-scale texture analysis, has the advantage of being directly applicable to original trabecular bone MR images by working in a unique computational environment[59,66]. These easier and faster versions of the method would speed up its application in both research and clinical practice.
OUTLINE OF THE METHOD
The method for TBA texture characterization by lacunarity analysis is made up on three main steps: image acquisition by MRI spin-echo multislice technique; image processing of vertebral axial images to produce a binary representation of trabecular bone structure; quantitative lacunarity analysis of trabecular network by hyperbola model function approximation of the graph produced by gliding box technique (Figure 1). The method has been developed on real images to avoid surprises during clinical application. In fact, often it happens that methods built up on simulation models fail when applied to real clinical images. The procedure adopted to set up the method, described in Zaia et al[57,58] is reported and discussed below.
MR images acquisition
Lumbar vertebras images were obtained by means of a standard clinical MRI system (1.5 T whole-body imaging system-Gyroscan Intera; Philips Medical System, ACR-Nema 1.0). Spin-echo multislice technique was applied to acquired a set of 9 axial vertebral images (512 × 512 pixels, pixel dimension equal to 0.469 mm) once the fourth lumbar vertebra (L4) is spotted (Figure 6). Clinical application of spin-echo multislice technique is still quite rare; therefore, vertebral axial section images were obtained on purpose from female patients underwent MR assessment for injuries of the column. Twenty-five subjects were enrolled in this study and classified as follows: 4 healthy young (mean age 33.7 years, age range 31-37); 10 perimenopausal (mean age 51.7 years, age range 43-65); 10 osteoporotic, on the base of classical BMD-DXA diagnosis (mean age 68.4 years, age range 64-74); 1 “healthy” 80 years old woman, in this case without clinical signs of osteoporotic pathology.
Figure 6 Set of images from magnetic resonance imaging spin-echo.
Nine axial images of the fourth lumbar vertebra obtained by magnetic resonance imaging multislice technique. Acquisition of the set of axial images of vertebral body to visualize the inner trabecular bone portion is performed from bottom (A: lower base) to top (I: upper base).
Image processing and analysis of MR images
MRI spin-echo data set of lumbar vertebras was transferred to a personal computer for image texture analysis. The middle axial image of L4 was considered to build the method. Two main considerations solicited the choice of the fifth axial section image. A certain degree of both inter and intra-vertebral heterogeneity had to be expected; therefore, only one L4 section, namely the fifth one was considered to set up the method. This choice makes easier to compare results from different subjects. In addition, limiting analysis procedures to one image only would save time and facilitate technological transfer of the method into clinical practice. The goodness of this choice is supported by results from a pilot study: a certain degree of variability among the sections was observed in any data set produced by the application of the method to the whole set of axial images of different vertebras. It is worth noting that the fifth section of any vertebra considered generally displayed intermediate values.
The algorithm for TBA characterization was developed on images pre-processed as described in Boutry et al[51] in a different computational environment by using ImageJ (version 1.26t), a software for image processing and image analysis available for free on the web. This pre-processing step, briefly described in Figure 7, generates a binary representation of TBA where the black regions represent bone marrow and the white ones represent bone. This binarization scheme has the advantage of keeping and/or excluding the “false-positive finding” of trabeculae in the black regions as well as of removing residual small artefacts.
Figure 7 Schematic representation of the pre-processing step of magnetic resonance images with image J program used in trabecular bone architecture lacunarity analysis method-version1.
Reverse grey-scale display of original images is used to ease visualization of trabecular network, a rectangular ROI is chosen as large as possible to fit the inner perimeter of trabecular bone. The boundary between cortical bone and soft-tissue background is defined by using an automatic contour detection algorithm. Trabecular bone is then extracted from bone marrow by using a Gaussian filter. This scheme generates a binary representation of the trabecular structure in which the dark regions represent bone marrow and the light regions represent bone. The Gaussian filter is applied first to the original ROI image and then to the thresholding one. The two resulting binary images are finally combined. ROI: Region of interest.
Estimate of fractal lacunarity
Gliding box algorithm, GBA, is the most popular among the definitions and calculating procedures proposed to estimate lacunarity[75]. It is based on the analysis of mass distribution in the set and involves the variance of box mass, M, at each step, wherein the box is moved one by one space unit. A frequency distribution of box masses, n(M,b), where b is the size of the gliding box, is produced by recounting the box mass throughout the whole region.
It is possible to consider only a discrete frequency distribution n(Mj,b), j=1,2,…,μ(b) by assuming that, for each b, only a finite number of masses Mj, j=1,2,…,μ(b) are found in the various gliding boxes of size b. It is worth noting that this assumption is true for binary images, where the mass of a generic box on the image is given by the number of pixels associated to the value 1 in the box, that is, the white pixels. From standard arguments the moments of order q of M, through the GBA, are given by:
Math 1
where the division by N(b) is necessary to convert n(Mj, b), j = 1, 2,…,μ(b) into a probability distribution. The definition of lacunarity function Λ is based on the first and the second moments of M only, that is
Math 2
The GBA method was implemented in a prototype software by using the MATLAB software package, version 6.1 (the MatWorks, Inc.). The program starts elaborating the binary image generated from the above mentioned binarization pre-processing step of MR lumbar vertebra images. It calculates lacunarity values, for each integer value of b comprised between bmin and bmax, where bmin, bmax are given integer multiples of the pixel size in the image considered. Once the lacunarity function Λ(b) is calculated, the program displays the results as a graph.
It was observed that the curvilinear behavior of lacunarity function resembles the hyperbola one (Figure 4, plot) for all the images analyzed; therefore, the following model function:
Math 3
was chosen to approximate the curvilinear plot of lacunarity where the three coefficients α, β, γ are potential parameters to quantitatively define TBA.
Note that the theoretical behavior of lacunarity function Λ for both ideal fractals and other complex random sets further stresses this observation. In addition, for self-similar fractals parameter α is related to the fractal dimension and parameter β characterizes the fractal lacunarity of the set[75]. The best interpretation of lacunarity Λ(b) by the hyperbola model function L(b) is computed as a least squares problem solution. The three coefficients α, β, γ are the independent variables and minimizer of this problem is a triplet of parameters, α* β* γ*, of the hyperbola model function that better describes mass density variation of pixels in any image considered[57-59].
Parameter β as an estimate of lacunarity
The method was applied to different types of bone trabecular structure to test the robustness of the method. The subjects included in this study were women chosen to be representative for age and physio-pathologic status: women in their thirties were analyzed as bone loss begins approximately at that age. Perimenopause subjects were also included in the study as they represent a condition characterized by a high risk of osteoporosis onset. Taking into account the possibility of a physiological aging, a healthy old woman (in this case without clinical signs of osteoporosis), was also considered. The pathologic patients were classified as osteoporotic on the basis of classical diagnosis based on BMD-DXA measurement. Application of the method to the same original image was repeated step by step for six times. Results from the fifth axial section image of L4 were obtained by fixing bmin equal to 5 pixels and bmax equal to 50% of ROI diagonal length, measured in pixels.
The three parameters (α* β* γ*) from our mathematical model define as univocal any lacunarity function, L(b) where α represents the order of convergence of L(b), β defines the concavity of the hyperbola, and γ is the translation term on X-axis. Parameters α and γ assume quite similar small values in all bone structures analyzed while parameter β ranges within a larger set of values. A low β value defines a hyperbola function with a slight concavity (high lacunarity), whereas a high β value defines a hyperbola function with a deeper concavity (low lacunarity)[57,58] (Figure 8).
Figure 8 Example of results from lacunarity analysis of three types of trabecular bone different for age and healthy status.
In the related plots, dotted lines represent the fitted hyperbola function. Goodness of the fitting represented by a almost complete overlapping of both theoric and experimental curve indicates that parametric lacunarity analysis by hyperbola model function is a proper choice. A: Young, 39s; B: Menopause, 53s; C: Osteoporosis, 74s.
Improvements of the method are in progress to validate the goodness of lacunarity parameter β as a marker of TBA deterioration in both aging and pathology. In particular, we have been focusing on both thresholding and segmentation procedures. ROI size and shape are also under consideration to improve robustness and reliability of the method before approaching a clinical study on a large size sample of subjects at risk of bone fragility fracture.
PECULIARITY OF THE METHOD
Different potential benefits can be expected from clinical application of the above described method to assess osteoporotic fracture risk. First of all, it has been built on vertebral images. From comparative studies dealing with structural analysis in different sites, such as distal radius, calcaneus, and spine, or even in minor sites such as wrist or finger, spine represents the best site to early predict bone fragility fracture status[76]. Second, the method uses MRI system to acquire axial section images of the fourth lumbar vertebra, described as the best site to assess osteoporosis. MRI has been reported as ideally appropriate for TB texture analysis. In fact, MRI spin-echo scans of TB, where bone marrow appears with uniform hyperintense signal while bone has low background intensity, generate a binary tomographic system that simplifies TB network thresholding. Furthermore, using spin-echo technique to obtain axial sections of TB allows avoiding imaging artefacts of trabecular enlargement as it occurs with gradient-echo[51].
Fractal lacunarity analysis of TBA represents the most intriguing benefit of the proposed method. In fact, fractal analysis has been becoming a common tool to quantify TBA complexity; however, only fractal dimension is applied and generally measured by the Hmean parameter[26,28,51]. As mentioned above, fractal lacunarity appears the most suitable tool to analyze TB texture as it has the potential to describe both bone network discontinuity and sizes of trabecular spaces, changes of which are an index of increased bone fragility[58]. In addition, dealing with lacunarity analysis as a more general approach to analyze complex patterns with or without fractal properties[75,77,78], it allows overcoming biomedical limits of fractal analysis, responsible for misinterpretation if neglected. Lacunarity analysis of TB in human vertebras has been previously described only for micro-radiography and CT images of bioptic specimens[79,80]. Fractal lacunarity analysis of TB in MR spin-echo images of the vertebra, therefore, represents an original application of a non-ionizing, non-invasive method to assess TBA[57,58]. A MRI-based method to in vivo assess osteoporosis in men has been also proposed[51]. It is based on the gradient-echo technique, less precise than spin-echo[38,39], to acquire sagittal images of the calcaneus. Dealing with clinical application of medical imaging techniques to assess osteoporosis, the choice of mirror sites, such as the calcaneus, is justified only for invasive/ionizing methods. Furthermore, histomorphometric methods are used for texture analysis of TB in the cited study[51]: the result is a set of 20 parameters, 13 of which reported as significantly different in osteoporotic patients. In the last decades, a lot of studies have been performed in both central and mirror sites with different image acquisition techniques. Most of them always propose histomorphometric-based methods alone or combined with other methods for TB texture analysis, among are anisotropy, co-occurrence, gradient matrices, gray level histogram, and runlength[81,82]. All of them provide quite a few parameters to be checked and analyzed thus limiting their application in clinical practice. As a matter of fact, 8 parameters of 32 calculated in[29] and 9 of 24 in[52] have been reported as featuring TBA changes.
This observation highlights the last, but not least, advantage of the method, that is the result: a triplet of values, α*, β*, γ*, corresponding to the three coefficients (α, β, γ) of the hyperbola model function which univocally characterize each single TBA pattern considered[57,58]. Among the three coefficients, parameter β, being highly sensitive to TBA changes with aging and pathology, is potentially candidate as a standard for TBA assessment and a parametric index useful in early diagnosis of osteoporosis and therapeutic assessment to prevent bone fracture risk[57-59,65].
DISCUSSION
Several clinical studies support the evidence that MRI-derived measurements of TBA are able to describe changes in both aging and pathologic status, and can discriminate patients with vertebral or hip fractures from fracture free individuals; the best featuring of patients with or without bone fragility fractures is achieved by combining structural parameters and BMD measurements. Nevertheless, very few data are available in literature on treatment-related changes of TBA, and no prospective studies can be found dealing with fracture risk prediction. Prospective studies on osteoporotic fractures on large scale therapeutic trials are necessary to give insight into the role of in vivo assessment of TBA. These studies are beneficial to define a set of diagnostic markers able to complement or improve fracture risk diagnosis based on BMD-DXA[83-85]. They would also give insights into the search for therapeutical approaches, such as antiresorptive drugs or fracture healing agents, more effective in preventing and/or curing osteoporotic bone fractures[24,86-89]. It is worth noting that Canadian Guidelines for the assessment of fracture risk in osteoporosis have been recently updated by introducing bone quality estimate[90] and last generation DXA devices start to be equipped with software for bone quality assessment[91-94]. Main obstacles to reach this goal are represented by limited technology dissemination in healthcare centers, a minimum standardization of protocols for image acquisition and processing, the number of parameters to characterize TBA. As a matter of fact, one most recent longitudinal study on alendronate therapy[48] uses a 3T-MRI system for image acquisition of peripheral bone districts such as distal tibia, distal radius, and proximal femur. The parameters measured to characterize TBA include BV/TV, Tb.N, Tb.Sp, Tb.Th, and seven parameters from topologic geodesic analysis (GTA). Only four GTA parameters and apparent Tb.N result significantly modified after 24 mo treatment in distal tibia when compared to BMD.
Efforts have been done for calibration and standardization; however, comparative multicentric studies are necessary to lay the bases for future multicentric clinical trials and prospective studies. These goals appear ever more unreachable dealing with the proposed methodologies, ever more sophisticate and developed without taking into account feasibility of their application in clinical practice. As the result, in 2013 we still find clinical studies on TBA characterization performed on bioptic specimens and analyzed by quantitative histomorphometric method in micro-CT images[95].
CT techniques represent a precious contribution to the knowledge of bone mechanical properties useful for a more precise evaluation of fracture risk. Studies on CT applications in this context are in progress but they can be performed only in research centers as they require highly qualified personnel and sophisticated software and hardware systems. Diffusion of CT techniques into clinical practice is further limited by high dose irradiation when compared to DXA.
Many advantages can be obtained by widening MRI role in the assessment of bone fragility. MRI is a non-invasive/non-ionizing method that allows defining bone tissue structure in vivo. MRI data set acquisition can be performed in different arbitrary axes with an image acquisition time as quick as 10-15 min. Impressive advances in the development of MRI techniques to assess bone fragility have been done over the last few decades. Further progress has to be expected with improved image-processing and image-analysis methods. The advent of new high-resolution imaging techniques has introduced new stimuli to the study of bone tissue microarchitecture as it has been possible to overcome the limit of bi-dimensional analysis belonging to conventional techniques. The intrinsic three-dimensional nature of new instruments has allowed visualizing and analyzing bone specimens directly in 3D although this progress has not been associated with a proper standardization of investigated parameters. Frequently, classic histological indices have been simply mimicked or adapted, so that the new indices introduced and tested for new methodologies result inadequate for an exhaustive description of bone trabecular tissue structure and properties. New high-resolution MRI methods recently proposed to characterize TBA require ever more powerful and sophisticate instruments with high field strength modalities (3 T, 7 T, or 11 T) for both bi- and three-dimensional characterization of TBA and, once again, are based on classical histomorphometric analysis in peripheral bone districts[32,96,97]. It depends on the choice of micro-CT as the gold standard, thus alienating the introduction in clinical practice of non-invasive non-ionizing tools[98], even though promising.
The original and innovative method described, based on fractal lacunarity of vertebral TBA, appears particularly promising. The method uses images acquired by 1.5 T-MRI system, widely disseminated in most healthcare centers, and provides only one parameter highly sensitive to TBA changes thus representing a suitable method for an easy and fast applicability into both research and clinical practice. Currently, a study for clinical validation of the method is in progress in our Institute. The study has been designed to be observational, cross-sectional and prospective, and schedules the enrolment of women at risk of bone fragility fractures. From this study lacunarity parameter β could emerge as a standard for TBA characterization and as a marker candidate for osteoporotic fracture risk.
