1.

Introduction

Biophotonic diagnosis of soft tissues using purely and hybrid-optical techniques has witnessed significant advancement during the past three decades and has been recognized as a promising nonionizing tool for the early screening of diseases due to its high resolution and contrast regarding tissue modifications, malignancies, and many other physiological changes.^1,2 Purely optical methods, like optical coherence tomography, diffuse optical tomography, fluorescence and luminescence detections, and Raman and infrared (IR) spectroscopies, give reliable results only for peripheral tissues since they suffer loss of coherence for deep zone imaging due to extreme light scattering.^3–11 Hybrid methods, like photoacoustic tomography/microscopy, ultrasound-modulated optical tomography, and photothermal (PT) imaging, have well-proven advantages in analyzing deep subsurface features with optical contrast.^12–19 Their success is due to the energy conversion process that generates the signal, which results in greatly enhanced dynamic range compared with direct energy methods. Despite the well-known advantages of these approaches for the (early) screening of cancer, cardiovascular diseases, ophthalmic disorders, brain tumors and failure of cerebral functions, and cellular and genetic disorders, which are fundamentally related to soft tissues, bio-optic diagnosis of hard tissues—especially bones—remains a sparsely explored domain. Bone diseases, particularly osteoporosis (OP), have been identified as the second leading health care problem after cardiovascular diseases. Among the World’s population, one in every four women and every eight men over the age of 50 is suffering from OP-related health risks.²⁰ The inefficiency of existing diagnostic techniques for early-stage detection of OP causes significant delays in initiating medical intervention and thus contributes to the severity of the disease in a major way. Established OP diagnostic modalities, like dual-energy x-ray absorptiometry, x-ray micro-computed tomography ($\mu \mathrm{CT}$), mechanical-response tissue analysis, and quantitative ultrasound tomography, have disadvantages such as poor sensitivity or contrast and exposure to ionizing radiation hazards.^21–23

The complex bone geometry consisting of the hard cortical layer with the inner trabecular network and the presence of soft tissue and skin overlayers hinder the adoption of existing biophotonic techniques toward efficient bone diagnosis.²⁴ During the early stages of OP, demineralization mainly occurs in the trabecular region and hence the capability to resolve the trabeculae-generated signals is essential for any modality to satisfy the early detection criterion. The depth range of optical coherence tomography is too small (typically $<1\mathrm{mm}$) to cope with this challenge. Fluorescence and luminescence give a depth-integrated signature, and deconvoluting the contribution of bone trabeculae from the composite signal is a cumbersome process. Raman spectroscopy is sensitive to both organic and mineral phases of bones. However, its advantages are beneficial mostly to in vitro studies because the Raman spectral lines are inherently weak and are strongly scattered by both soft and hard tissues.^25,26 Although diffuse optical tomography yields depth-resolved information, its resolution is inadequate for the trabecular network structure, which usually lies in the submillimeter regime.²⁷ Photoacoustic tomography can image depth-resolved distribution of soft tissue features with high sensitivity and resolution; however, it performs poorly with bone tissues due to the multiple reflections and severe attenuation of ultrasound signals within the trabecular pores, especially at higher frequencies necessary for enhanced spatial and axial resolutions.^13,28–31

Recently, we demonstrated a PT imaging modality called truncated-correlation photothermal coherence tomography (TC-PCT) for the depth-resolved, three-dimensional (3-D) visualization of the thermo-optical features of light absorbing translucent/opaque media.³² Traditional PT techniques which follow the detection of thermal waves generated by modulated optical excitation yield a depth-integrated two-dimensional (2-D) picture of the sample thermophysical features. Using TC-PCT, we have achieved the highest degree of energy localization in a parabolic-wave governed thermal-diffusion-wave field with controllable axial/depth resolution and range. TC-PCT further enables layer-by-layer separation for the independent analysis of various stacks contained in the sample tomogram, reminiscent of a similar optical coherence tomography capability. Furthermore, we have introduced a 3-D PT marker called the thermal wave occupation index (TWOI), derived from TC-PCT, which exhibits an inverse correlation with the degree of bone demineralization or a positive correlation with bone mineral density (BMD). Obviously, TC-PCT analysis of bones has the distinct advantage of computational slicing the cortical and trabecular regions for independent investigations on the status of BMD and, hence, for screening OP in a nonionizing manner. In this article, we investigate the potential of TC-PCT for the detection of bone demineralization in chemically etched goat rib samples in a depth-resolved manner. Since mineral loss is confined to the trabecular region during the initial stage of OP, tracing the trabecular BMD through the cortical overlayer is an advantage in the context of early detection of this disease. We used $\mu \mathrm{CT}$ analysis, assuming it to be the gold standard, for the estimation of a few 2- and 3-D x-ray morphometric parameters with well-known demineralization dependency. The cortical and trabecular sections were analyzed separately and together for investigating the variation of these parameters with demineralization time. The observed correlations between the demineralization time-dependent changes in the $\mu \mathrm{CT}$ parameters and their TC-PCT counterparts validate the potential of the latter for monitoring bone mineral loss in a depth-resolved manner. Furthermore, a sensitivity comparison has been carried out for various $\mu \mathrm{CT}$ and TC-PCT parameters.

2.

Principle of TC-PCT

Contrary to the depth-integrated nature of traditional PT techniques, the recently introduced depth-selective approach involving the chirped-pulsed PT radar exhibits high depth-resolving and penetrating capabilities along with excellent sensitivity.³³ The PT radar is a prominent class of methodologies proposed to achieve the energy localization in parabolic diffusion-wave fields.^19,34 If $R(\tau )$ and $S(t+\tau )$ are the excitation (reference) and response (signal), respectively, the latter being measured at the delay time $t$ and then the cross-correlation,

