Fractal dimension from radiographs of peridental alveolar bone A possible diagnostic

indicator of osteoporosis

Urs E. Ruttimann, PhD,” Richard L. Webber, DDS, PhD,b and Jane B. Hazelrig, Bethesda, Md., Winston-Salem, N.C., and Birmingham, Ala.

PhD,’

NATIONAL

THE

UNIVERSITY

INSTITUTES OF ALABAMA

OF HEALTH,

THE

BOWMAN

GRAY

SCHOOL

OF MEDICINE,

AND

AT BIRMINGHAM

The purpose of this study was to investigate whether a radiographic estimate of osseous fractal dimension is useful in the characterization of structural changes in alveolar bone. Ten dry mandibular bone segments were radiographed from three controlled projection angles (-5, 0, f5 degrees), before and after acid-induced partial decalcification. Fractal dimension was estimated by regression analysis of power spectra computed by Fourier transform of selected regions of interest in digitized images of the radiographs. Repeated-measures ANOVA showed that fractal dimension so determined varied over anatomic locations (p < .Ol), but increased after acid-induced demineralization (p < .0005), irrespective of the radiographic projection angles (p > .99). In vivo fractal dimension was computed from randomly selected intraoral radiographs of six premenopausal (ages, 32.8 f 3.9) and six postmenopausal (ages, 62.5 + 4.1) women. A significantly (p < .Ol) higher fractal dimension was observed in the older group. (ORAL SURC ORAL MED ORAL PATHOL 1992;74:98-110)

0

steoporosis is responsible for at least 1.2 million bone fractures in the United States each year; 44% occur at vertebral, 19% at hip, and 14% at distal forearm sites.’ The sequelae of these fractures are fatal in 12% to 20% of cases,2 and the direct and indirect costs of this bone disease were estimated in 1984 as $6.1 billion annually.3 Osteoporosis is defined as “a decrease in bone mass and strength leading to an increase in fractures.“4 Although, according to this definition, a fracture must occur for the disease to become clinically apparent, the disorder must exist before a trauma reveals its presence. Consequently, several techniques have been developed to measure bone mass density in patients. Cross-sectional studies indicate an age-related bone loss of about 1% per year in normal postmenopausal women. 5-7Hence, bone mass measurements require a dActing Chief, Diagnostic Systems Branch, National Institute of Dental Research, National Institutes of Health. bProfessor, Departments of Dentistry and Radiology, The Bowman Gray School of Medicine. ?Associate Professor, Biostatistics and Biomathematics, Health Affairs, The University of Alabama at Birmingham. 7/16/37312

98

very high precision if, as clinically desirable, the results are to be meaningful for the follow-up of individual patients. Ott et al8 compared the bone mass measurements in 49 osteoporotic women obtained at baseline and again after 1 year of treatment with the use of four techniques considered the most precise today: total body calcium by neutron activation, single photon absorptiometry at the radius, dual photon absorptiometry at the lumbar spine, and quantitative computed tomography at the spine. None of the measurements obtained by one technique could be used to predict changes as measured by any other of the techniques (no significant correlations). This was found to be true despite the fact that many of the measurements showed longitudinal changes that were greater than the reported precision of these methods (2% to 5%). In fact, patients who showed substantial gains by one method showed losses by other methods, and vice versa. Consequently, the most sophisticated techniques for bone mass determination are not precise enough to measure short-term changes in individual patients, and may be used only in the study of large groups of patients. Intraoral (periapical) radiographs that are rou-

Volume Number

74 1

tinely obtained by dentists to assessthe status of teeth and supporting bone generally include x-ray projections of parts of the mandibular or maxillary bones and typically have a high image resolution as a result of the direct film-exposure (nonscreen) technique used. Hence these radiographs would provide a vast source of data that already exist and could potentially be used for the longitudinal assessment of bone mass loss. However, the inability of x-ray film to distinguish between mass- and energy-dependent differences in attenuation (required to distinguish radiolucencies caused by soft- and hard-tissue changes, respectively) assures that film-based radiography could never attain the precision of the specialized methods mentioned above, which themselves are only borderline acceptable. These limitations persist even if meticulous contrast calibration and projection standardization procedures that were developed for digital subtraction radiography 9* lo are used. Therefore, to foster the possibility of eventually making use of the large amount of data that exist in radiographs of persons who receive dental care, we propose to investigate an alternative way to characterize the impact of the disease. As a decrease in bone strength is included in the clinical definition of osteoporosis, bone strength may be considered an equally valid indicator of disease status. As an alternative to attempts to measure loss of bone mass, we seek to assessother bone changes likely to be correlated with strength. We submit that functional changes that influence strength are not only determined by mass but also are related to structure (i.e., how mass is spatially distributed). Thus it is postulated that a likely correlate of bone strength can be assessed indirectly from the trabecular structure of the jawbones as revealed by the high-resolution images provided by intraoral radiography. As an objective measure to quantify bone structure, we propose to use Mandelbrot’s fractal dimension.” This measure has certain theoretical properties that make it attractive for retrospective use of radiographs, which will be discussed in the next section. However, to be applicable, the bone structure must sufficiently approximate certain general geometric properties that govern the complexity of its three-dimensional shape. Also, because it would be advantageous to use nonstandardized radiographs, this measure should not be sensitive to moderate changes of the radiographic projection geometry that typically occur in routine clinical follow-up. Therefore the purpose of this investigation was to explore whether (1) the fractal model is an adequate descriptor of the peridental bone pattern represented on intraoral radiographs, (2) the measured fractal dimension is independent of

