Characterization of multi-biomarkers for bone health assessment based on photoacoustic physicochemical analysis method.

Photoacoustic (PA) techniques are potential alternatives to histopathology. The physicochemical spectrogram (PCS) generated by the PA measurement at multiple wavelengths can presents the morphology and chemical composition target at multi-biomarkers simultaneously. In this work, via multi-wavelength PA measurements performed on rabbit bone models, we investigated the feasibility of using PCSs for bone health assessment. A comprehensive analysis of the PCSs, termed PA physicochemical analysis (PAPCA), was conducted. The "slope" and "relative content" were used as the PAPCA-quantified parameters to characterize the changes in the physical and chemical properties of bone tissue, respectively. The findings are consistent well with the gold-standard imaging results. It demonstrated that the PAPCA can be used to characterize both the microstructure and content of multi-biomarkers which highly related with bone health. Considering the PA technique is noninvasive and radiation-free, it has great potential in the implementation and monitoring of bone diseases progression.

Introduction

Osteoporosis is a major risk factor for fragility fractures. It is characterized by a decrease in bone mineral density (BMD), alterations in bone microarchitecture (BMA), and changes in chemical composition of the metabolic skeletal system. Recent diagnostic methods are generally based on X-rays or ultrasound [1], [2], [3]. The dual-energy X-ray absorptiometry (DEXA) is the “gold standard” for clinical diagnosis of osteoporosis, which can provide the BMD information. However, it cannot measure the BMA, bone elasticity, or other important factors that determine the risk of fractures. In the quantitative ultrasound (QUS) bone assessment method, ultrasonic attenuation and speed of sound (SOS) in calcaneus (heel) bone in vivo has a strong correlation with bone mineral content, which is an important determinant of fracture risk [3], [4], [5]. Besides, ultrasonic backscatter measurements can provide the QUS parameters having significant correlations with bone density and structural features [6], [7], [8]. However, these techniques mostly focus on the physical information of the bone, including BMD or BMA, and they are limited when pathogenic bone diseases are determined based on chemical changes.

Bone health depends not only on the bone mass and structure of the non-organic mineral matrix but also on the composition of the organic matrix, including lipids, blood cells, collagen, and proteins. For example, the blood content and oxygen saturation of the marrow is highly associated with bone blood flow and cellular metabolism. In addition, previous studies have shown that, based on clinical observations, fat filled in the bone marrow is not a simple filling tissue, but is also considered an actor within the bone microenvironment; it is involved in bone metabolism and bone conversion and is highly associated with weak bone mass in osteoporosis [39]. Magnetic resonance imaging (MRI) and magnetic resonance spectroscopy have shown potential as noninvasive imaging techniques for quantifying the chemical composition of bone, such as the lipids and blood perfusion [9], [10]. In addition, optical spectroscopic techniques have been used to evaluate the contribution of bone composition alterations to bone quality changes related to aging, disease, or injury [11]. However, MRIs can only provide limited chemical information and are bulky as well as expensive, whereas traditional optical techniques suffer from limited spatial resolution and overwhelming optical scattering in biological tissues, thereby reducing the efficacy of in vivo skeletal imaging and hard to provide the information related with BMA. Moreover, the findings of current bone imaging and detection technologies make it difficult to quantify the physical microstructure and chemical composition information of different compositions simultaneously. To overcome this issue, easy-to-apply biomarkers and new diagnostic approaches are required.

The biomedical photoacoustic (PA) technique, which is based on light and sound, has the unique ability to detect highly sensitive light absorption contrast in deep biological tissues [12], [13], [14], [15], [16], [17]. The PA technique has the potential to measure not only the parameters related to BMD and BMA, but also the content and distribution of chemical constituents, such as collagens, hemoglobin, and lipids, all of which are highly relevant to bone health. Lashkari et al. assessed the changes in bone structure, density, and collagen content using dual backscatter ultrasound and PA radar systems [18], [19]. Gu et al. used photoacoustic Fourier transform infrared spectroscopy to measure the mineral stoichiometry of cortical bone [20]. Wood et al. used an optimized PA imaging technique to assess oxygen saturation in the bone marrow cavity through disease progression in a murine model of acute lymphoblastic leukemia [21]. Steinberg et al. used a dual-modality multispectral PA system to quantify the blood/fat ratio in the marrow in vivo, which has been correlated with molecular changes in the long bone. In addition, our group investigated the feasibility of using thermal PA, PA spectral analysis (PASA), and PA time-frequency analysis to evaluate the BMD and BMA of trabecular bone in rat models [22], [23], [24]. Recently, our group evaluated the relative content of collagen in rabbit bone model ex vivo by performing multi-wavelength PA (MWPA) measurements at the wavelength of 1300 nm–1800 nm [25]. Feng et al. also investigated the feasibility of PA assessment of human bone health in vivo by performing MWPA measurements [26], [27], [28]. Since the PA signals represent the spatial distribution of the optical absorption properties in biological tissue, past studies have indicated that this technique combines micrometer to centimeter morphologies and the chemical compositions of tissue specimens [29], [30], [31], [32], [33], [34], [35]. For example, based on the PA measurements for each examined tissue specimen, a physicochemical spectrogram (PCS) that integrates the power spectra of the PA signals along the full optical spectrum can be formulated. The PCS presents a unique physicochemical signature of any type of tissue [33].