The pulsed PT radar is capable of acceptable SNR, while using low enough optical pump intensity to meet the laser safety requirement for maximum permissible exposure (MPE) in biological media. The TC-PCT concept was derived from this pulsed PT radar. Although the physical foundation, implementation, features, and potential applications of TC-PCT have been detailed elsewhere,³² a brief discussion is furnished here. A linear frequency-modulated (LFM) train of short optical pulses ($E$) of pulsewidth $W_E$ is used to excite the sample [Fig. 1(a)]. Truncated reference chirps of pulsewidth $W_T(W_T<W_E)$ with continuously variable delay/initial phase ($D=nd$, where $n=0,1,2,\dots$, and $d$ is the delay step) are synthesized ($R_0-R_n$), and each of these chirps is cross-correlated with the PT image captured by an IR camera (radar pulse compression) to produce 2-D axially resolved layer-by-layer images of the sample. For example, cross-correlation of the PT image chirp ($S$) with reference $R_0$, which is in phase with $E$, will give a surface image of the sample. In the case of $R_1$, the resulting image will be that of a near-surface (subsurface) thin layer. As the delay increases, the depth of the layer below the surface also increases. Furthermore, the depth resolution or “thinness” of a truncated-correlation image layer will increase as the truncated pulsewidth shortens. Parameters controlling the maximum penetration depth are the starting LFM frequency and the period/length of the chirp. For a given chirp and truncated pulsewidth, the depth range is controlled by the total delay or the maximum value of $n$. According to the theory of TC-PCT based on conductive PT waves, the truncated-correlation (in-phase/quadrature) output is given by³²

Fig. 1

Signals, algorithm, and instrumentation for TC-PCT. (a) Excitation, photothermal relaxation and typical truncated-correlation chirps (not drawn on the actual time scale); (b) Frequency-domain TC-PCT algorithm implemented in LabView environment. Various modules are fast Fourier transformation (FFT), complex conjugate ($Z^{*}$), inverse fast Fourier transformation (IFFT), multiplication ($X$), and binary exclusive-OR (EX-OR); (c) Instrumental configuration; and (d) Depth-coded planar images sampled at 20-ms reference-delay intervals, obtained from a goat rib.

Here, the various parameters are ${\omega }_1$: Starting chirp angular frequency; $\delta$: Dirac delta function; $t$: Time; $I_0$: Intensity of the excitation pulse; $\xi$: Gaussian characteristic time constant defining the reference pulse shape; $R$: Surface reflectivity of the sample at the excitation wavelength; $K$, $L$, ${\mu }_a$, $\gamma$, $\eta$, and $\alpha$: Sample’s thermal conductivity, thickness, optical absorption coefficient, thermal-wave transmission coefficient at the irradiated surface-air interface, light-to-heat conversion efficiency, and thermal diffusivity, respectively; ${\tau }_t$: An optothermal time constant ($=1/\alpha {\mu }_a^2$), indicating heat conduction from a depth equal to the optical penetration length; and $k=({\omega }_2-{\omega }_1)/T$ is the sweep rate. ${\omega }_2$ and $T$ are the ending angular frequency and period of the chirp, respectively. $n=0,1,2,\dots$, and $p$ is the number of pulses constituting a chirp. erfc refers to the complementary error function, and $\otimes$ represents the convolution operation.

The absolute TC-PCT amplitude ($A_{\mathrm{TC}\text{-}\mathrm{PCT}}$) and phase (${\phi }_{\mathrm{TC}\text{-}\mathrm{PCT}}$) can be estimated as $A_{\mathrm{TC}\text{-}\mathrm{PCT}}=\sqrt{A_0^2+A_{90}^2}$ and ${\phi }_{\mathrm{TC}\text{-}\mathrm{PCT}}={\tan }^{-1}(A_{90}/A_0)$. Furthermore, two additional channels—delay time and zero-crossing phase—have been derived for this tomography. Among these four outputs, the amplitude channel has the key advantage of high-dynamic range and sensitivity based on a one-to-one correspondence with the sample thickness, thus paving the way for a unique depth reconstruction method over about four diffusion lengths. The TC-PCT algorithm has been implemented in the frequency domain for faster execution in LabView environment [Fig. 1(b)]: The recorded excitation chirp is first converted to a square-wave chirp and then passed through the delay-sweep module. The delay-incremented square wave chirp (${\mathrm{C}}_{1f,0}$) is then frequency doubled ($C_{2f,0}$). Subsequently, $C_{1f,0}$ and $C_{2f,0}$ are subjected to the binary exclusive-or (EX-OR) operation to generate a quadrature square chirp $C_{1f,90}$. Truncated references $R_0$ and $R_{90}$ are synthesized from $C_{1f,0}$ and $C_{1f,90}$, respectively. The complex conjugate of the Fourier-transformed truncated references, both in-phase and quadrature, is multiplied with the Fourier transform of the PT chirp, and the products are inverse-Fourier transformed to generate the in-phase and quadrature TC-PCT outputs.