Fractal

dimension of alveolar

bone

99

the x-ray projection angle, and (3) estimates of the fractal dimension of peridental bone determined from preexisting clinical radiographs are sensitive enough to differentiate between premenopausal and postmenopausal women. MATERIAL AND METHODS Theoretic rationale

The establishment of the geometric concept of fractals by Mandelbrot” has resulted in a multitude of applications across many scientific disciplines. In particular, fractal geometry has been found to be a useful descriptor of surfaces that range from a molecular scale12 to the topography of landscapes in the millimeter to kilometer range,13 all the way up to a characterization of continental faulting.14 The mathematic attribute that makes fractals so well suited to the description of many natural objects is the scaling property that states that a fractal object statistically looks the same under whatever magnification is used.15 Hence, a portion of an object scaled to its original size looks similar to the original object. This preservation of shape regardless of the observation scale is specified in the statistical sense: the spatial properties of a fractal object are similar in distribution. Or, more generally, if the object is described by different physical quantities (e.g., length and time), the probability distributions of the various quantities may follow different scaling laws, in which case the object is called statistically self-affine. Whereas these scaling laws hold for all sizes of the mathematic model, real objects are statistically selfaffine only over a finite range of observations. Despite this idealization, fractal models have been applied successfully to biologic objects such as the branching of the bronchial tree or the His-Purkinje system in the heart.t6 In medical imaging, fractal geometry has been used for the description of the dendritic tree in neurons,17 the analysis of nuclear medicine scans of the liverI and the lung,19 and the characterization of the parenchymal pattern in mammograms.20 More specific to bone imaging, a fractal texture analysis has been applied to radiographs of the OS calcis (heel bone)21 and intraoral radiographs.22 A promising theoretic basis that relates fractal geometry to transmission x-ray images of trabecular bone follows from the observation that the texture of bone tends to change from coarse to fine as one crosses the boundary from cortex to spongiosa. In flat bones such as the mandibles examined in this investigation, this transition is manifested largely in a direction parallel to the x-ray beam. Hence, the associated radiographic projection necessarily consists of a heterogeneous mixture of spatial frequencies characterized

100

Ruttimann,

Webber, and Hazelrig

by the underlying spectrum of trabecular constituents. To the extent that the power spectrum of this mixture is relatively continuous and constrained so that its curvature over a reasonable range is approximately log-linear, the projection that results over this range of spatial frequencies necessarily is well approximated by a single fractal dimension. (This conclusion follows directly from arguments presented in the Appendix.) Various methods have been proposed to determine the fractal dimension of a given data set.15 One technique that has been applied previously to estimate the fractal dimension of the trabecular bone pattern21 is based on Mandelbrot’s generalization of Brownian motion.23 However, the use of the autocorrelation function in that approach resulted in biased estimates when applied to noise-corrupted data,2’ and thus it was modified by the use of the power spectrum instead.22 Although the two methods are obviously related by the Wiener-Khintchin theorem,24 the latter approach offers the capability of partial noise removal, which will be discussed later. On the basis of the premise that the Fourier transform of f(x) exists, a simple heuristic derivation is given in the Appendix that relates the power spectrum to the fractal dimension. The decision to use the power spectrum instead of the original data for the estimation of the fractal dimension was made for two main reasons. First, a very efficient numeric method to compute power spectra is available that uses the fast Fourier transform (FFT) algorithm. 25 As a result, the computation burden to derive fractal dimensions from images becomes much less in the spectral than in the original data domain. Second, the spectral method provides the opportunity for some noise filtering. As indicated before, the scaling law (see equation 1 of the Appendix) holds for real objects only over a limited range of scales, r = l/f. For measurements made from radiographs, the smallest object detail that can be resolved is limited by the modulation transfer function (MTF) of the imaging system.26 Hence, below a certain scale, rmin, or above a spatial frequency, f,,, = l/rmin in the power spectrum, the film density fluctuations are no longer determined by the object but, rather, by the random variation of the x-ray quanta and the film grain, which establish a “white noise floor” or a constant power. Consequently, the negative slope-b-in the double logarithmic plot extends only up to a certain “corner” frequency, above which the spectral power asymptotically approaches a slope of zero. Therefore, all power values above that corner frequency are considered noise and can be eliminated before the estimation of the characteristic spectral slope.

ORAL

SURG

ORAL

MIZD

ORAL

PATHOL July 1992

The upper limit of the scaling law is determined by the selection of the length-L-of the data window. A consequence of the finite value of L is “spectral leakage,“25 which necessitates some corrective steps at the low end of the power spectrum. The leakage is the result of the convolution indicated in equations 4 and 6 of the Appendix and causes the substantial power contained in the continuous spectrum between the lowest (f = 0) and next lowest (f = l/L) FFT samples to spread into adjacent higher FFT samples. To avoid a possible bias because of this leakage, the power at f = 0 was forced to zero by subtracting the average gray value of the image within the data window from each pixel before applying the FFT. Furthermore, a few of the lowest FFT samples were routinely excluded from the fitting of the spectral slope. Experimental

methods