To date, although the use of PA physicochemical analysis (PAPCA) of PCSs has been reported for visualizing tissue microstructures and chemical composition simultaneously, it has never been applied in bone assessment to the best of our knowledge. In this study, we investigated the feasibility of the PAPCA technique in assessing the variations of the BMD, BMA, and chemical composition in bone through ex vivo experimental studies on rabbit bone models. Besides, in order to targeted at both the bone mineral matrix (mostly hydroxyapatite) and bone marrow (mostly blood and lipid) [22], [31], [32], we chose the laser wavelength at the range of 690 – 950 nm. The PCSs of the bones from the osteoporosis and control groups were obtained and compared. Subsequently, the quantified PAPCA parameters were obtained and compared with gold-standard DEXA, micro-computed tomography (CT), and MRI images, which present the BMD, BMA, and lipid ratio of bone tissue, respectively. Finally, the quantified parameters of the chemical and physical information were classified using support vector machine (SVM) learning.

Materials and methods

Animal disease models

In this study, well-established rabbit osteoporosis disease models and bone preservation were employed. Five-month-old female New Zealand white rabbits were divided randomly into two groups: a control group (N = 10) and an osteoporosis group (N = 10). Bilateral ovariectomy was used in the osteoporosis group to develop the symptoms of osteoporosis [36], and sham surgery was used in the control group to avoid other factors affecting the experimental results. Twenty weeks after the surgery, the rabbits were euthanized, and the distal end of the left femur was dissected and used for PA assessment. However, 2 of the rabbits died during the modeling period. Therefore, 18 specimens (N = 9 for the osteoporosis group, N = 9 for the control group) were used to study the initial feasibility. This study was approved by the ethics committee of Nanjing University of Science and Technology (No. 202100129).

The bone specimens of the right complete leg for each rabbit from the osteoporosis and control groups were scanned using both DEXA (InAlyzer) and micro-CT (SCANCO, vivaCT 80) imaging systems to examine the BMD and BMA conditions, respectively. The imaging and quantified results of micro-CT imaging are shown in Fig. 1(B) and (C). The bone volume/total volume (BV/TV) data were obtained from the micro-CT images and used for further analyses. The statistical analysis results shown in Fig. 1(C) verify that the BV/TV of the osteoporosis group decreased significantly, with an 8.2% reduction. Fig. 1(D) shows a photograph of BMD measurement by the DEXA system. The DEXA images shown in Fig. 1(E) and (F) verify that the average BMD of the bone specimen in the osteoporosis group diminished significantly, resulting in a 7.3% reduction compared with the specimens in the control group. These analysis results based on DEXA and micro-CT imaging confirmed the differences in the BMA and BMD between the two rabbit bone models used in this study.

Fig. 1

Gold-standard imaging results for the physical and chemical properties of bone tissue: micro-CT, DEXA, and MRI imaging for the bone microarchitecture (BMA), bone mineral density (BMD), and the fat fraction assessment in the bone, respectively. (A) Photograph of BMA measurement by micro-CT system. (B) Examples of micro-CT images from the osteoporosis and control groups. (C) Statistical analysis of the bone volume/total volume (BV/TV) for the rabbit bone specimens from the osteoporosis and control groups. (D) Photograph of BMD measurement by DEXA system. (E) DEXA images from the healthy and osteoporosis rabbits. (F) Statistical analysis of the BMD. (G) Photograph of lipid fraction measurement by MRI system. (H) MRI images from the healthy and osteoporosis rabbit bone specimens. (I) Statistical analysis of the fat fraction in the bone for the osteoporosis group compared with the control group. * indicates p < 0.05 from unpaired t-tests comparing the results from the two specimen groups.