As shown in Fig. 1(c), a semiconductor diode laser capable of generating short optical pulses as triggered by a custom-developed pulse synthesizer is the excitation source. The pulse synthesizer is triggered by a programmable sinusoidal LFM chirp generator (Agilent 33220A, Mississauga, Ontario). The laser source is thermally stabilized using a thermoelectric cooler and a temperature controller. The current through the laser diode is controlled by a laser driver. The laser output is delivered through an optical fiber and is intensity homogenized before irradiating the bone sample. An IR camera (Cedip 520M, Burlington, Ontario, 3.6- to $5.1\text{-}\mu \mathrm{m}$ spectral response) captures the evolution of laser-induced thermal signatures from the bone surface. The entire system is under the control of a computer, which executes the TC-PCT algorithm as well. For its deep-penetrating capability in tissues, an 808-nm (120 W maximum peak power, quasi-continuous-wave) laser diode (Jenoptic JOLD-120-QPXF-2P, Goeschwitzer, Jena) was chosen for excitation. Unless otherwise specified, an excitation pulsewidth of 10 ms was used in our experiments. The starting and ending frequencies for LFM were 0.2 and 0.6 Hz, respectively, with a chirp period/length of 12 s (henceforth, termed chirp-1 waveform). For 50-W peak power and $\sim 3{\mathrm{cm}}^2$ area of illumination, the pulse chirp energy density³³ is $\sim 0.8\mathrm{J}{\mathrm{cm}}^{-2}$. The MPE ceiling for soft tissue is $\sim 3\mathrm{J}{\mathrm{cm}}^{-2}$ and is even higher for bones.³⁶ The truncated chirp pulsewidth was 2 ms, the minimum value limited by the temporal resolution of our camera. This corresponds to an axial/depth resolution of $\sim 25\mu \mathrm{m}$ in bone. It has been observed that the truncated-correlation amplitude offers the highest sensitivity and dynamic range among various output channels and one can generate a sequence of 2-D amplitude images at increasing depths below the surface by incrementing the delay of the reference signal. In Fig. 1(d), a set of such depth-coded planar images (imaged area is $4.0\times 3.2{\mathrm{mm}}^2$) obtained for a goat rib sample ($\sim 50.0\times 15.0\times 4.3{\mathrm{mm}}^3$) is shown. The cortical thickness was $\sim 0.8\mathrm{mm}$, and the trabecular network extended down to $\sim 3.3\mathrm{mm}$ below the irradiated surface.²⁴ The corresponding truncated-correlation delay increment and the approximate depth below the surface are given for every image. These depth-coded slices can be stacked to form a 3-D amplitude image (amplitude tomogram). Furthermore, the amplitude data can be binarized such that those above and below a threshold can be designated as Boolean “1” and “0,” respectively, to depict the physical presence or absence of thermal-wave power at a point within the sample. This binarized amplitude tomogram is independent of the absolute PT amplitude and has advantages in the structural examination of bones, both demineralized and intact.³⁷

3.

Experiment and Results

In traditional x-ray tomography, 2-D morphometric parameters of bone are evaluated from binarized cross-sectional images, whereas surface-rendered volume reconstruction is used for the estimation of 3-D parameters.^38,39 For both 2-D and 3-D, there exist two types of analysis: (1) All objects in the region of interest (ROI) are analyzed together to yield an integrated output. For example, mean total cross-sectional area, cross-sectional thickness, and centroid in 2-D, and object volume, volume of closed pores, and fragmentation index in 3-D. (2) Discrete objects can be categorized, and separate results are generated for individual classes of objects. For example, object area, form factor, and aspect ratio in 2-D, and object volume, volume of pores, and number of pores in 3-D. In addition to the 3-D TWOI, we introduce a set of morphometric protocols by defining a 2-D parameter, as explained below, to evaluate the performance of TC-PCT as a potential bone diagnostic modality. Subsequently, an analogous set of μCT-estimated x-ray morphometric parameters will be used to validate the 2-D and 3-D TC-PCT markers, assuming the former as the gold standard. Throughout this article, we assume the Bone ASBMR (American Society for Bone and Mineral Research) definitions, units, and symbols for the x-ray morphometric parameters:^38,39

In order to formulate TC-PCT morphometric protocols, we made use of the results obtained for a goat rib sample ($\sim 50.0\times 15.0\times 4.3{\mathrm{mm}}^3$). This specimen was prepared by keeping Canadian goat ribs immersed in water for about 3 weeks for a complete decomposition of meat and blood. Clean ribs were rinsed with water and dried in sunlight. For artificial demineralization, the sample was dipped in ethylenediaminetetraacetic acid (EDTA) solution, which was prepared by diluting 0.5-M acid with equal volume of distilled water. The TC-PCT was recorded for the intact state and after 0.5-, 1-, 2.5-, 5-, 10-, and 20-h demineralization in EDTA. Figure 2 shows 2-D amplitude images and 3-D binarized stack sequences, both at different depths, as well as absolute amplitude and binarized full tomograms, all for the rib piece before and after 1-h demineralization. Images in Figs. 2(a1–10) and 2(c1–10) are the binarized TC-PCT amplitude slices ($\sim 300\text{-}\mu \mathrm{m}$ thick), vertically aligned as a function of depth ($z$-coordinate) below the irradiated surface before and after demineralization, respectively. Each slice is formed by stacking and binarizing 12 depth-resolved planar images, each of thickness $\sim 25\mu \mathrm{m}$. Figures 2(b1–10) and 2(d1–10) are the planar images sampled at each 24-ms delay interval before and after demineralization, respectively. They depict the distribution of absolute TC-PCT amplitude in the specified plane. In the figure, BT1 and BT2 are, respectively, the binarized amplitude tomogram for a depth of $\sim 3\mathrm{mm}$ below the surface for the rib sample before and after demineralization. AT1 and AT2 are the corresponding absolute amplitude tomograms before and after demineralization for the same depth. BT1 and BT2 are generated from AT1 and AT2, respectively, through binarization followed by surface-mesh formation for volume rendering. The imaged area is $4.0\times 3.2{\mathrm{mm}}^2$. The experiment was performed using chirp-1 excitation with $\sim 50\text{-}\mathrm{W}$ peak power. Both truncated-reference pulsewidth and delay-step were 2 ms. Maximum delay applied was 240 ms, which corresponds to a depth of $\sim 3\mathrm{mm}$ in bone. For the binarized amplitude tomograms (BT1 and BT2), the 3-D morphometric parameter TWOI measures the fractional volume occupied by the photothermally active bone components contained in a cube element (cuboid):