This research consisted of two experiments. The first is an analysis of radiographs obtained in vitro from two human mandibles. The bones were each sectioned into five segments that contained teeth as shown in Fig. 1. A special jig was used to position each segment reproducibly relative to the x-ray source. This enabled experimental control of the projection geometry of the radiographs to be made before and after artificial bone decalcification. Each cut specimen was radiographed at three different projection angles that spanned plus and minus 5 degrees with the use of 65 kVp x-rays filtered with 2.5 mm aluminum. The exposure for all specimens was identical (48 impulses at 10 mA) and the film used was nonscreen with an intermediate speed rating (Ultraspeed, Kodak, Rochester, N.Y.). The exposure sequence was repeated after decalcification of each specimen with 10% formic acid solution. Decalcification was accomplished with each specimen positioned so that one of the cut surfaces was oriented horizontally. The acid used to demineralize was dripped into the marrow space of this surface from above. The acid permeated the entire trabecular region through capillary action, which ultimately yielded a salt-laced filtrate that flowed out through the opposite cut surface. With a drip rate of approximately 16 drops per minute, this process required 35 minutes for each specimen. The entire filtrate was collected, and the calcium leached from each specimen was precipitated through the addition of ammonium oxalate. The amount of calcium removed from each specimen was determined by spectrophotometric analysis of the turbidity that resulted from the oxalate precipitation process. The analysis was calibrated through the use of control solutions that contained

Fractal

Volume 74 Number

dimension of alveolar bone

101

1

Specimen Preparation and Data Sampling Configuration

Incisor Region

3

Premolar Region

2 m

4

Premolar Region Molar

60

cm

Schematic Transaxial View of Typical Mandible Showing Relative Positions of Sectioning Cuts Dividing It into Five Segments

-

1

Reproducibly Positioned Bone Segment

=I

c X-Ray Film Packet

Example of Projection Geometry Used in Radiography of a Typical Jaw Segment

Fig. 1. Sampling configuration and projection geometry of in vitro specimens,

known quantities of calcium chloride that were chemically precipitated in a similar manner. The exposed films were processed in an automatic film processor (Model 08 I-PlO, Hope Industries, Willow Grove, Pa.) and digitized with a video camera with a charge-coupled device as the input transducer (Javelin JE2362, Atlantic Video, Birmingham, Ala.). The digitized video output of one image frame was stored as a 512 pixel X 512 pixel X 8 bit data matrix with a commercially available frame buffer (PCVision Plus, Imaging Technology, Woburn, Mass.). The optical image that reached the camera was produced with a macro lens (MVL-25752, Atlantic Video, Birmingham, Ala.), adjusted to a magnification sufficient to permit isolation of a small square portion of the radiographic image of interradicular trabecular bone that measured 7.5 mm X 7.5 mm with a resolution of 128 pixels X 128 pixels (Fig. 2). The window was positioned in a relatively uniform region made up predominantly of uniformly textured bone. The par-

ticular location was chosen to minimize the encroachment of the window on the projections of teeth and other anisotropic structures to the maximum extent possible as determined by visual inspection. An identical window was placed in near absolute registration over the radiographic images of the same tissues produced after demineralization. This was accomplished with registration techniques developed for subtraction radiography. lo The original digitized image was contrast-inverted and copied into the display memory as a reference. The postdemineralization radiograph was then placed into the field of view of the video camera and electronically superimposed over the reference by mixing the video signals that arose from the camera and the display memory, which accomplished a real-time subtraction. Spatial registration was achieved by manipulation of the radiograph while the mixed video signal was viewed on the monitor until recognizable structures were completely cancelled. Placement of the data window at identical pixel coor-

102

Ruttimann,

Webber, and Hazelrig

ORAL

Fig. 2. In vitro radiograph of typical specimen shows square window in peridental sures 128 X 128 pixels used to compute fractal dimension.

dinates in the corresponding digital images then yielded comparable registration with an estimated precision of f I pixel. To investigate whether this precise window registration procedure can be somewhat relaxed in future clinical applications, a secondseriesof data was generated with only visually matched subregions of the bones selected without the aid of the subtraction technique. Also, no attempt was made to circumvent or to calibrate the uncontrolled effects of the automatic gain control mechanismdesignedinto the camera. Hence, average exposure and large-area contrast gradients produced in the digitized images necessarily varied from the image of one specimento the next and produced differences in the power spectra at very low frequencies. This lack of control was intentional, to assessthe sufficiency of suppressingthis experimental variation with the elimination of the lowestfrequency components in the digitally computed power spectra. The average value of the digitized data within each window was subtracted and the result converted by the two-dimensional FFT algorithm into two-dimensional periodograms. One-dimensional power spectra were then produced with the expression of the spectral data in polar coordinates and the determination of the average of the power values acrossall anglesfor each of the 64 radial increments of the spatial frequency (i.e., the average was computed over circles with radii that corresponded to the discrete spatial

SURG

ORAL

Mm

ORAL

PATHOL July 1992

alveolar bone that mea-

frequencies). This yielded a smoothed, isotropic, onedimensional approximation to the spectral power distribution associated with the trabecular image and enabled characterization by a single slope-parameter. According to equation 9 in the Appendix, this parameter-b-was computed by linear regressionanalysis of the double-logarithmic plot of the power spectrum. The least-squaresfit to the slope was obtained after the lowest four frequency data were deleted to suppresspossibleartifacts due to large-area contrast gradients as described above. The MTF of the imaging process that yielded radiographs of the in vitro specimens was sufficiently broad so that all data points in the power spectrum up to the highest (Nyquist) frequency could be included in the fit. The effects of anatomic location, radiographic projection angle, and experimental decalcification on the measured fractal dimension were investigated by repeated-measuresANOVA (program BMDP 4V).27 In this analysis, the location and angle were considered as grouping (between) factors, and the data before and after calcification asrepeated measurements (within factor). The secondexperiment assessedthe clinical feasibility of the method. Of interest was whether the quality of intraoral radiographs was sufficient to characterize the fractal dimension and reveal structural differences between women of two age groups. Nonstandardized periapical radiographs were obtained retrospectively from dental records acquired