To validate the findings of the differences in the fat ratio obtained via the PA measurements between the control and osteoporosis groups, the bone specimens from the right leg of each rabbit were scanned using a 9.4 T MRI system (Bruker Biospec 94/20 USR), as shown in Fig. 1(G). For the acquisition of lipid and water separation images, the measurements were performed with the MRI system using the 3-point Dixon technique with water-lipid shifts of (0, π, 2π) [9]. The measured data were processed using the iterative decomposition of water and lipids with echo asymmetry and a least-squares estimation algorithm. The resulting water and lipid images were used to calculate the lipid fraction maps, as shown in Fig. 1(H). The basic parameters of the measurement were as follows: repetition time of 2500 ms; echo time of 24 ms; echo spacing of 8.0 ms; slice thickness of 2.0 mm; resolution of 0.313 × 0.313 mm; matrix of 128 × 128; field of view of 40 × 40 mm; and bandwidth of 100 kHz. The quantified MRI results shown in Fig. 1(I) verify that the average lipid fraction of the bone specimen in the osteoporosis group increases by approximately 10.1% compared with that in the control group.

Experimental setup

The PA signals of each bone specimen from the two groups were measured in the range of 690–950 nm at 10-nm intervals. The experimental setup is illustrated in Fig. 2(A). The light beam generated by a Nd:YAG laser-pumped optical parametric oscillator (Vibrant B, Opotek) was divided into two parts by a beam splitter. In one part, 10% of the laser energy was projected onto a black rubber component and the generated PA signal. The PA signal was recorded by an ultrasound transducer (V310-SU, Olympus) and used as the reference signal for the calibration of the laser energy. The remaining 90% of the laser energy was illuminated from the top to the bone surface directly. The diameter of the beam was approximately 8 mm, and the energy of the light fluency was controlled within 20 mJ/cm2, which is within the ANSI safety margin. After the PA signal generation, the signal was received by a needle hydrophone (HNC-1500, Onda Co., Sunnyvale, CA, USA) perpendicular to the laser beam with a broad bandwidth of 0–10 MHz, as shown in Fig. 2(A) inset. Then, the PA signal was enlarged by a pre-amplifier connected with the hydrophone to enhance the signal-to-noise ratio (SNR). In addition, the PA signals was amplified by 25-dB gain using an amplifier (5072PR, Olympus Corp., Tokyo, Japan) with a 1-MHz high-pass filter. Next, the PA signal was digitized and recorded by a digital oscilloscope (HDO6000, oscilloscope, Teledyne ap, USA). To improve the SNR, the PA signal was averaged over 50 laser pulses. Fig. 2(B) shows a photograph of a bone specimen and Fig. 2(C) shows the corresponding PA signal generated.

Fig. 2

Experimental setup. (A) Schematic of the experimental setup. (B) Photograph of the rabbit bone specimen. (C) An example of the PA signal generated by the bone specimen. (D) PA physicochemical analysis (PAPCA) process.

Signal processing

The PA signal generated by the bone specimens was calibrated using the peak-to-peak value of the PA signal generated by the black rubber component. Based on the calibrated PA signals for each tissue specimen, the PCS that integrates all the power spectra of the PA signals at the optical spectrum (690–950 nm) can be formulated and analyzed, as shown in Fig. 2(D).

The chemical information was then analyzed. First, the amplitude of the PA signal in the frequency domain at 1–8 MHz was summarized for each wavelength. The cutoff frequency at 1 MHz and 8 MHz was used to avoid the effect from large low-frequency (< 1 MHz) and high-frequency (> 8 MHz) noises. Next, a PA energy curve was obtained over the wavelength range of 690–950 nm, which represents the optical absorption spectrum of each bone specimen. As was measured by the PA technique, it is also referred to as the PA absorption curve in this paper. The PA spectrum curves for each bone sample obtained in three different directions were then averaged for further analysis to reduce measurement errors. Finally, spectral unmixing was conducted based on the least-squares regression method and the PA absorption spectrum of each bone specimen was decomposed. Then, the decomposed contribution of each chemical component to the PA absorption spectrum was used as the proportion of the corresponding chemical component.

To analyze the physical information, the ultrasonic frequency spectrum for each bone specimen at each wavelength was first quantified by linear fitting. Then, the slope of the fitted line was obtained and used as the quantified parameter. A PA slope curve at the range of wavelength 690–950 nm was obtained for each bone specimens. Next, the PA slope curves for each bone sample obtained in three different directions were averaged for further analysis to reduce measurement errors. Finally, the PA slopes from the two groups at several typical targeted wavelengths for hydroxyapatite (690 nm), blood (800 nm), and lipids (930 nm) were compared.

Experimental results