Fig. 2

Binarized amplitude segments (a, c) and layer-by-layer planar TC-PCT amplitude images (b, d) obtained for the intact and 1-h demineralized conditions of a goat rib. The 10 binarized sequences, $\sim 300\text{-}\mu \mathrm{m}$ thick each, compose of tomograms BT1 and BT2. The 120 planar sequences, $\sim 25\text{-}\mu \mathrm{m}$ thick each, constitute of the absolute amplitude tomograms AT1 and AT2. Demineralization-induced reduction in the binarized volume can be seen for c1–10 compared with its intact state depicted in a1–10. The mean PT signal for every plane in d1–10 is less than that for its counterpart in b1–10.

To investigate the potential of these markers for monitoring the degree of demineralization, we considered four goat rib samples, labeled as RB10 to RB13, and demineralized in EDTA solution for 10, 20, 30 and 40 h, respectively. They were subjected to chirp-1 TC-PCT analysis before and after demineralization. The results obtained for MPA before and after demineralization are shown in Fig. 3. A general trend observed in Fig. 3 for MPA is the progressive decrease as the depth below the surface increases. This is attributed to the attenuation of optical intensity with depth in the sample, due to absorption and scattering. In contrast to this steady behavior, there appears an anomalous response near the surface. One potential reason for this is the presence of peripheral contamination and its random reactions with EDTA. In addition, the temporal uncertainty in the camera frame trigger ($\sim 2\mathrm{ms}$) is likely to cause a lapse in synchronizing the first frame with the start of the thermal relaxation. This may affect the surface/near-surface data contained in the early transient. Another possibility is the thermal confinement near the surface, which tends to raise the mean temperature in the very near subsurface region. Possible mechanisms responsible for the changes in MPA following demineralization are presented in Sec. 4. Figure 4 shows the variation of TWOI with the duration of demineralization for cortical and trabecular regions independently and together. Four ROIs ($4.0\times 3.2{\mathrm{mm}}^2$) were selected randomly on the sample surface for averaging. For individual analysis, the cortical and trabecular regions were reconstructed separately by suitably selecting the number of depth-resolved planar images, and the ratio of sequences being determined by visually examining the onset of the trabecular signature. The common trend observed is that the TWOI follows an inverse correlation with demineralization time for 10 and 20 h. For very high demineralization (30 and 40 h), the response is unpredictable. The significance of these trends and the mechanisms responsible for them will be considered in detail in Sec. 4. In Figs. 5(a)–5(c), binarized amplitude tomograms of intact goat rib sample RB10 are shown. They were used for TWOI estimation and were reconstructed individually and together for cortical and trabecular regions. The corresponding absolute TC-PCT amplitude tomograms are shown in Figs. 5(d)–5(f). In Figs. 5(a) and 5(d), $z=0$ corresponds to the irradiated cortical surface. These samples were further analyzed using a $\mu \mathrm{CT}$ scanner (Skyscan 1173, Bruker, Kontich, Brussels) over the same ROI before and after demineralization.

Fig. 3

Variation of MPA with depth below the irradiated sample surface before and after demineralization for 10, 20, 30, and 40 h. For low-demineralization levels, usually up to 10 to 20 h for a goat rib in 50% EDTA solution, the PT signal decreases with the time of acid treatment. For longer durations, the response is unpredictable due to the trabecular collapse and bone volume compression, especially along the thickness.

Fig. 4

Variation of TWOI with duration of demineralization estimated separately (a and b) and together (c) for the cortical and trabecular regions. As in the case of MPA variations, the TWOI response maintains a negative correlation with the degree of mineral loss for low-demineralization levels. It becomes erratic for very strong mineral loss. The decrease in the absorption and scattering coefficients with demineralization is the main reason for the fall in TWOI for 10- and 20-h samples. For chirp-1 excitation, uncertainty in the estimation of TWOI is $\sim 1.3\%$.

Fig. 5

Binarized amplitude TC-PCT images (a–c) and absolute amplitude tomograms (d–f) for the intact sample RB10. (a and d) Cortical alone; (b and e) Trabecular alone; and (c and f) Cortical-trabecular combined. Imaged area is $4.0\times 3.2{\mathrm{mm}}^2$ and depth is $\sim 2.8\mathrm{mm}$. For tomograms (a, d, c, and f), the plane $z=0$ coincides with the irradiated cortical surface. For images (b and e), $z=0$ lies on the cortical-trabecular interface. Tomogram in (c) is the view from the trabecular bottom toward the irradiated cortical layer, which can be used to assess the resolution of TC-PCT. Geometric measurements show that the trabeculae thinner than 100 μm have been resolved.