Fractal

Volume 74 Number 1

dimension of alveolar bone

103

Fig. 3. Portion of typical in vivo dental radiograph shows window measuring 64 X 64 pixels in peridental alveolar bone used to compute fractal dimension.

Table I. Comparison

of two random samples of premenopausal

and postmenopausal

Premenopausal Age

Mean SD

37 30 32 32 28 38 32.8 3.9

Postmenopausal

b 1.919 1.883 1.545 1.433 1.495 1.598 1.645* 0.206

.fc

D

Age

0.44 0.50 0.53 0.56 0.44 0.56 0.51 0.06

2.540 2.558 2.727 2.183 2.152 2.701 2.671 0.103

63 67 51 58 65 65 62.5 4.1

b = Negative slope of log(power) versus log(frequency). f, = Upper limit of the frequency components usedfor slope D = Fractal dimension computed from spectral slope. *p = 0.008, two-sample t test.

estimation,expressed asa fractionof

Institute on Aging as part of an ongoing survey of the health needs of aging Americans (Baltimore Longitudinal Study on Aging). Six premenopausal* (28 to 38 years of age) and six postmenopausal* (57 to 67 years of age) women were selected randomly (Table I). The radiographs were converted to 512 X 5 12 X 8 bit digital frames by a solid-state video camera (SC series, Applied Intelligent Systems, Inc., Ann Arbor, Mich.) that was interfaced with an image-processing system (IP 6400, by the National

*All participants had completemedicalhistories and hencemenopausal status was established

independent

women

of age.

b 1.364 1.317 1.394 1.364 1.290 1.428 1.360* 0.050

fc

D

0.53 0.50 0.44 0.53 0.56 0.53 0.52 0.04

2.818 2.841 2.803 2.818 2.855 2.786 2.820 0.025

the total band width.

Gould DeAnza, Fremont, Calif.). The optical magnification was adjusted so as to place a 64 X 64 pixel data window over the image area of the interradicular bone between the maxillary premolar and first molar teeth (Fig. 3). This resulted in sampling areas of 2.5 X 2.5 to 3.5 X 3.5 mm2. The subsequent data processing steps were the same as for the first experiment except for the fit of the spectral slope. Since the length of the data window was one half the size used in the first experiment, only the two lowest-frequency components were deleted before the fit. Also, the MTF of the imaging process in the case of these clinical radiographs was sufficiently limited so that noise

104

Ruttimann,

Webber, and Hazelrig

ORAL SURG ORAL MED

ORAL PATHOL July 1992

(Molar region, Calcified specimen) 6

1,

y = 5.7879 - 2.7052x r = 0.992 Estimated Fractal Dimension I 2.14

5-

4-

0

I 0.3

0.0

'

I 0.6

'

I 0.9

log Spatial

'

I 1.2

'

1 1.5

'

I 1.8

'

Frequency

Fig. 4. Representativeplot of power spectrumassociatedwith radiographobtainedfrom typical in vitro specimenshowsregressionline from which fractal dimensionwascomputed. Table II. Repeated-measuresANOVA Source

of fractal dimension, spatially registered images Mean Square

4

F

P value

Between Location (L) Angle (A) LxA Error

2 2

0.1315

6.11

0.0081*

0.000 1

4

0.0002 0.0215

0.00 0.01

0.9970 0.9997

21

Within Calcific. CXL CxA CxLxA Error

(C)

1

0.0504

17.40

0.0004*

2 2 4

0.0033 0.0015 0.0004

1.14 0.50 0.12

0.3373 0.6119 0.9723

21

0.0029

*Signifies statistical significance at p < 0.05.

became dominant above a certain “corner frequency”-f,-as discussedin an earlier section of this article. The corner frequency was determined by visual inspection of the spectral plots, and data points above its value were excluded from the slope fit. The decisionfor exclusion wasmade before the regressionused to estimate the slope, so that its impact on the final fractal computations was unknown at the time. The difference between the two age groups with respect to the spectral slopesor, equivalently, fractal dimension, was analyzed by the two-sample t test.

RESULTS In vitro data Fig. 4 shows a representative plot of data used to compute fractal dimension from the power spectrum associated with a typical specimen. The regression line ignoresthe first four points that correspond to the lowest spatial frequencies. All regressionscomputed in this way produced high correlation coefficients, r > 0.97. A repeated-measures ANOVA of the in vitro data showed that significant differences existed between experimental estimates of fractal dimension

Fractal dimension of alveolar bone 105