PCS results

Fig. 3(A) and (B) show the PCSs of the osteoporosis and normal control bone specimens at all optical wavelengths from 690 nm to 950 nm, respectively, as obtained through ex vivo measurements in the rabbit models. The X axis shows the optical wavelengths, whereas the Y axis shows the ultrasound frequencies. Considering the wavelength axis (X axis), the different amplitudes of the PA signal correspond to different optical absorptions of the major chemical components of the bone tissue samples. Over the spectral range of 690–950 nm, the main absorbed components in the bone are minerals (mostly hydroxyapatite), lipids, oxygenated hemoglobin, and deoxygenated hemoglobin. For the frequency axis (Y axis), the PA frequency domain power distribution at different wavelengths correlates to the microstructures of different tissue clusters in the bone specimens. For example, smaller thicknesses of the trabeculae or smaller sizes of the lipid cluster produce higher ultrasonic frequency components in the PCS.

Fig. 3

Examples of photoacoustic physicochemical spectrogram (PCS). (A) Averaged PCS spectrogram for the osteoporosis group. (B) Averaged PCS spectrogram for the control group.

Owing to the differences in chemical components and microstructures, each bone specimen with different conditions leads to a unique PCS signature as a combination of organic and non-organic tissues. A comparison of the PCS of the osteoporosis bone specimens with that of normal control shows the following differences: (1) higher magnitudes of the hemoglobin (Hb) at approximately 690 and 760 nm are observed in the former group (osteoporosis group); (2) higher magnitudes of the lipid fingerprints at approximately 930 nm are observed in the former group; (3) a further extension along the frequency axis toward a higher ultrasonic frequency at approximately 700 nm in the former group is observed.

The PCS is also sensitive to the ultrasonic attenuation in the bone. For the ultrasonic attenuation in bone tissue, it contains the ultrasonic attenuation in porous trabecular bone and damping effect of bone marrow. Comparing with the ultrasonic attenuation in porous trabecular bone, the damping effect of lipid or marrow is a minor factor to the ultrasonic attenuation in bone tissue. The reported research shows that major sources of attenuation in the bone include absorption, scattering, and phase cancellation [3]. Among those factors, scattering is a significant contribution to ultrasonic attenuation in bone tissue. Especially for the longitudinal-shear scattering, which due to the mode conversion at scatterer interfaces of the solid tissue (mineralized trabeculae) and soft tissue (marrow), it is a major contributor to the ultrasonic attenuation in bone. It means that higher proportion of scatterers (porous trabecular bone) corresponds to higher ultrasonic attenuation in the bone. Therefore, compared with normal bone, the lower proportion of trabecular bone results in lower ultrasonic attenuation of the PA signal in osteoporosis bone, leading to a higher value of PA frequency power distribution in the PCS of the osteoporosis specimen.

Chemical information analysis of the bone specimens from the PCS data

The multi-wavelength PA optical absorption curve can be obtained by summarizing the intensities of the PA frequency powers from 1 to 8 MHz for each PCS, and it can be used for further analysis. Fig. 4 shows the chemical analysis results of the osteoporosis (N = 9) and control (N = 9) groups. As shown in Fig. 4(A), the two solid lines were the respective averaged PA optical absorption spectrum measured from the osteoporosis and control groups normalized at 800 nm, and the shaded area along the solid lines was the standard deviation for each group. It is evident that the absorption curve exhibits peaks at approximately 690, 760, and 930 nm. Fig. 4(B) shows the optical absorption spectra of the major chemical components in the bone in the wavelength range of 690–950 nm. By comparing the PA absorption spectra shown in Fig. 4(A) and the optical absorption spectra shown in Fig. 4 (B), it can be seen that the strong relative absorption value around 690 nm shown in Fig. 4(A) is mainly caused by the absorption of both deoxygenated hemoglobin (Hb) and hydroxyapatite. The peak at approximately 760 nm in the PA spectrum was attributed to the absorption of deoxygenated Hb and lipids. The absorption peak at 930 nm was mainly a result of the light absorption of the lipids. Therefore, it is reasonable and preparatory to estimate that the Hb/blood and lipid/blood ratios in the osteoporosis group are higher than those in the control group. However, owing to the mixed optical absorptions at 690 and 760 nm, it is difficult to obtain the relative contents of hydroxyapatite and blood from the PA absorption curves directly.

Fig. 4