Since $\mu \mathrm{CT}$ parameters follow well-established correspondence with the extent of mineral loss, their primary role in this study is to understand the dependence of TC-PCT markers on the demineralization process. Although TC-PCT is at a very early stage of development compared with $\mu \mathrm{CT}$, the tandem investigations using these complementary modalities are helpful to establish clear distinctions between their strengths and weaknesses in the context of tracking demineralization. A parallel investigation using TC-PCT and $\mu \mathrm{CT}$ was carried out for a broad range of mineral loss (up to 40 h in EDTA, as already mentioned) that is much larger than the natural range observed in osteoporotic cases. However, the number of samples was restricted to four (10-, 20-, 30-, and 40-h demineralization) due to the limited availability of the $\mu \mathrm{CT}$ analyzer during this work. For this reason, an assessment of the absolute sensitivity of the two methods was not carried out. All the samples were subjected to TC-PCT (Toronto) and then $\mu \mathrm{CT}$ (San Antonio) analysis in their intact state. They were then demineralized in separate containers such that sample 1 (RB10), sample 2 (RB11), sample 3 (RB12), and sample 4 (RB13) remained in EDTA for 10, 20, 30, and 40 h, respectively. TC-PCT and $\mu \mathrm{CT}$ were recorded for the demineralized samples also. The ideal approach would have been to demineralize a sample sequentially (say, for $0,0.5,1.0,1.5,\dots \mathrm{h}$) and to record its TC-PCT and $\mu \mathrm{CT}$ after each demineralization stage, but that was not possible due to the physical distance of the two facilities. A statistical analysis on the behavior of a set of such sequentially demineralized samples could have led to a more detailed understanding of the interplay among various bone parameters responsible for the observed signal variations. Leaving this as a future initiative, a large set of 2-D and 3-D x-ray morphometric parameters of samples RB10 to RB13 was estimated for cortical and trabecular regions separately and together as the reference markers for validating the TC-PCT results. Boundaries for different ROIs were selected as shown in Figs. 6(a)–6(c): for purely cortical region, the irradiated surface and cortical-trabecular interface were the boundaries (a). Cortical-trabecular (upper) and trabecular-cortical (lower) interfaces defined the periphery for purely trabecular analysis (b). Irradiated cortical surface and lower trabecular-cortical interface formed the border for combined cortical-trabecular measurements (c). In all cases, regions with bluish highlight are the bone voxels of interest, and the red-filled spaces are nonbone voxels (air) excluded from measurements. In Fig. 7, the pre- and post-demineralized $\mu \mathrm{CT}$ images of slices normal to the bone thickness (a1–d1 and a3–d3) and along the length (a2–d2 and a4–d4) are shown. With the decrease in mineral content, the image contrast falls as x-ray attenuation in the organic phase is much lower than that in the mineral phase.

Fig. 6

Region of interest (ROI) boundary selection for $\mu \mathrm{CT}$ analysis. Regions with bluish highlight are the bone voxels of interest, and the red-filled spaces are nonbone voxels (air) excluded from the measurements. (a) For purely cortical region, the irradiated surface and cortical-trabecular interface were the boundaries. (b) Cortical-trabecular (upper) and trabecular-cortical (lower) interfaces defined the periphery for purely trabecular analysis. (c) Irradiated cortical surface and lower trabecular-cortical interface formed the border for combined cortical-trabecular measurements.

Fig. 7

Planar $\mu \mathrm{CT}$ images of samples RB10to RB13 before and after demineralization: (a1–d1and a3–d3) Plane normal to the bone thickness before and after demineralization, respectively. (a2–d2 and a4–d4) Plane along the length for the intact and demineralized conditions, respectively. Even for demineralization time as low as 10 h, considerable mineral loss occurred at the trabecular region. This has greatly reduced the trabecular contrast in post-demineralized images. It can also be seen that very little contrast change occurs for this region, as one moves through a3–b3–c3 or a4–b4–c4. Perceptible trabecular mineral traces exist in d3 and d4. This can be either due to the natural immunity of this sample to demineralization or weak EDTA penetration into the trabecular network. However, this mineral residue does not make any significant change to the general demineralization response of the trabeculae.

In Fig. 8, pairs of typical 2-D and 3-D $\mu \mathrm{CT}$ parameters are depicted for separate and combined cortical and trabecular regions. In general, these parameters maintain an inverse correlation with the degree of demineralization time up to 30 h. However, an increase compared with the result for 30-h is observed for the 40-h sample. Due to a failure in $\mu \mathrm{CT}$ protocol, a pair of trabecular thickness (Tb.Th) data (20 and 30 h) is missing.

Fig. 8

Various 2-D and 3-D $\mu \mathrm{CT}$ morphometric parameters (left and right columns, respectively) estimated for the rib bones, RB10–13, before and after demineralization. These parameters were evaluated separately and together for the cortical and trabecular regions. Since $\mu \mathrm{CT}$ estimates the absolute value of these parameters, their variation with demineralization time can be used to validate the potential of TC-PCT markers for monitoring mineral loss (see text for details). The uncertainty in the measurement of all $\mu \mathrm{CT}$ parameters is $<1\%$.

4.

Discussion

4.3.

Demineralization Dependency of $\mu \mathit{CT}$ and TC-PCT Parameters

The ideal procedure for understanding the morphometric variations with mineral loss is to sequentially demineralize a sample by incrementing the time of EDTA treatment and to record its TC-PCT and $\mu \mathrm{CT}$ after each treatment. As already mentioned (Sec. 3), due to the physical distance between these facilities, we considered only four samples with different EDTA treatment times (10, 20, 30, and 40 h) and their independent analysis. It was intentional to select a broad range of this kind (intact to 40 h) with a reasonably long separation (10 h): due to the highly inhomogeneous nature, bone properties (mineral distribution, trabecular dimensions, cortical thickness, surface area, density, etc.) may vary from sample to sample even if they belong to the same rib. So, if a shorter range (say, 1, 2, 3, and 4 h) of demineralization is considered, then the resulting changes for consecutive demineralization times need not remain resolvable always, i.e., overlap may occur. For the logical analysis of the demineralization dependency of various parameters (both TC-PCT and $\mu \mathrm{CT}$) based on this multiple sample approach, a normalization process has been carried out for every parameter by defining a fractional change $[\delta (\text{parameter})]$ for the intact and demineralized states:

Eq. (7)

$$\delta (\text{Parameter})=\frac{{\text{Parameter}}_{\text{Intact}}-{\text{Parameter}}_{\text{Demineralized}}}{{\text{Parameter}}_{\text{Intact}}}.$$

4.3.1.

Demineralization dependency of 2-D parameters: T.Ar, Tb.Th, and MPA

The T.Ar is the mean of the cross-sectional area for all slices in the selected volume of interest. The ROI periphery assumed here is the outer boundary of the object, without including any internal spaces or structure in the measurement. For the goat rib sample, the outer surface and the cortical-trabecular interface are considered as the ROI periphery for evaluating the cortical T.Ar [Fig. 6(a)]. Any pores or inclusions within the cortical slab will be discarded. In general, the cortical T.Ar reduces on demineralization, as shown in Fig. 8(a). Up to 30 h, the decrease in this parameter follows a positive correlation with the demineralization time. This trend is natural as the removal of mineral would reduce the volume of bone material, and all non-bone voxels including newly created voids at the locations of mineral loss are excluded from T.Ar estimation. This in turn reduces the cross-sectional area. But as one moves from the 30- to 40-h sample, the fractional change in this parameter decreases from 0.51 to 0.37. Corresponding $\mu \mathrm{CT}$ images, Figs. (7c2 and c4 and d2 and d4, respectively), support this measurement as the top cortical thickness in Fig. 7(d4) appears as a larger fraction of its intact state compared with Fig. 7(c4). This can be due to the higher biochemical immunity of the 40-h sample to demineralization compared with its 10-, 20-, and 30-h counterparts. Another possibility is the severe trabecular collapse on 40-h EDTA treatment. In this case, thin trabecular whiskers collapsing on volume compression following very long demineralization may remain crowded near the cortical-trabecular interface. This is likely to happen since compression occurs predominantly along the bone thickness. Obvious increase in the $\mu \mathrm{CT}$-estimated cortical thickness and hence T.Ar is thus possible. A nearly comparable behavior also is observed for the combined cortical-trabecular T.Ar [Fig. 8(e)]. Unlike cortical-alone results, T.Ar does not differ substantially for 10- and 20-h demineralizations. For 30- and 40-h samples, the response is similar to that of the cortical T.Ar. For a better understanding of the trabecular T.Ar response, we analyzed the Tb.Th, another trabecular 2-D parameter [Fig. 8(c)]. The estimation of this parameter for 20- and 30-h EDTA treatments was not possible due to a failure in the $\mu \mathrm{CT}$ protocol. Similar breakdown was observed for other cortical 2-D parameters like trabecular number and trabecular separation (not accounted in this article). This kind of protocol failure may arise while analyzing samples of broad BMD variations, since the definitions made for locating the ROI periphery may not be equally valid for low- and high-demineralization levels. These results lead to a conclusion that the 2-D $\mu \mathrm{CT}$ parameters need not always yield reliable morphometric information, especially for the trabecular region, if the demineralization range is broad enough to cause the protocol malfunction. Hence, a truthful validation of the 2-D TC-PCT parameter (MPA) response may not be possible using 2-D $\mu \mathrm{CT}$ parameters if the range of demineralization is very broad.

As a consequence of optical attenuation resulting from absorption and scattering in bone, the MPA decays with the depth below the irradiated cortical surface for all demineralization conditions, as seen in Fig. 3. In Figs. 9(a)–9(c), the mean MPAs estimated using the data constituting Fig. 3 for the prominently cortical (0.3 to 1.2 mm), purely trabecular (1.2 to 3 mm), and combined cortical-trabecular (0.3 to 3 mm) regions have been shown for the intact and demineralized conditions of all four samples. Owing to uncertainty mentioned in Sec. 3, near-surface (0 to 0.3 mm) data have been omitted in all cases. For all three regions, the post-demineralized MPA maintains an inverse correlation [positive correlation for $\delta (\mathrm{MPA})$] with mineral loss for 10- and 20-h samples. This feature can be attributed to the fall in absorption and scattering coefficients with the degree of mineral loss. Another possibility is the increase in bone thermal diffusivity due to the decrease in density on mineral loss. Since diffusivity is a measure of the rate of heat flow/removal from the point of generation, it would be reasonable to assume that a drop in PT amplitude is possible with the loss in density. Since the exact variations of specific heat capacity and thermal conductivity with demineralization are not known, the role of diffusivity and volumetric specific heat cannot be precisely interpreted at this stage. For the prominently cortical region, the post-demineralization MPA exceeds its intact value by about 1% and 7%, respectively, for the 30- and 40-h samples. Since absorption and scattering coefficients keep on decreasing with mineral loss, thus reducing the TC-PCT amplitude and hence MPA, it is logical to assume that one or more mechanisms capable of boosting PT signal compensating the influence of falling optical coefficients come into play for very high demineralization levels. Both the PT amplitude [Eq. (6)] and diffusivity follow an inverse correlation with $\rho c$, the volumetric specific heat capacity. If there exists a decrease in thermal conductivity with mineral loss at a rate larger than the rate of descend in volumetric specific heat, then a strong signal enhancement may occur. The possibility of such a substantial reduction in thermal conductivity can be justified using the void fraction hypothesis described in Sec. 3: excessive removal of bone material on prolonged demineralization may enlarge the void fraction within both the cortical and trabecular regions. Given the poor thermal conductivity of air (assuming that the pores are air filled), the increased void volume may have impeded the conductive loss rate, thereby sustaining a stronger thermal-wave field within the cortical layer and trabecular network. Generally, in all three cases (cortical and trabecular regions independently and together), $\delta (\mathrm{MPA})$ follows comparable variations with the time of EDTA treatment. The highest impact of demineralization, as estimated by $\delta (\mathrm{MPA})$, is on the purely trabecular section and the lowest on the prominently cortical segment. This is possible, since, for a given bone volume, the trabecular region maintains a larger contact area with EDTA compared with the cortical layer.