Volume 74 Number 1

2.5

Effect of Decalcification on Fractal Dimension of Spatially Registered Images Calcified c] Decalcified

II 5 !j o 5 5 o s 5 0 5 5 0 5 5 0 5 5 0 5 5 0 5 5 0 5 5 0 5 Angles lnasor V rremow prernoler rremotar Premolarr rvblar Molar r Location #5 #l #2 #3 #4 #6 #7 #a #9 # 10 Specimen Fig. 5. Bar chart showsrelationshipsthat exist betweenestimatedfractal dimensionthat existsbeforeand after acid decalcificationof in vitro specimens asfunctions of samplelocation and projection angle.

obtained from different anatomic locations (incisor, premolar, and molar regions), and also between those determined before and after acid decalcification. On the other hand, projection angle changeshad no significant effect on the final estimate of fractal dimension (Table II). To determine whether the angle independence demonstrated in this investigation dependson precise spatial registration of the subregionsselected for digitization and subsequent analysis, data from the independently selected but otherwise identical series of only visually matched subregionswere analyzed in the sameway. The relative angular independencewas retained under these conditions, but the lack of positional control did slightly reduce the significance level for the effects of the acid decalcification from 0.0004 to 0.0012. In other respects the data were comparable to those shown in Table II. The effects of decalcification per se are graphically shown in Fig. 5, in which the data are paired so that the dimensions measured before and after decalcification in each specimenare related. Notice that in 25 out of 30 casesthe computed fractal dimensionsafter decalcification are larger than before (p < 0.00001, by the sign test). A plot of fractal dimension versus data averaged

Table III. Relationship between lost calcium and associatedchange in fractal dimension Specimen No. 1

2 3 4 5 6 7 8 9 10

Location

(me)

Average fractal change

Incisor Incisor Premolar Premolar Premolar Premolar Molar Molar Molar Molar

58.08 64.75 29.32 32.96 23.94 33.93 67.78 41.71 28.94 36.13

0.0167 0.0433 0.2133 0.0534 0.0000 -0.0033 0.0800 0.1400 0.0467 0.0834

Lost Ca++

within anatomic locations is shown in Fig. 6. It appearsthat the structural differences between mandibular alveolar bone located near incisors and premolars are smaller than their differences relative to bone in the molar region. A quadratic regression was used to calibrate the spectrophotometric assay scheme from which the amount of calcium leached from each specimenwas determined. The resultant curve fit was excellent

106

Ruttimann,

Webber, and Hazelrig

Effect of Location of Bone on Fractal Dimension

ORAL

SUKC

ORAL

MED

ORAL

PATHOL July 1992

two groups. Similarly, in both groups, about 50% of the spectral components provided by the FFT could be used for slope fitting. The average slopes of the regression lines were significantly (p < .Ol) steeper for the premenopausal than for the postmenopausal women and yielded fractal dimensions of 2.68 and 2.82, respectively. Hence, the data indicate a difference of the bone structure between the two age groups. DISCUSSION In vitro data

Incisor

Premolar

Molar

Data averaged According to Location Fig. 6. Bar chart showsthe effectsof in vitro data averagedaccordingto anatomiclocation on estimatedfractal dimension.

(R2 > .999) and assured that the technique provided a reasonable quantitative estimate of the total calcium lost from each specimen. The data that related this aggregate loss per specimen to the change in fractal dimension produced in the associated subregion are shown in Table III. The correlation between these measurements was statistically not significant. In vivo data

Fig. 7 shows examples of the power spectrum plots obtained from intraoral radiographs of premenopausal and postmenopausal women. Data points excluded before fitting the regression lines are indicated by round symbols. The results of these regression analyses and the corresponding fractal dimensions are summarized in Table I. The correlation coefficients of the regressions ranged between 0.984 and 0.997, with no systematic difference between the

The excellent fit seen in Fig. 4 assures the validity of the self-affine assumption that underlies the application of a fractal model as a meaningful basis for data interpretation. As a result, it follows from equations 1 and 10 in the Appendix that the fractal dimension-D-provides a scale-independent measure of bone structure. This scale independence of the measure was indeed our theoretical motivation for its selection, with the strong practical consequence that it is by definition independent of the radiographic magnification associated with intraoral radiographs. The results from the ANOVA show that fractal dimension discriminates between anatomic location of the trabecular bone (p < 0.01) and the presence or absence of induced demineralization (JJ < O.OOl), irrespective of the projection angles used to produce the radiographs. Furthermore, the lack of significant interactions indicates relative independence of these experimental factors. The experimentally confirmed angular independence within a range of IO” is particularly interesting because, combined with the inherent scale independence, it indicates that fractal dimension measured in this way is relatively insensitive to variations in projection geometry that underlie routine radiographic examinations performed in a clinical environment. Inasmuch as such variations are difficult to eliminate in practical situations and tend to limit the precision of current absorptiometric methods for the measurement of osteoporotic changes, this approach provides an attractive alternative to the

status quo. That these results were obtained despite uncontrolled variations in contrast, mean density, and other sensimetric properties of radiographs and common video systems is encouraging and indicates that the expedient of excluding the first few Fourier coefficients from the spectral fit is sufficient. The fact that experimental decalcification produced structural changes that were nearly all in the same direction (Fig. 5) suggests that the developed method might serve as the basisfor a diagnostic test for osteoporosis or other bone loss-related diseases.