Chemical information extracted from the PCS of the bone specimens from the osteoporosis (N = 9) and control (N = 9) groups. (A) Relative optical absorption spectra of the bone specimens for the osteoporosis and control groups calculated from the PCSs. The solid curves show the average values of the bone specimens in each group, whereas the shaded area shows the standard deviation at each wavelength. (B) Optical absorption spectra of the major chemical components in the bone specimens, including the deoxy-hemoglobin (Hb), oxy-hemoglobin (HbO2), lipids, and hydroxyapatite. (C) Examples of fitted optical absorption spectra using the least-squares regression method for spectral unmixing. The quantified correlations, R2, were 0.94 and 0.97, respectively. The mean and standard deviation of R2 is 0.7 ± 0.1709 for the osteoporosis group, and 0.8 ± 0.1774 for the control group. (D) Results of the relative contents of minerals, lipids, and blood (Hb + HbO2) obtained from spectral unmixing in the two groups. * indicates p < 0.05, ns indicates p > 0.05.

By performing PA spectral unmixing with the least-squares regression method, the relative contribution of each major chemical component in the bone to the PA absorption spectrum was derived. Fig. 4(C) shows examples of PA and fitted spectra from the osteoporosis and control groups. After spectral unmixing, the relative contents of the major chemical components in the bone were obtained: hydroxyapatite, lipids, deoxygenated hemoglobin, oxygenated hemoglobin, and blood. The relative contents of the minerals, lipids, and blood—which represents the total deoxygenated and oxygenated hemoglobin—from the two groups of bone samples are shown in Fig. 4(D). To assess whether there is statistically significant differences in the chemical composition between the osteoporosis and control groups, an unpaired two-tailed independent samples t-test (with Welch’s correction in cases of unequal variances) was conducted using GraphPad Prism 7.0 software. Compared to the control group, the changes of chemical composition in the osteoporosis group showed statistical significance, including decreased mineral (hydroxyapatite) content and increased lipid content. This is consistent with the gold-standard DEXA and MRI results, as shown in Fig. 1. However, there was no statistically significant difference in blood changes between the control and osteoporosis groups. This may be because a portion of the blood was lost during the ex vivo bone specimen experiments.

Physical information analysis of the bone specimens from the PCS data

In addition to chemical analysis, the PCS data contain rich microstructural information for different chemical components, such as lipid clusters (rich in lipids), marrow (rich in both blood and lipid), and trabeculae (rich in hydroxyapatite), as shown in Fig. 5(A). Based on the quantitative analysis of the PCS for characterizing the tissue microstructure, the PA power spectral densities (PSDs) derived from the PCSs between the osteoporosis and control bones were compared and quantified.

Fig. 5

PA power spectral density (PSD) analysis for the osteoporosis and control groups. (A) The target PSDs marked in the PCS spectrogram. (B) Examples of the PA PSDs obtained at the wavelength of 800 nm and fitted lines for the osteoporosis and normal bones. (C) Quantified PA parameter (slope) over the wavelength range of 690–950 nm for the osteoporosis and control groups. (D) Statistical analysis results for the wavelengths of 690, 800, and 930 nm. * indicates p < 0.05, * * indicates p < 0.01, ns indicates p > 0.05.

As an example, Fig. 5(B) shows the PSD curves of the bone specimens from the osteoporosis and control groups. It was found that the PSD of the osteoporosis bone specimens contains more relatively low-frequency components than that in the control group. To characterize the difference in the PSD curves from the osteoporosis group and control group, the PSD curves were quantified using linear fitting in the frequency range of 1–8 MHz, as described in Section 2.3, and then the quantified parameter slope of the linear fitting line was obtained. Next, the parameter slope of the PSD was quantified at all wavelengths in the range of 690–950 nm, and then the slope spectra for each bone specimen were obtained. The slope spectra that present the slope as a function of the wavelength are shown in Fig. 5(C). The two solid lines in this figure show the averaged PA slope spectra measured from the two groups of bone specimens, and the standard deviation is shown by the shaded area along each curve. It is obviously that the average PA slope for the osteoporosis group is lower than that for the control group at most of all the wavelengths in the range of 690–950 nm.

The quantified PA parameter slope is affected by both the microstructure of absorbers and the ultrasonic attenuation in bone tissue. The major absorbers in the bone in the wavelength range of 690–950 nm are the minerals, blood, and lipids. Among these three types of absorbers, the contributions of the blood (~40%) and lipids (~40%) to the PA signal are much greater than that of the minerals (~20%) in the wavelength range of 690–950 nm, as shown in Fig. 4(D). Therefore, compared with microstructure of trabecular bone, the PA parameter slope is more dependent on the size of the marrow which is mostly combined by blood and lipid. Then, the lower slope of the osteoporosis group can be explained by the fact that the osteoporosis bone has a larger porosity (as shown in Fig. 1(C)) that is filled with bone marrow, which exhibits a high optical absorption at the wavelength targeted in the blood and lipid. Additionally, large-scale spatially distributed marrow generates higher PA signals in low frequencies components. Furthermore, compared with control bone, the lower proportion of trabecular bone results in lower ultrasonic attenuation in osteoporosis bone. Those factors lead to a high value of the relatively low ultrasonic frequency components and a lower slope in osteoporosis bone, as shown in Fig. 5(B) and (C), respectively.