Fig. 9

Mean MPAs, estimated using the data constituting Fig. 3, for the prominently cortical (a), purely trabecular (b) and combined cortical-trabecular (c) regions. Near-surface (0 to 0.3 mm) data have been omitted in all cases due to the uncertainties arising from synchronization issues, surface contamination, and near-surface thermal confinement, which is randomly influenced by the ambient air currents. For all three regions, the post-demineralized MPA maintains an inverse correlation with mineral loss for 10- and 20-h samples, the same trend observed for TWOI (Fig. 4). Uncertainty in the estimation of MPA is $\sim 0.9\%$ for chirp-1 excitation.

4.3.2.

Demineralization dependency of 3-D parameters: TMD, BMD, and TWOI

By definition, TMD is calculated from the average x-ray attenuation in the bone tissue only without including the contributions of nonbone voxels like biofluids and soft tissue. On the other hand, BMD estimation accounts for the collective attenuation offered by both the bone and nonbone voxels.^38,39,46 Since our rib bone samples are deprived of nonbone bioconstituents, both TMD and BMD estimate the attenuation by bone component only. The TMD estimation is considered for the purely cortical and combined cortical-trabecular sections because of the dominance of bone volume fraction over the air volume fraction. The TMD, as shown in Fig. 8(b), trails a monotonic increase for its fractional change $[\delta (\mathrm{TMD})]$ with demineralization time up to 30 h. For the 40-h sample, $\delta (\mathrm{TMD})$ makes a significant drop (47%) from its value for 30-h demineralization (79%). A monotonic increase in the fractional change or $\delta (\text{parameter})$ up to 30-h and then a decrease for 40-h demineralization are a characteristic feature of all three sections. The continuing increase in $\delta$ (TMD and BMD) for samples RB10 to RB12 is an immediate consequence of the fall in mineral content with EDTA treatment time. However, for RB13, very high demineralization may lead to significant bone volume compression. This can result in a higher value for the estimated mineral density compared with what is expected under a volume conserved demineralization process.

The 3-D TC-PCT marker, TWOI, yields a marginally increased $\delta$ value between 10 and 20 h [Fig. 4(a)]. For 30 h, the trend is toward a decrease but with high uncertainty, which is also true for 40 h. For purely trabecular and cortical-trabecular combination, similar responses are observed for TWOI [Figs. 4(b) and 4(c)]. As in the case of MPA, the justifications based on the reduction in optical absorption and scattering coefficients with mineral loss are applicable here also for the 10- and 20-h samples. For 30- and 40-h EDTA treatment, the void fraction increase and resulting fall in thermal conductivity should be responsible for reversing the trend.

4.4.

Relative Sensitivity of $\mu \mathit{CT}$ and TC-PCT Parameters

For both $\mu \mathrm{CT}$ and TC-PCT, there exists a common characteristic: in the case of all parameters, the variation of $\delta$ with demineralization time follows a change in direction at a stage that is early for TC-PCT compared with $\mu \mathrm{CT}$. Due to the highly inhomogeneous nature of bones, this change-over time should vary from sample to sample and be dependent on the demineralization conditions. Since naturally occurring mineral loss may not be as high as what was observed for very long EDTA treatment (say, 30 to 40 h), we confine our major interest to the first two responses as far as OP diagnosis is considered as the prime target. The broad range variation will be significant while developing the general purpose TC-PCT protocols and allied calibration procedures. Assessment of the relative sensitivity of various parameters to the changes in cortical and trabecular sections, independently and together, is important prior to arriving at vital conclusions of this investigation. For this, we assumed the 10-h $\delta$ value as the baseline/offset, for all parameters, and made it “0.” Such a baseline subtraction, although not strictly accurate to scale a nonlinear response, is helpful to judge the relative capability of different parameters to trace the effects of EDTA treatment. The offset-corrected $\delta$ parameters are plotted against respective demoralization time in Figs. 10(a)–10(c). For sensitivity comparison, we consider 10- and 20-h demineralized samples, the practically interesting cases. For the prominently cortical layer, $\mu \mathrm{CT}$ parameters are more sensitive than TC-PCT markers. Since a clear idea is not yet available on the variation of different thermophysical parameters (especially, specific heat capacity and thermal conductivity) with cortical demineralization, no vital conclusion can be made on the inferior TC-PCT sensitivity to cortical mineral loss. For purely trabecular and cortical-trabecular combination, TC-PCT offers better sensitivity than $\mu \mathrm{CT}$. One can see in Fig. 7 that the $\mu \mathrm{CT}$ trabecular contrast is very low even for the least demineralized sample (RB10) and very little contrast change is observable as one moves to higher levels of mineral loss (RB11 to RB13) for this porous region. X-ray attenuation is mainly due to the mineral phase, most of which are removed during the 10-h demineralization itself. Also, one can see that the absolute trabecular mineral density is much smaller than the cortical value for all samples [Figs. 8(b) and 8(d)]. The poor $\mu \mathrm{CT}$ sensitivity to trabecular mineral loss is a dynamic range limitation, in which the large changes in cortical mineral content override minor trabecular changes. Unlike μCT, the TC-PCT sensitivity to changes in both mineral and organic phases is a benefit in the case of trabecular analysis. Since osteoporotic demineralization first occurs at the trabecular region, the enhanced trabecular sensitivity of TC-PCT is a remarkable advantage in the early detection scenario. For a better picture on the dynamic range features of these modalities, we analyzed the binarized amplitude TC-PCT and $\mu \mathrm{CT}$ images of the 20-h treated rib before and after demineralization (Fig. 11). The tomogram is a view from below, i.e., from the trabecular surface toward the cortical layer. It can be seen that the TC-PCT offers much better dynamic range which allows the trabecular structure of the demineralized sample to be clearly visible along with the cortical portion, justifying our previous argument.