However, the analysis shown in Table II and illus-

Fractal dimension of alveolar bone 107

Volume 74 Number I

y = 0.891 - 1.883x. r = 0.993 Estimated Fractal Dimension I 2.58 n fiied l not fiied

Age 63 y I 0.531 - 1.364x , r - 0.990 Estimated Fractal Dimension = 2.82 CI fitted 0 not fitted

t 6 CL B -1

n

Log Power Spectra of Representative Pre- and Post-menopausal Patients 1

-2 !

I 0.2

I 0.4

I

I

log Sp&~i Freqkcy

I 1.0

I 1.2

1

plots of power spectra associated with clinical dental radiographs obtained from typical premenopausal and postmenopausal patients showsregression lines from which respective fractal dimensions were computed. Fig. 7. Representative

trated in Fig. 6 indicates that estimates of D, and thus trabecular bone structure, are location dependent. Although this is not a surprise, it implies that efforts must be made in comparative data analyses to place the data-sampling window over corresponding anatomic regions. In the case of cross-sectional studies, present data provide only coarse guidelines in terms of gross anatomic sites to meet this requirement. However, the higher statistical significance attained with precise registration of the data window compared with only a visually guided placement does indicate the potential importance of a reproducible window-placement procedure for individual patient follow-up. Furthermore, the location dependence of the trabecular structure suggests that multiple anatomic sites should be sampled in future clinical studies. Then, if the pathologic condition to be diagnosed causes global changes of bone strength and structure, an effect on fractal dimension is expected to occur across all sampled sites. The combination of all data in a multivariate statistical approach as illustrated in Table II will then substantially increase the statistical power and thus the specificity of the method. Although this experiment was designed only to induce bone structure changes and not to assess a pos-

sible association between the changes in mass and structure, the lack of correlation displayed in Table II warrants some comments. Whereas, in theory, changes in one without the other are possible (e.g., a structural change can be achieved without a change of mass simply by a spatial redistribution of the mass), this was probably not the case for the in vitro decalcification procedure used. The lack of correlation is most likely due to the fact that the bone region sampled by the data window was only a very low fraction of the bone specimen and the flow of the acid solution was in no way controlled to achieve uniform decalcification; consequently, the integrated calcium loss from the entire specimen was a poor indicator of the calcium lost in the region imaged by the window. However, the contraints in biologic systems impose some interdependence between mass and structure in vivo, which may furthermore change under different pathologic conditions. It may indeed become of clinical relevance to assess the correlation between mass and structural changes that occur under biologic conditions of calcium loss. In vivo data

The significant difference between fractal dimensions computed from the radiographs of premeno-

108

Ruttimann,

Webber, and Hazelrig

pausal and postmenopausal women shown in Table I is most reassuring for several reasons. First, it demonstrates that age-related changes in the trabecular bone pattern produce an increase in fractal dimension, which is in agreement with the expected higher bone resorption in postmenopausal women as well as with the results of the experimentally induced decalcification of the bone specimens. Although the acidinduced bone loss cannot be expected to be a satisfactory model of the in vivo bone demineralization, both processes resulted in loss of spatial correlation in the trabecular structure. This loss of organized fine structure caused more erratic and abrupt density changes in the radiographic images and increased the relative amount of spectral energy in the higher-frequency bands and thus the fractal dimension. Second, the radiographs were retrieved from commonly maintained dental records and thus were made under routine clinical conditions, the data windows applied were of slightly different sizes, and their relative anatomic sites were not rigidly defined. This, coupled with the fact that no data preprocessing other than removal of the average gray level within the windows was applied, supports the premise that fractal analysis may be a robust method to assess structural changes associated with bone demineralization. Third, although the large age difference between the two groups of women studied was deliberately chosen, the realization that statistic significance was attained with such a small sample does suggest a potential sensitivity of the fractal approach to detect biologically relevant, age-related bone structure changes from readily available clinical radiographs. On the other hand, many questions about the ultimate use of this approach remain unanswered. This investigation was limited to observable changes in the alveolar process. Whether these changes are predictive of similar changes in the columnar vertebrae and/or other weight-bearing bones has yet to be determined. Furthermore, whereas the method revealed changes in bone structure between premenopausal and postmenopausal women, the relevance of this information with respect to the development of osteoporosis is not clear. Several uncontrolled factors could be responsible for the observed fractal differences that might have little to do with osteoporosis or age-related bone loss per se. For example, the older women may have had fewer teeth or more nonfunctional teeth, which conceivably might have influenced the local distribution of trabecular bone in some as yet unrecognized fashion. Although one of the purposes of this investigation was to determine whether the fractal model could be applied successfully to peridental bone patterns as recorded in routine dental radiographs, the approach

ORAL

SURC ORAL

MED

ORAL

PATHOL July 1992