To evaluate whether these two average PA slope curves have a statistically significant difference, a two-way ANOVA test was conducted using GraphPad Prism 7.0 software with p < 0.05. The results indicated that the PA slope spectra can be used to distinguish the bone directly. To further study the microstructure of different clusters in the bone, the slopes at the wavelengths of 690, 800, and 930 nm were compared, as shown in Fig. 5(D). An unpaired two-tailed independent samples t-test (with Welch’s correction in cases of unequal variances) was conducted for the slope results by using GraphPad Prism 7.0 software. For the wavelength of 690 nm, it targeted both hydroxyapatite and Hb. For the wavelength of 800 nm, it targeted the blood that is distributed uniformly in the marrow and can be used to characterize the size of the bone porosity that is filled by the bone marrow. For the wavelength of 930 nm, it targeted the lipid and can be used to characterize the size of the lipid clusters in the bone marrow.

The statistical analysis results for 800 nm and 930 nm indicated p < 0.05, demonstrating that the PA spectral parameters (slope) acquired at those two laser wavelengths (800 and 930 nm) can be used to differentiate the control group from the osteoporosis group by assessing the changes in the bone microarchitecture. This can be explained by the fact that the existence of large-scale spatially distributed marrow in osteoporosis bone results in high PA signals at low frequencies, and then, leads to a lower slope for the osteoporosis group at the wavelength of 800 nm. Furthermore, the lower slopes in the osteoporosis group at the wavelength of 930 nm indicate that the lipid clusters in the osteoporosis group is larger than that in the control group, which is consistent with the results of past studies [37], [38], [39], [40]. Besides, compared with control bone, the ultrasonic attenuation in osteoporosis bone is lower, which leads to much more energy and feature information of the PA frequency spectrum obtained at each wavelength are preserved in osteoporosis group.

It was found that there were no significant differences in the slopes obtained at the wavelength of 690 nm. At this wavelength, both the hydroxyapatite and Hb contribute to the quantified PA parameter slope. The hydroxyapatite is the major components of trabecular bone, while the Hb is one of the major components of the bone marrow. Therefore, the slope at the wavelength of 690 nm is not only affected by the microstructure of trabecular bone but also bone marrow filled in the pore of bone tissue. In addition, the PA parameter slope is also very sensitive to the ultrasonic attenuation in the bone. Those comprehensive factors lead to no significant differences results of quantified PA parameter slope at the wavelength of 690 nm.

Statistical study and categorization of bone conditions by SVM classification

In this study, each of the bone specimens in the control group (N = 9) and osteoporosis group (N = 9) was scanned, therefore, a total of 18 PCS were generated. For each PCS, several parameters used for characterization of chemical and physical properties were obtained, as shown in Fig. 4 (D) and Fig. 5 (D). For those quantified PCS parameters, the content and microstructure of mineral, blood, and lipid were used as features to classify the bone conditions. The content was represented by the relative content of mineral, blood, and lipid, while the microstructure of mineral, blood, and lipid was represented by the slope at the wavelength of 690 nm, 800 nm, 930 nm, respectively. SVM, which is a multi-variant analysis MATLAB tool, was used to examine the feasibility of characterization of the bone conditions based on those quantified PCS parameters.

To investigate the categorization of the bone conditions by employing the SVM, a three-fold cross-validation approach was used in this study to avoid the overfitting of the SVM model. The nine sets of data acquired from each bone condition were randomly divided into three groups with three datasets. Two of the three groups for each bone condition were used for training the SVM, i.e., a total of 12 datasets (2 groups × 3 datasets × 2 bone conditions), and the remainder was used for testing, i.e., a total of 6 datasets (1 groups × 3 datasets × 2 bone conditions). A three-fold cross-validation was performed using all the datasets for training and testing the SVM model. In this study, a linear kernel was used in SVM algorithm. Essentially, the linear kernel is used to create a linear hyperplane which can separates the datasets into two classes. Besides, the linear kernel SVM model is less prone to overfitting than non-linear kernel, even for the small amount of data. The liner kernel is defined as:where x is the variable matrix, i and j are the indices of the matrix. To better examine the feasibility of characterization the bone conditions based on physical parameters and chemical parameters, we conducted the categorization of bone conditions based on 3 kinds of different parameters combinations. First, only the contents (including the relative contents of mineral, blood and lipid) were used to perform the categorization. Second, only slope was used to perform the categorization, including the slope at the wavelength of 690 nm, 800 nm and 930 nm. Finally, both the content and slope were used to perform the categorization of bone conditions. For the categorization with each parameter combination, a three-fold cross-validation approach was used to avoid the overfitting of the SVM model. The testing accuracy was considered as the categorization accuracy. The SVM categorization results are shown in Fig. 6 and the SVM categorization based on 3-fold cross-validation are shown in Table 1.