Fig. 10

Sensitivity of various $\mu \mathrm{CT}$ and TC-PCT parameters to (a) purely cortical, (b) purely trabecular, and (c) combined cortical–trabecular demineralizations. Assuming as the baseline/offset, the 10-h $\delta$ value was subtracted from corresponding higher demineralization values for all four regions; a coarse calibration procedure. For purely trabecular and combined cortical-trabecular regions, TC-PCT offers better sensitivity than $\mu \mathrm{CT}$. For the cortical-dominant region, $\mu \mathrm{CT}$ is more sensitive. Since osteoporotic demineralization first occurs at the trabecular region, the enhanced trabecular sensitivity of TC-PCT is a remarkable advantage in the early detection scenario. Trabecular $\delta$ (T.Ar) values have been discarded due to a failure in the measurement protocol. The legend is common for all three plots.

Fig. 11

Binarized amplitude TC-PCT and $\mu \mathrm{CT}$ images of the goat rib sample before and after demineralization. $\mu \mathrm{CT}$ images of the whole sample ($\sim 38\text{-}\mathrm{mm}$ long) are shown in the central column. The ROI for TC-PCT ($4.0\times 3.2{\mathrm{mm}}^2$) is outlined in green color. The TC-PCT image and zoomed-in version of $\mu \mathrm{CT}$ image corresponding to the ROI are shown in the left and right columns, respectively. Due to the high-dynamic range advantage of TC-PCT, both cortical and trabecular regions are resolved with high contrast for the post-demineralized sample.

X-ray based methods measure the absorption/attenuation of photons (mostly) by the mineral content in bones, and tomographic reconstruction is done using back-projected planar images. Since scattering and diffraction effects are negligible for bone x-rays, a well-defined structural imaging is possible. In the case of TC-PCT, depth-coded planar images primarily depict the (photo)thermal-wave distribution associated with optical absorbers. Not only the absorber distribution, but also the spatial variation of thermophysical parameters (e.g., density, effusivity, etc.) influences the thermal-wave distribution within and around an object. This means that the TC-PCT generates a functional image of the thermo-optical features of the sample. If the separation between adjacent objects (say, trabeculae) is large enough to be axially depth-resolved locally without overlapping thermal-wave contributions from other depths, then structural imaging is possible with TC-PCT on a spatial resolution scale of the thermal wavelength. Since the thermal diffusion length is inversely correlated with the square root of frequency, higher resolution for improved structural imaging would be possible at higher frequencies that demand shorter excitation pulsewidth. However, this will limit the depth range. With chirp-1 excitation, the lateral resolution observed in bone is $\sim 100\mathrm{\mu m}$.

5.

Conclusions

The current investigation has been focused on the nonionizing functional imaging and diagnosis of artificial demineralization in animal bones. Specific attention has been paid for resolving the trabecular structure through the cortical layer, quantifying the extent of mineral loss using the recently introduced TC-PCT technique, and validating the results using $\mu \mathrm{CT}$ imaging. The TC-PCT computational slicing advantage has enabled the independent analysis of cortical and trabecular regions and the comparison of these results with the respective $\mu \mathrm{CT}$ images. Analogous to x-ray morphometric parameters, a 3-D marker (TWOI) and a new 2-D parameter (MPA), both derived from the TC-PCT output, have been defined and estimated for goat rib bones demineralized at intervals of 10 to 40 h using EDTA solution. For shorter intervals, 10 and 20 h, both TC-PCT and $\mu \mathrm{CT}$ markers follow the same trend. This is consistent with our previous report,⁴⁵ which investigated the early mineral loss detection capability of TC-PCT. For higher extent of demineralization, 30 and 40 h, both the $\mu \mathrm{CT}$ and TC-PCT results show anomalous behavior. This is not a matter of prime concern, since such severe demineralization is unlikely to occur naturally in living bones. A direct comparison of tomograms derived using TC-PCT and $\mu \mathrm{CT}$ has revealed the superior dynamic range of the former. As a consequence, TC-PCT offers better sensitivity to trabecular and combined cortical-trabecular demineralizations compared with $\mu \mathrm{CT}$. However, comparison of TC-PCT with other modalities with respect to the structural imaging must be made carefully: depending on the resolution of the IR camera, distribution of thermophysical properties of the sample, and the frequency bandwidth of the TC-PCT image, there will be a limit on the possible lateral resolution ($\sim 100\mu \mathrm{m}$ at present) for a given threshold of binary TC-PCT reconstruction. This limit clearly characterizes the conditions under which purely structural imaging is possible for bone using this modality. As far as the application of TC-PCT for in vivo bone diagnosis is concerned, an immediate possibility is the case of femur (or similar bones) of small animals, like mice, for which the soft tissue overlayer thickness is not large (1 to 2 mm).