taken may not necessarily be limited to this particular imaging technology. To make the method more generally applicable, the impact of the MTF of the radiographic imaging and subsequent image-transduction system on the computed fractal dimension needs to be established. As the observed power spectrum of the image is the product of the power spectrum of the object and the MTF,26 the shape of the latter inlluences the slope of the observed power spectrum and thus fractal dimension. In the case of the clinical radiographs in this study, the effect of the MTF was clearly evident by the limited useful band width of the data (Fig. 7). That band width was similar for the two groups of women (fc in Table I), so that the observed difference of the fractal dimension was not likely due to a systematic difference of the effective MTFs. However, the actual values of the estimated fractal dimensions must be considered biased. It remains to be seen whether, with knowledge of the effective MTF as well as the level of noise corruption, it is possible to back-correct to the desired power spectrum of the object to remove some of this imaging system-induced bias. Clearly, at this stage of the methodologic development, fractal dimensions can be compared only when derived from images obtained with the same imaging methods. If relatively high spatial and contrast resolution images such as those associated with dental nonscreen-type x-ray film turn out to be essential, then the method will necessarily remain limited to images produced outside the mainstream of current medical radiographic resources. SUMMARY

AND CONCLUSIONS

This investigation explored whether the trabecular bone patterns recorded in conventional intraoral dental radiographs can be characterized by a fractal descriptor that may provide a diagnostic indicator of bone mineralization state. In vitro data were produced from 10 mandibular osseous jaw segments radiographed from three controlled projection angles that spanned a range of 10 degrees. Two series of images made before and after demineralization of the bone specimens with 10% formic acid were used. Fractal dimension was estimated from the regression slope of thedouble-logarithmic plot of radially averaged power spectra versus spatial frequency, produced from digitized video images of the radiographs. Repeatedmeasures ANOVA showed that fractal dimension so determined was different for different anatomic locations on the mandibular bone but increased after acid-induced demineralization, irrespective of the projection angles used to produce the radiographs. Computation of fractal dimension by this method applied to intraoral radiographs from two groups of six randomly selected premenopausal and postmeno-

Fractal dimension of alveolar bone 109

Volume 74 Number 1

pausal women showed an increase in the older group. This result and the insensitivity to radiographic projection geometry suggest that this method can be applied to nonstandardized intraoral radiographs and may find potential application in the clinical diagnosis of osteoporotic or other age-related changes in bone structure. We expressour sincerestappreciationto Dr. J.R. Pate1 of the Departmentof DiagnosticSciences,Schoolof Dentistry, University of Alabama at Birmingham,whowasprimarily responsiblefor creating the in vitro database;Dr. Jack Lemmonsof the Departmentof Biomaterials,School of Dentistry, University of Alabama at Birmingham,who assistedin the designof the in vitro experimentand in the implementationof the in vitro assay;andDr. JonathanShip, Clinical Investigationsand Patient Care Branch, National Institute of Dental Research,National Institutesof Health, who provided the intraoral radiographsusedfor the computation of fractal dimensionfrom in vivo data. REFERENCES 1. Riggs BL, Melton LJ III. Involutional osteoporosis. N Engl J Med 1986;314:1676-84. 2. Cummings SR, Kelsey JL, Nevitt MC, O’Dowd KJ. Epidemiology of osteoporosis and osteoporotic fractures. Epidemiol Rev 1985;7:178-208. 3. Holbrook TL, Grazier K, Kelsey JL, Stauffer RN. The frequency and occurrence, impact and cost of selected musculoskeletal conditions in the United States. Chicago: American Academy of Orthopedic Surgeons, 1984. 4. Raisz LG. Local and systemic factors in the pathogenesis of osteoporosis. N Engl J Med 1988;318:818-28. 5. Riggs BL, Walner HW, Seeman E, et al. Changes in bone mineral density of the proximal femur and spine with aging. J Clin Invest 1982;70:716-23. 6. Mazess RB. On aging bone loss. Clin Orthop 1982;165:239-52. 7. Cohn SH. Vaswani A. Zanzi I. Aloia JF Roainski MS. Ellis KJ. Changes in body chemical composition wiyh age measured by total-body neutron activation. Metabolism 1976;25:85-95. 8. Ott SM, Kilcoyne RF, Chesnut CH III. Longitudinal changes in bone mass after one year as measured by different techniques in patients with osteoporosis. Calcif Tissue Int 1986;39:133-8. 9. Ruttimann UE, Webber RL, Schmidt EF. A robust digital method for film contrast correction in subtraction radiography. J Periodont Res 1986;21:486-95.