Fig. 6

SVM categorization results. (a) SVM categorization based on the relative contents of minerals, lipids, and blood. (b) SVM categorization based on the slopes at 690, 800, and 930 nm. “▲” stands for training data, “•” stands for testing data, “⚪” stands for mis-classification data in testing datasets.

Table 1

SVM categorization based on 3-fold cross-validation. All the data sets are randomly divided into 3 groups, each group are in turns used for testing, while the other two groups are used for training. Cate. = Categorization. Samp. = Sample.

Test sample condition	Cate. #	Test samp.#	Cycle 1			Cycle 2			Cycle 3
			Relative content	Slope	Both	Relative content	Slope	Both	Relative content	Slope	Both
Osteoporosis	1	1	1	2	2	1	1	1	2	1	1
	1	2	1	1	1	2	2	2	2	2	2
	1	3	1	1	1	2	1	1	1	1	1
Control	2	4	2	2	2	1	1	1	2	2	2
	2	5	1	2	2	1	2	2	2	2	2
	2	6	2	2	2	2	2	2	2	1	2
Training accuracy			66.7%(8/12)	66.7%(8/12)	91.7%(11/12)	83.3%(10/12)	83.3%(10/12)	91.7%(11/12)	75.0%(9/12)	83.3%(10/12)	83.3%(10/12)
Categorization accuracy (Testing accuracy)			83.3%(5/6)	83.3%(5/6)	83.3%(5/6)	33.3%(2/6)	66.7%(4/6)	66.7%(4/6)	66.7%(4/6)	66.7%(4/6)	83.3%(5/6)

As shown in Table 2, the SVM trained based on only the PA parameter of relative contents and only the PA parameter of slopes achieved average accuracies of 61.1% and 72.2%, respectively, in terms of categorizing the test datasets. When both spectral parameters were used, the accuracy reached 77.8%. It indicates that the chemical information combined with the physical information can provide a better classification accuracy.

Table 2

The training and testing accuracy by 3-fold cross-validation SVM.

3-fold SVM (the ratio between training set and testing set was 2:1)
	Relative Content	Slope	Both
Mean training accuracy	75.0%	77.8%	88.9%
Mean testing accuracy	61.1%	72.2%	77.8%

To evaluate the sensitivity of the SVM model used in this study, the receiver operating characteristic (ROC) curves of SVM models based on different training data were analyzed, as shown in Fig. 7. It shows that the area under curve (AUC) values for all testing groups are larger than 0.70. Especially for the SVM categorization based on the PA parameters of both relative content and slope, they were larger than 0.85 for all cycles as shown in Fig. 7 (C) (F) (I). It means that the SVM model with the linear kernel we used based on the parameters of chemical information and physical information acquired from PAPCA can provide a better classification accuracy.

Fig. 7

ROC curves and associated AUC values of the SVM models based on different parameters and circles. (A) SVM categorization model of relative content in cycle 1. (B) SVM categorization model of slope in cycle 1. (C) SVM categorization model of both relative content and slope in cycle 1. (D) SVM categorization model of relative content in cycle 2. (E) SVM categorization model of slope in cycle 2. (F) SVM categorization model of both relative content and slope in cycle 2. (G) SVM categorization model of relative content in cycle 3. (H) SVM categorization model of slope in cycle 3. (I) SVM categorization model of both relative content and slope in cycle 3.

Discussion and conclusion

This study demonstrated that the PCS obtained from the PA measurements of bone tissue can be used to characterize the microstructure and chemical contents of the organic matrix and non-organic tissue in bone, which are associated with the pathology information of bone tissue. The results showed that the quantified PA parameters, namely the “slope” and “relative content,” from the PCS were highly correlated with the gold-standard DEXA, micro-CT, and MRI results, which can provide rich information on the BMD, BMA, and lipid content, respectively. Furthermore, based on the assessment of both the microstructure and chemical composition, the PA technique is better able to differentiate between the osteoporosis and control groups.