10. Bolander ME, Gardiner J, Coyle J, et al. Detection and measurement of simulated early rheumatoid lesions of the hand using digital subtraction radiography. Invest Radio1 1990; 25:708-13. 11. Mandelbrot BB. The fractal geometry of nature. New York: WH Freeman, 1983:15, 353-65. 12. Avnir D, Farin D, Pfeifer P. Molecular fractal surfaces. Nature 1984;308:261-3. 13. Burrough PA. Fractal dimensions of landscapes and other environmental data. Nature 1981;294:240-2. 14. Davy Ph, SornetteA, SornetteD. Someconsequences of a proposed fractal nature of continental faulting [Letter]. Nature 1990;348:56-8. 15. Feder J. Fractals. New York: Plenum Press, 1988:184-9. 16. West BJ, Goldberger AL. Physiology in fractal dimensions. Am Scientist 1987;75:354-65. 17. Smith TG Jr, Marks WB, Lange GD, Sheriff WH Jr, Neale EA. A fractal analysis of cell images. J Neurosci Methods 1989; 27:173-80. 18. Cargill EB, Barrett HH, Fiete RD, Ker M, Patton DD, Seeley GW. Fractal physiology and nuclear-medicine scans. SPIE Proc Medical Imaging II 1988;914:355-61. 19. Cargill EB, Donohoe K, Kolodny G, Parker JA, Zimmerman RE. Analysis of lung scans using fractals. SPIE Proc Medical Imaging III 1989;1092:2-9. Caldwell CB, Stapleton SJ, Holdsworth DW, et al. Characterisation of mammographic parenchymal pattern by fractal dimension. Phys Med Biol 1990;35:235-47.. 21. Lundahl T. Ohlev WJ. Kav SM. Siffert R. Fractional Brownian motion: a maximum likelihood estimator and its application to image texture. IEEE Trans Med Im 1986;MI-5: 152-6 1. 22. Kuklinski WS, Chandra K, Ruttimann UE, Webber RL. Application of fractal texture analysis to segmentation of dental radiographs.SPIEProcMedicalImagingIII 1989;1092:11 l-7. 23. Mandelbrot BB, Van Ness JW. Fractional Brownian motions, fractional noises and applications. SIAM Rev 1968;10:422-37. 24. Bracewell RN. The Fourier transform and its applications. 2nd ed. New York: McGraw-Hill, 1986:52-3, 115. 25 Bergland GD. A guided tour of the fast Fourier transform. IEEE Spectrum 1969;6:41-52. 26 Dainty JC, Shaw R. Image science. London: Academic Press, 1974:21 l-5. 37 Dixon WT, ed. BMDP statistical software. Berkeley: University of California Press, 1985. .

L,.

Reprint

requests:

Richard L. Webber, DDS, PhD Departments of Dentistry and Radiology The Bowman Gray School of Medicine Medical Center Blvd. Winston-Salem, NC 27 157-1093

Appendix

If f(x) denotes an attribute of a fractal object (e.g., the height of a mountain range as a function of distance from some base), then f(x) must be statistically self-affine. Hence, if the scale of x is changed by the factor r, then the values of the function change by a factor r H. Formally, for any r > 0 and 0 < H < 1, f(x) and the properly resealed function g(x) = l/rHf(rx) (1) must have identical probability distribution and power spectrum. The power spectrum, lF(s)12,is obtained as the Fourier transform of the autocorrelation function, Rf(u), (Wiener-Khintchin) IF(

= s Rr(u) e-i2Kusdu.

(2)

110

ORAL SURGORAL MED ORAL PATHOL July 1992

Ruttimann, Webber, and Hazelrig

The autocorrelation

function of f(x), estimated from a finite observation interval, L, is given by w

Rr(u;L) = l/L

s

f(x)f(x + u)dx

-L/2 03

= l/L

s

P(x/L)P((x

+ u)/L)f(x)f(x

+ u)dx

= l/L Wc/L)Rr(u),

(3)

where P(x/L) is the rectangular window function of width L, and W(u/L) is the triangular window function that results from the self-convolution of P(x/L).*~ The Fourier transform of the product of functions in equation 3 yields the convolution of the corresponding Fourier-transformed terms /F(s;L)]* = sinc*(Ls)*]F(~)]~,

(4)

where sine*(x) = (sin?rx/?rx)2. If the autocorrelation function in equation 3 exists in the limit for L+oo, equation 4 approaches equation 2, as limL-co [sinc*(Ls)]+(s), the Dirac impulse function. Analogous to equation 3, the autocorrelation of the resealed function g(x) in equation 1 is cc

then

Rs(ru;rL)

= 1/(rZH L)

s

P(rx/rL)P((rx

co = l/(r2H rL) P(t/rL) s

+ ru)/rL)f(rx)f(rx

P((t + ru)/rL)f(t)f(t

+ ru)dx

+ ru)dt

= 1/(r2H+1 LjooW(ru/rL)Rr(ru;rL),

(5)

where in the second integral the substitutions rx = t and dx = dt/r were made. In order to apply equation 2 to equation 5, the substitutions ru = v and du = dv/r are required, yielding the power spectrum of g(x) ]G(s;L)12 = 1/(r2H+1) sinc2(Ls)*]F(s/r)i2.

(6)

For statistic self-affinity, the power spectra in equations 4 and 6 must coincide, and thus IF(

= 1/(rZHfl) ]F(s/r)]*.

(7)

Formally setting in equation 7 s = 1 and r = l/f results in IF(

= ]F(1)]2/f2H + l cx l/l?

Hence, the power spectrum of f(x), when plotted as the logarithm quency, f, must have a negative slope of magnitude b=2H+l.

of the power versus the logarithm

(8)

of the fre(9)

H is called the Hurst exponent and, as discussed by Mandelbrot,” is the fractal codimension of the space spanned by the set of points on f(x), that leads to the following heuristic: A statistically self-affine function fh,x*, . . . ,xn) in E + 1 Euclidian dimensions is a surface with fractal dimension D = E + 1-H = E+(3-b)/2,

(10)

where equation 9 was substituted into the last equation. Hence, with 0 < H < 1, and equations 9 and 10, it follows that 1 < b < 3 and E < D < E + 1. Referring to the example of a mountain range with f(x,y) displaying the altitude at position x,y, the mountain surface has the fractal dimension 2 < D = 3 - H < 3, and the intersectionof avertical planewith that surfacegenerates a self-affine profile with fractaldimension 1 < D = 2 - H < 2.