In this study, PAPCA is an easy means to characterize bone pathology. However, the PCSs used for PAPCA are affected by the detection mode. Therefore, the PA absorption spectrum and PA power density spectrum obtained from the PCS should be corrected before quantification depending on the detection mode used in the study. For example, during detection in the reflection mode, the PA signal contained a large reverberation from the tissue surface. In the transmission mode, the PA signal was less affected by the reverberation from the surface; however, light and ultrasound were attenuated in the bone. Therefore, the PCS was affected by “spectral coloring” and “frequency-dependent ultrasonic attenuation,” which should be corrected before analysis. “Spectral coloring” can be eliminated by using model-based or machine-learning-based methods in future work, while the frequency-dependent ultrasound attenuation can be partly reduced by compensation. It is worth noting that the BUA obtained from the QUS system can reflect the relationship between the frequency and ultrasonic attenuation; therefore, the “frequency-dependent ultrasonic attenuation” in PCSs can be corrected by using a combined QUS and PA system. Thus, it is reasonable to expect that the PA signals from bone could carry richer and more direct information regarding the microstructure and metabolism of bone tissue.

However, this study has several limitations. First, the bone tissue was scanned ex vivo after it was sacrificed, which led to the inability to observe oxygenated and deoxygenated hemoglobin individually. Therefore, in vivo studies should be conducted in the future. Second, because of the relatively low optical absorption of hydroxyapatite compared with that of blood, the contributions of bone minerals to the PA signal were relatively lower than those of blood. This led to reduced mineral contrast in the osteoporosis and control bones. Therefore, the algorithm for decomposition should be improved to achieve higher precision and improved contrast in future studies. Third, for the in vivo case, the individual soft tissue overlying the bone tissue could result in uncertain optical and ultrasonic attenuations of the PA signal. However, because the bone specimens used in this study were relatively small and thin, the effects of light and ultrasonic attenuation in the soft and bone tissues were not considered. In a previous work, we indicated that the PA signal from the overlying soft tissue in vivo can be partly removed by using ultrasound-guided detection and segmentation [41]. By using the ultrasound-guided method, the thickness of soft tissue can be evaluated, and the PA signal acquired by using transmission mode and its reverberations can be partly eliminated. In addition, we investigated the light attenuation in both bone and overlying soft tissue for in vivo bone models [28]. Therefore, light attenuation in the overlying soft tissue can be partly compensated by using a model-based method. The effects of the overlying soft tissue on the PCS should be further investigated for in vivo studies in future work. Forth, the samples were limited for categorization of bone conditions by SVM classification in this study. In order to solve this issue, we used the SVM algorithm with linear kernel which is suitable for classification of small amount of data. Essentially, the linear kernel is used to create a linear hyperplane which can separates the datasets into two classes. Besides, the linear kernel SVM model is less prone to overfitting than non-linear kernel, even for the small amount of data. Furthermore, to prevent the overfitting based on the small number of training datasets, the 3-fold cross-validation was used. Therefore, the SVM algorithm based on linear kernel and 3-fold cross-validation is applicable for small amount of data. In future work, large number of training data sets should be included to further investigate the feasibility of the PCS combined with SVM algorithm in the classification of bone conditions.

In clinical applications, DEXA imaging, micro-CT imaging, MRI, and QUS are already being used for bone assessment. The gold standard is DEXA imaging, which is the most commonly used clinical method; it focuses on the BMD information only and provides less information on bone metabolism. Further, micro-CT imaging can determine the BMD and BMA with a higher resolution than DEXA, but it does not provide metabolism information. In addition, MRI can provide only limited information on bone marrow, therefore, it is not commonly used for clinic diagnosis of osteoporosis. Moreover, the QUS bone assessment techniques can provide BUA and SOS information, which is highly associated with the BMD and BMA; however, QUS provides limited metabolism information. In the PA technique for bone assessment, based on the endogenous contrasts of the chemical components in the bone tissue, the PCS can simultaneously present the microstructure and metabolic information of bone without multiple detection methods. Compared to the conventional DEXA and MRI methods, PAPCA technique for bone assessment, based on the endogenous contrasts of the chemical components in the bone tissue, it can simultaneously present the structure and metabolic information of bone without multiple detection methods. The PA method is advantageous as it is noninvasive and nonionizing, and enables sufficient penetration into calcified and non-calcified tissue. Hence, the PA technique shows considerable potential for clinical applications.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.