Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Real-time and accurate estimation ex vivo of four basic optical properties from thin tissue based on a cascade forward neural network

Open Access Open Access

Abstract

Double integrating sphere measurements obtained from thin ex vivo tissues provides more spectral information and hence allows full estimation of all basic optical properties (OPs) theoretically. However, the ill-conditioned nature of the OP determination increases excessively with the reduction in tissue thickness. Therefore, it is crucial to develop a model for thin ex vivo tissues that is robust to noise. Herein, we present a deep learning solution to precisely extract four basic OPs in real-time from thin ex vivo tissues, leveraging a dedicated cascade forward neural network (CFNN) for each OP with an additional introduced input of the refractive index of the cuvette holder. The results show that the CFNN-based model enables accurate and fast evaluation of OPs, as well as robustness to noise. Our proposed method overcomes the highly ill-conditioned restriction of OP evaluation and can distinguish the effects of slight changes in measurable quantities without any a priori knowledge.

© 2023 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

1. Introduction

The optical properties (OPs) of biological tissues are characterized by wavelength-dependent parameters such as the absorption coefficient ${\mu _a}$, scattering coefficient ${\mu _s}$, anisotropy factor g, and refractive index n [1]. Accurate acquisition of OPs is a prerequisite for describing the transmission and distribution of light in tissues, and is also the basis of diffuse optical tomography [2] and laser surgery [3]. The double integrating sphere (DIS) system is an ex vivo measurement technique that is widely used for measuring the OPs of tissues. It provides a fast and simple method of characterizing the spectral information of biological tissues such as the brain and skin [4,5]. In a typical DIS system, visible to near-infrared light is used to illuminate the tissue, and multiple reflected light inside the integrated sphere produce steady illumination. By detecting the power on the inner surface of the sphere, the diffuse reflectance ${R_d}$, diffuse transmittance ${T_d}$, and collimated transmittance ${T_c}$ of the tissue placed between the two spheres can be obtained, which in turn derives the OPs [5]. In contrast to other in vivo measurement techniques, such as spatially resolved reflectance and spatial frequency domain imaging, the DIS system enables an external spectrophotometer connecting to simultaneously estimate the OPs of tissues within a wide spectral range. In addition, the system is sensitive to the layered structure of biological tissue and thus has become the most used and accepted measurement method [6].

However, the accurate determination of OPs depends not only on the DIS system but also on the efficiency of the reconstruction models that establish the intrinsic inverse mapping relationship between a set of spectral measurements and their corresponding OPs. Several methods have been proposed to evaluate OPs from DIS systems, such as the inverse adding-doubling (IAD) and inverse Monte Carlo (IMC) methods. Both methods model the forward process of light propagation in biological tissues by the adding-doubling (AD) and Monte Carlo (MC) methods, respectively. Combined with the spectral information measured by the DIS system as the input, forward models are solved iteratively using specific optimization algorithms until the errors between the calculated and measured values converge to the target accuracy [79]. The IAD method is fast but only applicable to tissues with a known refractive index [7,11]. The IMC method was developed as the gold standard due to its high accuracy. However, because the MC method relies on repeated random sampling, its application in real-time computation and clinical applications is inevitably limited. Although many accelerated methods have been proposed, such as scaling MC, perturbation MC, hybrid MC, variance reduction, and parallel-computed MC [10], the inverse determination of OPs remains time-consuming. Furthermore, the numbers of OPs extracted by the IMC and IAD methods are limited by the measurable quantity of tissues. When the sample is so thick that the ${T_c}$ and ${T_d}$ cannot be detected, the ${\mu _s}$ and g are simplified to a single parameter by introducing the reduced scattering coefficient $\mu _s^{\prime}$ or simply assuming g is a constant [7,12]. This simplifies the ill-posedness of the inverse problem so that solutions can be determined. However, for strongly scattering tissues, both the direction and distance that photons travel between every two scattering events have a considerable effect on the distribution of light in tissues; therefore, the separate determination of ${\mu _s}$ and g is essential. Meanwhile, in a continuous wave measurement system, the refractive indexes of the tissue sample ${n_s}$ and cuvette holder ${n_c}$ are affected significantly by slight changes in wavelength. Therefore, treating the refractive index as a fixed value limits the accuracy of estimated OPs and the consequent diagnosis of disease [13,14].

Deep learning, which benefits from its powerful multi-layer nonlinear learning capability, can automatically fit and approximate complex relationships in high-dimensional space [15], thus eliminating the need for computationally demanding inverse algorithms. Hence, deep learning has been extensively used in various frameworks for extracting OPs [1517]. The first use of neural networks to infer OPs from DIS measurements was demonstrated by Li et al. [18], who constructed a back-propagation neural network that was able to retrieve ${\mu _a}$ and ${\mu _s}^{\prime}$ from ${R_d}$ and ${T_d}$ precisely. Nishimura et al. [19] developed a multilayer perceptron (MLP) model that was designed to estimate ${\mu _a}$ and ${\mu _s}$ at specific wavelengths using ${R_d}$ and ${T_d}$ combined with the pre-fitted wavelength-dependent g and ${n_s}$ of human dermal tissue. Although both Li et al. [18] and Nishimura et al. [19] achieved an inverse evaluation of OPs through neural networks, only part of the OPs were estimated because slight changes in ${R_d}$ and ${T_d}$ led to significant changes in the estimated OPs when the physical thickness of the tissue was relatively thin [14]. Accordingly, in the previous studies the thicknesses of tissue samples were usually thick that resulted in the ${T_c}$ becoming too weak to be detected by the DIS system. As a result, only ${R_d}$ and ${T_d}$ were used as inputs that limited the number of valid features that neural networks could extract. Additionally, the accuracy saturates and then rapidly deteriorates as the network depth increases [20]. The creation of direct short paths from early layers to later layers provides a solution to overcome the degradation of deep neural networks and has been broadly used in various neural networks, such as ResNet and DenseNet [2024]. The cascade forward neural network (CFNN), as a variant of an MLP model, is unique in that each layer (except the input layer) is connected to all the previous layers through a weighted link. Profiting from such a network framework, the CFNN can successfully avoid the degradation of the deep neural network and enhance the network's ability to explore the relationship between the input and output layers [20,22].

The purpose of this study is to propose a CFNN-based inverse model for the accurate extraction of the ${\mu _a}$, ${\mu _s}$, g and ${n_s}$ of biological tissues at any wavelength from fitted ${n_c}$ and ${R_d}$, ${T_d}$, ${T_c}$ of thin ex vivo tissues measured by the DIS system. GPU-based MC modeling software of light transport in multi-layered tissues, CUDAMCML, was adopted to quickly and efficiently build the database required for training. Tens of thousands of MC simulation datasets were used to train the CFNN. The OPs were obtained upon completion of training by combining the refractive index of the cuvette holder calculated using the Sellmeier equation, thus testing the robustness and generalization of the CFNN. To the best of our knowledge, the simultaneous determination of OPs such as ${\mu _a}$, ${\mu _s}$, g and the successful evaluation of the refractive index of the tissue sample ${n_s}$ from the combination of the DIS system and neural networks, has not yet been reported. Our new method overcomes the highly ill-conditioned restriction of the evaluation of OPs from thin ex vivo tissues and can distinguish the effects of slight changes in measurable quantities on each OP without any a priori knowledge and complex data processing. Furthermore, the relationship between the wavelength and OPs is fully considered, enabling real-time and highly accurate estimation of all basic OPs within a wide spectral range.

2. Method

2.1 Cascade forward neural network

Similar in construction to MLP, CFNN have direct connections from the input to the hidden and output layers. These additional connections enable the CFNN to perform nonlinear mapping from input to output while maintaining the linear relationship between the input and output layers [23]. The relation through the inputs to the first hidden layer of the CFNN can be described as [22]:

$$\varphi _m^{(1)} = b_m^{(1)} + \sum\nolimits_{i = 1}^j {w_{mi}^{(1)}{x_i}}$$
$$O_m^{(1)} = f_m^{(1)}(\varphi _m^{(1)})$$
where ${x_i}$ is the input signal of the ${i^{th}}$ dimension, j is the number of dimensions, $w_{mi}^{(h)}$ denotes the weight of node m for the input i of layer h, $b_m^{(h)}$ is the weight from the bias node m of layer h, $f_m^{(h)}$ is the activation function of the ${m^{th}}$ node of layer h, and $O_m^{(h)}$ is the output node m for layer h. The superscript h represents the ${h^{th}}$ hidden layer and n is the number of neurons.

Each subsequent layer of the CFNN is connected to the input layer and all previous hidden layers via weighted links. The connections between these layers are determined by the following equations:

$$\varphi _m^{(h)} = b_m^{(h)} + \sum\nolimits_{p = 1}^{h - 1} {\sum\nolimits_{k = 1}^n {O_{mk}^{(p)}w_{mk}^{(p)} + \sum\nolimits_{i = 1}^j {w_{mi}^{(h)}{x_i}} } }$$
$$O_m^{(h)} = f_m^{(h)}(\varphi _m^{(h)})$$

In this study, the Neural Network Toolbox (The MathWorks, Inc.) in MATLAB was used to establish a dedicated CFNN inverse model for each OP. Thus, all networks consisted of only one node in the output layer. This strategy enabled us to eliminate the impact of potential crosstalk that might occur when two or more OPs are simultaneously estimated by the same network [15]. Moreover, compared with broadening the width of the neural network, extending the depth of the neural network makes it easier for it to have better nonlinear approximation capabilities [25]. Therefore, each network comprises one input layer, one output layer, and eight hidden layers, each of which has 10 neurons. The inputs to the CFNN include the diffuse reflectance ${R_d}$, diffuse transmittance ${T_d}$, collimated transmittance ${T_c}$ measured by the DIS system, and refractive index of the cuvette holder ${n_c}$ fitted by the Sellmeier equation. The structure of the CFNN model and the workflow of OP mapping using the model are shown in Fig. 1. To avoid possible risks from different dimensions of the training data, inputs and outputs were scaled to [-1, 1] by min-max normalization before training. A linear, rectified linear unit (RELU) activation function was used in the output layer and in each hidden layer. The update of the weights and biases of the CFNN was based on the Levenberg-Marquardt optimization with Bayesian regularization using a mean square error (MSE) loss function and a learning rate of ${10^{ - 4}}$. Moreover, hyperparameters were tuned via grid search on a log scale and the training was stopped early when the value of loss function no longer improved for 10 successive epochs.

 figure: Fig. 1.

Fig. 1. Structure of CFNN model and workflow of OP mapping.

Download Full Size | PDF

2.2 Wavelength-dependent refractive index

To stabilize the sample and minimize the influence of the rough surface on the measurement results in the DIS system, the tissue sample is usually held in a quartz cuvette [7], since the Fresnel equation is adopted in both the MC and AD methods to cope with the mismatched boundaries [26]. Consequently, the spectral information measured by the DIS system contains contributions from the refractive index of the cuvette holder ${n_c}$ and tissue sample ${n_s}$. Owing to the heterogeneity of biological tissues, the refractive index of most tissues cannot be directly used. In contrast, the refractive index of the cuvette as a holder is already available when the wavelength is known, providing potential information for the inverse determination of the OPs of tissues. The relationship between the refractive index and wavelength of the quartz cuvette can be derived from the Sellmeier equation [27]:

$${n^2}(\lambda ) - \textrm{1 = }\frac{{{B_1} \cdot {\lambda ^2}}}{{{\lambda ^2} - {C_1}}} + \frac{{{B_2} \cdot {\lambda ^2}}}{{{\lambda ^2} - {C_2}}} + \frac{{{B_3} \cdot {\lambda ^2}}}{{{\lambda ^2} - {C_3}}}$$
where n is the refractive index, $\lambda $ is the wavelength in $\mathrm{\mu m}$, and ${B_1}$, ${B_2}$, ${B_3}$, ${C_1}$, ${C_2}$ and ${C_3}$ are the Sellmeier coefficients. Once the refractive index of the cuvette holder at a given wavelength is obtained, combined with the spectral information measured by the DIS system, the OPs and refractive index of the tissue at the corresponding wavelength can be immediately obtained through the CFNN model.

2.3 Collimated transmittance

In a DIS system, the collimated transmittance is often too weak to detect because of the high total attenuation coefficient ${\mu _t}$ or the high physical thickness of the tissue. Nonetheless, accurate measurement of the collimated transmittance can be achieved by reducing the thickness of the tissue sample, using a sensitive detector combined with a stronger source, or collecting the collimated transmitted light using a third integrating sphere [7]. Previous studies neglected the unique advantage of ex vivo measurements over in vivo measurements (i.e., the ability to separately analyze the layered structure of tissues); therefore, the samples used in the experiments were frequently thick. Consequently, the collimated transmittance measured by the DIS system was not used as an effective input to the neural network. However, the collimated transmittance contains significant information about the OPs, that is [26]:

$${T_c} = \frac{{(1 - {r_1})(1 - {r_2})\exp ( - \tau )}}{{1 - {r_1}{r_2}\exp ( - 2\tau )}}$$
where $\tau $ is the optical thickness, ${r_1}$ and ${r_2}$ represent the unpolarized Fresnel reflections for light that are normally incident on the front and back surfaces of the tissue, respectively. In this study, by reducing the physical thickness of tissue samples, the collimated transmittance was maintained within the detectable range; thus, it may be used as an effective input of the neural network and participate in the estimation of OPs.

2.4 Integrating sphere correction

We consider a setup with a double integrating sphere (DIS) system, which is used to measure the diffuse transmittance ${T_d}^m$, diffuse reflectance ${R_d}^m$, and collimated transmittance ${T_c}^m$. The following definitions are used:

$${T_d}^m = \frac{{{P_t}^s}}{{{P_t}^t}},\textrm{ }{R_d}^m = \frac{{{P_r}^s}}{{{P_r}^r}}$$
where ${P_t}^s$ and ${P_t}^t$ are the detected signals in the transmission sphere for the sample and the standard, respectively. ${P_r}^s$ and ${P_r}^r$ are the detected signals in the reflection sphere for the sample and the standard, respectively. The superscripts t, r and s represent the transmission standard, the reflection standard, and the sample, respectively.

In the experiment setup, it is reasonable to assume that the beam diameter is small (far smaller than the sample diameter) with the illumination incident on the sample center perpendicularly, since the beam diameter had only a small effect on the OPs estimation [13]. Besides, a baffle is usually installed between the sample port and the detector port to avoid the collection of direct light from the sample. The equations for the detected signals obtained from two identical integrating spheres with baffles when the target is illuminated by a collimated incident light can be derived according to [7,28]:

$${P_r}^i = \delta \frac{{mR_c^i(1 - R_{dd}^is) + R_d^im\alpha + \rho T_{dd}^i(T_d^i + mT_c^i)m\alpha }}{{[1 - (m\alpha + R_{dd}^is)][1 - {{(\rho T_{dd}^i)}^2}]}}{P_{in}},\textrm{ }i = r,s$$
$${P_t}^j = \delta \frac{{mT_c^j(1 - R_{dd}^js) + T_d^jm\alpha + \rho T_{dd}^j(R_d^j + mR_c^j)m\alpha }}{{[1 - (m\alpha + R_{dd}^js)][1 - {{(\rho T_{dd}^i)}^2}]}}{P_{in}},\textrm{ }j = t,s$$
where $R_{dd}^i$ and $R_d^i$ are the diffuse reflection factor of the object i with the diffuse and collimated incident light, respectively; $T_{dd}^j$ and $T_d^j$ are the diffuse transmission factor of the object j with the diffuse and collimated incident light, respectively; $R_c^i$ ($T_c^i$) is the specular reflection (collimated transmission) factor of the object i with collimated incident light; ${P_{in}}$ is the incident power; $\rho $ is the ratio:
$$\rho = \frac{s}{{1 - (m\alpha + {R_{dd}}s)}}$$

Other parameters and their corresponding values are summarized in Table 1. Parameters used for DIS correction were carefully selected based on the published literature [5,7]. To cast Eq. (8) and (9) in the closed form we used the approximation proposed in [29], namely:

$$R_{dd}^s = R_d^s(1 - \Delta r_c^s) + \Delta r_c^s$$
$$T_{dd}^s = T_d^s(1 - \Delta r_c^s)$$
where $\Delta r_c^s = r_{cd}^s - r_{cc}^s$, $r_{cd}^s$ and $r_{cc}^s$ are Fresnel reflection coefficients of the sample for diffuse light and collimated light at right angle incidence, respectively.

Tables Icon

Table 1. Parameters of the integrating sphere

For the reference power, a diffuse reflectance standard plate with 99.9% diffuse reflectance is used as a reflectance standard and air is used as a transmission standard.

2.5 Database preparation

The basic framework of MC modeling of light transport in multi-layered tissues (MCML) was initially developed by Wang et al. [30]. Because MC simulation depends on repeated random sampling to model the radiative transfer process, it inevitably requires substantial computation time to produce reliable results. To enhance the computation speed, CUDAMCML [31], which is based on Compute Unified Device Architecture (CUDA) to realize parallel simulation of multiple photonic packets on a graphics processing unit (GPU), was used to calculate the diffuse reflectance ${R_d}$, diffuse transmittance ${T_d}$, and collimated transmittance ${T_c}$ of the tissue for different combinations of OPs, thereby establishing the training database of the CFNN. The ranges of absorption coefficient ${\mu _a}$, scattering coefficient ${\mu _s}$, anisotropy factor g, and refractive index ${n_s}$ were 0.1–2 cm-1, 10–400 cm-1, 0.65–0.95, and 1.35–1.45, respectively. The ranges used for simulations were carefully selected based on the published literature [4,32,33], and the increments for ${\mu _a}$, ${\mu _s}$, g and ${n_s}$ were set as 0.2 cm-1, 20 cm-1, 0.02 and 0.01, respectively. In each CUDAMCML simulation, 1 × 106 photon packets were launched [34], and the Henry–Greenstein phase function was used. We considered a 5-layer geometry (air/glass/tissue/glass/air) with uniform scattering and absorbing properties in each layer. For the non-absorbing glass layer, the absorption coefficient ${\mu _a}$ and scattering coefficient ${\mu _s}$ were set to 0, and the anisotropy factor g was set to 1 to ensure that light propagates in the glass along a straight way. To ensure that the collimated transmittance of the tissue could be accurately detected by the DIS system, the physical thickness of the sample was set to 100 µm and the radial width was infinite, thus the sample could be considered as an infinite plane-parallel slab of finite thickness. In addition, because the tissue was sandwiched by a quartz cuvette, the influence of the mismatched boundary caused by the cuvette holder (1.48$\le {n_c} \le $1.55, with an increment of 0.01) was also considered in the simulations.

2.6 Performance metric

The relative error (RE) and mean relative error (MRE) reflect the degree to which the measured value deviates from the true value, and the R-square (R2) is used to evaluate the fitting level of the regression line to the estimated value. Accordingly, RE, MRE, and R2 were used to evaluate the performance of the CFNN model in OP mapping. They were calculated using the following:

$$RE = \frac{{|{{X_{estimated}} - {X_{expected}}} |}}{{{X_{expected}}}} \times 100\%$$
$$MRE = \frac{1}{n}\sum\nolimits_{i = 1}^n {\frac{{|{{X_{i,estimated}} - {X_{i,expected}}} |}}{{{X_{i,expected}}}}} \times 100\%$$
$${X_{avg}} = \frac{1}{n}\sum\nolimits_{i = 1}^n {{X_{i,expected}}}$$
$${R^2} = 1 - \frac{{\sum\nolimits_{i = 1}^n {{{({X_{i,estimated}} - {X_{i,expected}})}^2}} }}{{{{\sum\nolimits_{i = 1}^n {({X_{i,estimated}} - {X_{avg}})} }^2}}}$$
where X can represent ${\mu _a}$, ${\mu _s}$, g or ${n_s}$, ${X_{estimated}}$ and ${X_{expected}}$ are the estimated and expected values of OPs, respectively. Smaller values of RE, MRE, and the R2 value close to one mean better performance of the CFNN model for determining optical property mapping.

3. Results and discussion

3.1 Impact of thickness

In the DIS system, there is no defined criterion for the selection of the physical thickness of the tissue sample. However, studies have shown that reducing the thickness of the tissue can efficiently eliminate the potential risk of light losses, and thereby significantly improve the accuracy of the estimated ${\mu _a}$ [14,35]. To evaluate the effect of the physical thickness of the tissue sample on the inverse determination of the OPs, CUDAMCML was used to calculate the diffuse reflectance ${R_d}$ and diffuse transmittance ${T_d}$ using 31,200 sets of OP combinations of ${\mu _a}$, ${\mu _s}$, and g in the ranges 0.1–3 cm-1, 10–400 cm-1, and 0.7–0.95, respectively, in increments of 0.1 cm-1, 10 cm-1, and 0.01, respectively, for four different thicknesses. A total of 2500 sets of OP combinations were randomly selected for the convenience of display. Figure 2(a) shows the relationship between ${R_d}$ and ${\mu _a}$, ${\mu _s}$, and g for different physical thicknesses ($d$ = 50, 100, 300, 500 µm). Each dot represents a value for ${R_d}$ through its color. The differences between the ${R_d}$ of tissues in different OP combinations gradually decreases as the physical thickness of the tissue reduces, resulting in the distribution of colored dots becoming more concentrated. Meanwhile, the distribution of ${T_d}$ under the same OP combinations is depicted in Fig. 2(b). Similar to ${R_d}$, the range of ${T_d}$ rapidly declines and the degree of dispersion of the data drops sharply with the reduction in the physical thickness of the tissue.

 figure: Fig. 2.

Fig. 2. Distribution of (a) the diffuse reflectance ${R_d}$ and (b) the diffuse transmittance ${T_d}$ as a function of the absorption coefficient ${\mu _a}$, scattering coefficient ${\mu _s}$, and anisotropy factor g for different physical thicknesses ($d$ = 50, 100, 300, 500 µm), represented by colored dots reflecting the respective values.

Download Full Size | PDF

To further explore the impact of the physical thickness of the tissue sample on the accuracy of estimating ${\mu _a}$, ${\mu _s}$, g and ${n_s}$, respectively, sensitivity was introduced as an evaluation index, thereby helping determine the response of the change of the OPs to perturbations in the values of the spectral information. Moreover, it can assist in identifying if unique solutions for evaluating all OPs exist as well as those parameters that can be derived when it is not possible to uniquely obtain all parameters from the measurements [36]. CUDAMCML was used to calculate the mapping between ${R_d}$, ${T_d}$ and each OP of the tissue sample for different physical thicknesses with the benchmark settings being ${\mu _a}$= 1 cm-1, ${\mu _s}$= 200 cm-1, $g$ = 0.9, ${n_s}$ = 1.4 (i.e., the values of OPs when one of the OPs changes). The results are shown in Fig. 3. As the physical thickness of the tissue sample decreases, the curves (Figs. 3(a)–3(c)) of ${R_d}$ and ${T_d}$ for each OP except ${n_s}$ suddenly grow flat, indicating that the sensitivity (i.e., ${{\Delta {\mu _a}} / {\Delta {R_d}}}$ and ${{\Delta {\mu _a}} / {\Delta {T_d}}}$, etc.) of the OPs to the changes in ${R_d}$ and ${T_d}$ rise sharply when the tissue sample is relatively thin. In other words, a small change in spectral information caused by noise will lead to a great change in the results of obtained OPs. In addition, the OPs are more sensitive to variations in ${R_d}$ and ${T_d}$ in the region of high absorption, high scattering, and low anisotropy factor when thickness of the tissue sample is known, thus affecting the accuracy of determinated

 figure: Fig. 3.

Fig. 3. Mapping between the diffuse reflectance ${R_d}$ and diffuse transmittance ${T_d}$, and (a) absorption coefficient ${\mu _a}$, (b) scattering coefficient ${\mu _s}$, (c) anisotropy factor g, and (d) refractive index of tissue sample ${n_s}$, for different physical thicknesses d.

Download Full Size | PDF

OPs in the above range. Moreover, the curves (Fig. 3(d)) of ${R_d}$ and ${T_d}$ for ${n_s}$ no longer satisfy the monotonicity with the reduction in thickness, and ${n_s}$ may not be uniquely obtained from a set of ${R_d}$ and ${T_d}$. This is reasonable due to the mismatched boundary of the cuvette holder and tissue sample, and hence photons undergoing multiple scattering may deviate from their original direction by refraction or total reflection. To sum up, the reduction in the physical thickness of the tissue sample significantly increases the ill-conditioned nature of this inverse problem, resulting in a great challenge in the high-precision estimation of OPs.

3.2 Training results of CFNN

To train the neural network and ensure network convergence, the pre-computed MC database was randomly divided into three subsets: the training set (70%) for optimizing the model parameters, the test set (15%), and the validation set (15%) for adjusting the hyperparameters and preventing potential overfitting, only the test set was used to evaluate the network performance. As the IAD and traditional IMC methods are only applicable to tissue samples with known refractive index [37], the prediction results of OPs other than ${n_s}$ in the test set of the CFNN model were compared with those of the IMC, IAD, and MLP methods without ${T_c}$ and ${n_c}$ after training. The error thresholds of the IMC and IAD methods were both set to 0.1% [8]. Figure 4 shows the distribution of the prediction errors of the IAD, IMC, MLP, and CFNN methods for each OP, where the abscissa is the common logarithm of RE. The results indicated that the IAD method fail to extract ${\mu _a}$ from thin ex vivo tissues, and its relative errors were mostly above 100%. In contrast, the CFNN achieved the highest accuracy in the prediction of the three remaining OPs with relative errors below 0.01% in most cases. Moreover, compared with the MLP method, a significant improvement in the performance of the CFNN model was achieved by introducing ${T_c}$ and ${n_c}$ as additional parameters. More extractable features make it possible to eliminate the bottleneck of high sensitivity in the determination of OPs from thin ex vivo tissues, while simplifying the ill-posedness of the inverse problem so that solutions can be determined.

 figure: Fig. 4.

Fig. 4. Distribution of estimation REs of (a) the absorption coefficient ${\mu _a}$, (b) scattering coefficient ${\mu _s}$, and (c) anisotropy factor g, using the IAD, IMC, MLP, and CFNN model.

Download Full Size | PDF

To evaluate the generalization ability of the CFNN model, a new test set of 2000 OP combinations and the corresponding diffuse reflectance ${R_d}$, diffuse transmittance ${T_d}$, and collimated transmittance ${T_c}$ values were randomly generated by CUDAMCML. The estimated versus the expected OP values are shown in Fig. 5. Each black dot represents an OP data point and the red line represents the regression line. The MRE, R2 and linear regression equations are shown in each subplot. It was observed that the MRE in the mappings were 0.2247%, 0.0731%, and 0.0699% for ${\mu _a}$, ${\mu _s}$, and g, respectively. Meanwhile, the R2 is as high as 0.9999, and the determination is highly correlated in every case, further indicating that the CFNN model proposed in this study is capable of accurately and efficiently predicting OP mapping from thin ex vivo tissues.

 figure: Fig. 5.

Fig. 5. Estimated versus expected for the CFNN model: (a) absorption coefficient, (b) scattering coefficient, and (c) anisotropy factor. The black dots indicate the randomly generated OPs, and the red lines are the linear fit.

Download Full Size | PDF

To explore the estimation ability of the CFNN model for the refractive index of the tissue ${n_s}$, the prediction results of the test set were compared to those obtained using a customized IMC. The customized IMC uses CUDAMCML as a forward model of light propagation in biological tissues, and the Nelder-Mead simplex method detailed in [7] was adopted as the optimization algorithm with an error threshold of 0.1%. The prediction results of IMC versus CFNN for the refractive index ${n_s}$ are shown in Fig. 6(a), where the abscissa is the common logarithm of RE. Comparable to the gold standard IMC model in terms of estimation accuracy, the CFNN model achieved an RE below 1% in the vast majority of cases, as well as an MRE of 0.6645%. Meanwhile, 2000 sets of randomly generated CUDAMCML simulation results were again used to test the generalization ability of the refractive index ${n_s}$ derived from the CFNN model. Figures 6(b)–6(c) show the results. The red line in Fig. 6(b) represents the fitting curve of RE, and each blue dot represents the estimated error value. The results demonstrate that the CFNN model can successfully extract the refractive index of tissue samples, with the prediction error of most ${n_s}$ being less than 3% and the expected MRE being 2.3608%.

 figure: Fig. 6.

Fig. 6. (a) Distribution of estimation REs of the refractive index of the tissue sample ${n_s}$ obtained using the IMC and CFNN. (b) Prediction error for ${n_s}$ obtained using the CFNN. The blue dots indicate the randomly generated OPs, and the red line is the fitting curve of RE. (c) Boxplot of estimation error for ${n_s}$ obtained using the CFNN.

Download Full Size | PDF

3.3 Robustness of CFNN

Noise is ubiquitous in electronic components, networks, and systems. To explore the impact of noise on OP mapping determined by the CFNN model, different levels of Gaussian noise were added to the DIS measurements, and MRE and R2 were calculated to assess the model performance, with the results shown in Fig. 7. It can be seen that the value of R2 decreased as the noise level increased for all OPs, while the MRE showed an opposite trend. At a Gaussian noise level of 5%, the R2 values of the individual OPs absorption coefficient ${\mu _a}$, scattering coefficient ${\mu _s}$, and anisotropy factor g are all greater than 0.99, with the MRE of 5.5846%, 0.5047%, and 0.2028%, respectively, indicate that the CFNN network is robust to noise with no signs of overfitting. Additionally, it is worth noting that the prediction results for high ${\mu _s}$ and low g are inferior to those for low ${\mu _s}$ and high g, which is consistent with the sensitivity analysis results shown in Fig. 3. However, overall the MRE is as low as 2.8592% and 1.0880%, and there is still a strong correlation between the estimated results and the expected values even when the noise level was as high as 10%. Moreover, it was also found that the estimations for accuracy of ${\mu _s}$ and g modestly outperformed those for ${\mu _a}$ at a noise level of 10%, which is reasonable since the ${\mu _s}$ of biological tissue tends to be 1–2 orders of magnitude higher than ${\mu _a}$. Consequently, the CFNN model was better able to resist noise for determining ${\mu _s}$ than ${\mu _a}$ due to the fact that the contribution of ${\mu _a}$ to the DIS measurements is far less than that of ${\mu _s}$, which is consistent with the results in previous studies [18,34]. However, the possible noise of the DIS system and others is only 1% [18], and hence the CFNN network is fully competent for the high-precision prediction tasks of the four basic OPs from the DIS system.

 figure: Fig. 7.

Fig. 7. Estimated values using CFNN model and different noise levels: (a) and (b) absorption coefficient ${\mu _a}$, (c) and (d) scattering coefficient ${\mu _s}$, and (e) and (f) anisotropy factor g.

Download Full Size | PDF

3.4 Computational efficiency

Because the MC simulation relies on repeated random sampling, the data required to train the CFNN model took an entire day to run on a computer with an Intel Core i7-12700 H CPU and NVIDIA GeForce RTX 3070 Ti Laptop GPU. However, owing to the simplicity of the CFNN model, once the required training set was obtained, it only took less than five minutes to complete training on the same computer. The estimation of a single OP can be executed very quickly by a well-trained CFNN model, taking only 16.1 ms to calculate a batch size of 10,000 (average execution time 1.61 µs / parameter). Table 2 summarizes the optimization algorithms, computer configuration, and time required to predict a set of OPs for the CFNN, IMC, and IAD models used in this study. Although the time for a single MC simulation is reduced to 0.12 s by CUDAMCML, it is still excessively time-consuming to solve the inverse problem iteratively. However, the CFNN model is approximately five orders of magnitude faster than the IAD model and has an equivalent or improved accuracy compared to IMC, thus enabling a highly accurate and fast evaluation of four basic OPs of biological tissue within a wide spectral range.

Tables Icon

Table 2. Computational Efficiency of Inverse Algorithms

4. Conclusion

In summary, we investigated the potential of the CFNN model for the fast and accurate estimation of the basic OPs comprising the absorption coefficient ${\mu _a}$, scattering coefficient ${\mu _s}$, anisotropy factor g, and refractive index of tissues sample ${n_s}$. This is the first attempt to fully extract all basic OPs from a combination of the DIS system and the neural network. To accurately determine the OPs, the impact of the physical thickness of the tissue sample on the prediction accuracy of the OPs was analyzed using an MC simulation. The sensitivity analysis results show that the ill-conditioned nature of this inverse problem increases sharply with the reduction in tissue thickness, resulting in an excessive increase in errors when using available inversion algorithms. However, thin ex vivo tissues have substantially greater homogeneity as well as more measurable quantities, when compared to thick ex vivo tissues, thus emerging as one of the most promising alternatives to the full estimation of the entire set of basic OPs. To fully estimate all basic OPs from thin ex vivo tissues, a dedicated CFNN was built for each OP. Moreover, the refractive index of the cuvette holder ${n_c}$ fitted by the Sellmeier equation was additionally introduced as an efficient input to the CFNN model based on the diffuse reflectance ${R_d}$, diffuse transmittance ${T_d}$, and collimated transmittance ${T_c}$ measured by the DIS system. Our proposed CFNN model overcomes the serious ill-conditioned limitations in the evaluation of OPs from thin ex vivo tissues with higher accuracy and speed than those used in previous studies. The mean relative error values of the absorption coefficient ${\mu _a}$, scattering coefficient ${\mu _s}$, anisotropy factor g, and refractive index of the tissue sample ${n_s}$ determined from the CFNN model were 0.2247%, 0.0731%, 0.0699%, and 2.3608%, respectively. Moreover, the estimation of OPs using the CFNN model can be accomplished with reasonable values in real time, even at a level of 10% Gaussian noise. In contrast to the traditional IMC and IAD methods, the CFNN model can estimate the OPs approximately five orders of magnitude faster without loss of accuracy. This study provides a new implementation for highly accurate and fast prediction of four basic OPs of biological tissues within a wide spectral range and is promising for the optical diagnosis of diseases and quantitative analysis of laser surgery. In addition, the CFNN-based method has shown great potential in reflectance-based diffuse optical tomography (DOT), which was used to non-invasively visualize the breast tissues and their anomalies [24]. In future work, we will extend the CFNN model to other reflective geometry [38,39] to validate its performance of resolving multiple layers in vivo.

Funding

STI2030-Major Projects (2021ZD0204402); National Natural Science Foundation of China (61775066); Independent Innovation Fund of WNLO.

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

References

1. S. L. Jacques, “Optical properties of biological tissues: a review,” Phys. Med. Biol. 58(11), R37–R61 (2013). [CrossRef]  

2. K. Liu, Y. Jiang, W. Li, H. Chen, Q. Luo, and Y. Deng, “High-fidelity mesoscopic fluorescence molecular tomography based on SSB-Net,” Opt. Lett. 48(2), 199–202 (2023). [CrossRef]  

3. K. Sudo, Y. Shimojo, T. Nishimura, and K. Awazu, “Three-dimensional transient simulation of CO2 laser tissue vaporization and experimental evaluation with a hydrogel phantom,” J. Innov. Opt. Health Sci. 15(03), 2250016 (2022). [CrossRef]  

4. E. A. Genina, A. N. Bashkatov, D. K. Tuchina, P. A. Dyachenko, N. Navolokin, A. Shirokov, A. Khorovodov, A. Terskov, M. Klimova, A. Mamedova, I. Blokhina, I Agranovich, E Zinchenko, O. V. Semyachkina-Glushkovskaya, and V. V. Tuchin, “Optical properties of brain tissues at the different stages of glioma development in rats: pilot study,” Biomed. Opt. Express 10(10), 5182–5197 (2019). [CrossRef]  

5. A. ul Rehman, I. Ahmad, and S. A. Qureshi, “Biomedical applications of integrating sphere: a review,” Photodiagn. Photodyn. Ther. 31, 101712 (2020). [CrossRef]  

6. R. Lu, R. Van Beers, W. Saeys, C. Li, and H. Cen, “Measurement of optical properties of fruits and vegetables: a review,” Postharvest Biol. Technol. 159, 111003 (2020). [CrossRef]  

7. S. A. Prahl, “Everything I think you should know about inverse adding-doubling,” Oregon Medical Laser Center, St. Vincent Hospital 1, 1–74 (2011).

8. M. Friebel, A. Roggan, G. J. Müller, and M. C. Meinke, “Determination of optical properties of human blood in the spectral range 250 to 1100 nm using Monte Carlo simulations with hematocrit-dependent effective scattering phase functions,” J. Biomed. Opt. 11(3), 034021 (2006). [CrossRef]  

9. B. Morales Cruzado, S. V. y Montiel, and J. A. D. Atencio, “Genetic algorithms and MCML program for recovery of optical properties of homogeneous turbid media,” Biomed. Opt. Express 4(3), 433–446 (2013). [CrossRef]  

10. D. Hu, T. Sun, L. Yao, Z. Yang, A. Wang, and Y. Ying, “Monte Carlo: A flexible and accurate technique for modeling light transport in food and agricultural products,” Trends Food Sci. Technol. 102, 280–290 (2020). [CrossRef]  

11. B. Quistián-Vázquez, B. Morales-Cruzado, E. Sarmiento-Gómez, and F. G. Pérez-Gutiérrez, “Retrieval of absorption or scattering coefficient spectrum (RASCS) program: a tool to monitor optical properties in real time,” Lasers Surg. Med. 52, 552–559 (2020). [CrossRef]  

12. M. Hohmann, M. Späth, D. Ni, D. Dörner, B. Lengenfelder, F. Klämpfl, and M. Schmidt, “Random laser as a potential tool for the determination of the scattering coefficient,” Biomed. Opt. Express 12(9), 5439–5451 (2021). [CrossRef]  

13. D. Hu, R. Lu, Y. Huang, Y. Ying, and X. Fu, “Effects of optical variables in a single integrating sphere system on estimation of scattering properties of turbid media,” Biosyst. Eng. 194, 82–98 (2020). [CrossRef]  

14. F. Foschum, F. Bergmann, and A. Kienle, “Precise determination of the optical properties of turbid media using an optimized integrating sphere and advanced Monte Carlo simulations. Part 1: Theory,” Appl. Opt. 59(10), 3203–3215 (2020). [CrossRef]  

15. Y. Zhao, Y. Deng, F. Bao, H. Peterson, R. Istfan, and D. Roblyer, “Deep learning model for ultrafast multifrequency optical property extractions for spatial frequency domain imaging,” Opt. Lett. 43(22), 5669–5672 (2018). [CrossRef]  

16. M. Ivančič, P. Naglič, F. Pernuš, B. Likar, and M. Bürmen, “Efficient estimation of subdiffusive optical parameters in real time from spatially resolved reflectance by artificial neural networks,” Opt. Lett. 43(12), 2901–2904 (2018). [CrossRef]  

17. B. H. Hokr and J. N. Bixler, “Machine learning estimation of tissue optical properties,” Sci. Rep. 11, 1–7 (2021). [CrossRef]  

18. C. Li, H. Zhao, Q. Wang, and K. Xu, “Artificial neural network method for determining optical properties from double-integrating-spheres measurements,” Chin. Opt. Lett. 8(2), 173–176 (2010). [CrossRef]  

19. T. Nishimura, Y. Takai, Y. Shimojo, H. Hazama, and K. Awazu, “Determination of optical properties in double integrating sphere measurement by artificial neural network based method,” Opt. Rev. 28(1), 42–47 (2021). [CrossRef]  

20. K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” in Proc. IEEE Comput. Soc. Conf.770–778 (2016).

21. G. Huang, Z. Liu, L. van der Maaten, and K. Q. Weinberger, “Densely connected convolutional networks,” in Proc. IEEE Comput. Soc. Conf.4700–4708 (2017).

22. M. Alzayed, H. Chaoui, and Y. Farajpour, “Maximum power tracking for a wind energy conversion system using cascade-forward neural networks,” IEEE Trans. Sustain. Energy 12(4), 2367–2377 (2021). [CrossRef]  

23. M. H. Esfe and D. Toghraie, “Cascade forward artificial neural network to estimate thermal conductivity of functionalized graphene-water nanofluids,” Case Stud. Therm. Eng. 26, 101194 (2021). [CrossRef]  

24. G. M. Balasubramaniam and S. Arnon, “Regression-based neural network for improving image reconstruction in diffuse optical tomography,” Biomed. Opt. Express 13(4), 2006–2017 (2022). [CrossRef]  

25. R. Eldan and O. Shamir, “The power of depth for feedforward neural networks,” in Conference on Learning Theory, (PMLR, 2016), 907–940.

26. S. A. Prahl, M. J. van Gemert, and A. J. Welch, “Determining the optical properties of turbid media by using the adding–doubling method,” Appl. Opt. 32(4), 559–568 (1993). [CrossRef]  

27. Y. Arosa and R. de la Fuente, “Refractive index spectroscopy and material dispersion in fused silica glass,” Opt. Lett. 45(15), 4268–4271 (2020). [CrossRef]  

28. M. Mõttus, A. Hovi, and M. Rautiainen, “Theoretical algorithm and application of a double-integrating sphere system for measuring leaf transmittance and reflectance spectra,” Appl. Opt. 56(3), 563–571 (2017). [CrossRef]  

29. I. V. Yaroslavsky, A. N. Yaroslavsky, T. Goldbach, and H. J. Schwarzmaier, “Inverse hybrid technique for determining the optical properties of turbid media from integrating-sphere measurements,” Appl. Opt. 35(34), 6797–6809 (1996). [CrossRef]  

30. L. Wang, S. L. Jacques, and L. Zheng, “MCML—Monte Carlo modeling of light transport in multi-layered tissues,” Comput. Methods Programs Biomed. 47(2), 131–146 (1995). [CrossRef]  

31. E. Alerstam, T. Svensson, and S. Andersson-Engels, “Parallel computing with graphics processing units for high-speed Monte Carlo simulation of photon migration,” J. Biomed. Opt. 13(6), 060504 (2008). [CrossRef]  

32. S. J. Madsen, Optical Methods and Instrumentation in Brain Imaging and Therapy, (Springer Science & Business Media, 2012), Vol. 3.

33. M. Azimipour, R. Baumgartner, Y. Liu, S. L. Jacques, K. W. Eliceiri, and R. Pashaie, “Extraction of optical properties and prediction of light distribution in rat brain tissue,” J. Biomed. Opt. 19(7), 075001 (2014). [CrossRef]  

34. Z. Sun, L. Xie, D. Hu, and Y. Ying, “An artificial neural network model for accurate and efficient optical property mapping from spatial-frequency domain images,” Comput. Electron. Agric. 188, 106340 (2021). [CrossRef]  

35. D. Zhu, W. Lu, S. Zeng, and Q. Luo, “Effect of light losses of sample between two integrating spheres on optical properties estimation,” J. Biomed. Opt. 12(6), 064004 (2007). [CrossRef]  

36. H. Cen and R. Lu, “Quantification of the optical properties of two-layer turbid materials using a hyperspectral imaging-based spatially-resolved technique,” Appl. Opt. 48(29), 5612–5623 (2009). [CrossRef]  

37. D. Fukutomi, K. Ishii, and K. Awazu, “Determination of the scattering coefficient of biological tissue considering the wavelength and absorption dependence of the anisotropy factor,” Opt. Rev. 23(2), 291–298 (2016). [CrossRef]  

38. A. Kienle, M. S. Patterson, N. Dögnitz, R. Bays, G. Wagnières, and H. van den Bergh, “Noninvasive determination of the optical properties of two-layered turbid media,” Appl. Opt. 37(4), 779–791 (1998). [CrossRef]  

39. T. J. Farrell, M. S. Patterson, and M. Essenpreis, “Influence of layered tissue architecture on estimates of tissue optical properties obtained from spatially resolved diffuse reflectometry,” Appl. Opt. 37(10), 1958–1972 (1998). [CrossRef]  

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (7)

Fig. 1.
Fig. 1. Structure of CFNN model and workflow of OP mapping.
Fig. 2.
Fig. 2. Distribution of (a) the diffuse reflectance ${R_d}$ and (b) the diffuse transmittance ${T_d}$ as a function of the absorption coefficient ${\mu _a}$ , scattering coefficient ${\mu _s}$ , and anisotropy factor g for different physical thicknesses ( $d$  = 50, 100, 300, 500 µm), represented by colored dots reflecting the respective values.
Fig. 3.
Fig. 3. Mapping between the diffuse reflectance ${R_d}$ and diffuse transmittance ${T_d}$ , and (a) absorption coefficient ${\mu _a}$ , (b) scattering coefficient ${\mu _s}$ , (c) anisotropy factor g, and (d) refractive index of tissue sample ${n_s}$ , for different physical thicknesses d.
Fig. 4.
Fig. 4. Distribution of estimation REs of (a) the absorption coefficient ${\mu _a}$ , (b) scattering coefficient ${\mu _s}$ , and (c) anisotropy factor g, using the IAD, IMC, MLP, and CFNN model.
Fig. 5.
Fig. 5. Estimated versus expected for the CFNN model: (a) absorption coefficient, (b) scattering coefficient, and (c) anisotropy factor. The black dots indicate the randomly generated OPs, and the red lines are the linear fit.
Fig. 6.
Fig. 6. (a) Distribution of estimation REs of the refractive index of the tissue sample ${n_s}$ obtained using the IMC and CFNN. (b) Prediction error for ${n_s}$ obtained using the CFNN. The blue dots indicate the randomly generated OPs, and the red line is the fitting curve of RE. (c) Boxplot of estimation error for ${n_s}$ obtained using the CFNN.
Fig. 7.
Fig. 7. Estimated values using CFNN model and different noise levels: (a) and (b) absorption coefficient ${\mu _a}$ , (c) and (d) scattering coefficient ${\mu _s}$ , and (e) and (f) anisotropy factor g.

Tables (2)

Tables Icon

Table 1. Parameters of the integrating sphere

Tables Icon

Table 2. Computational Efficiency of Inverse Algorithms

Equations (16)

Equations on this page are rendered with MathJax. Learn more.

φ m ( 1 ) = b m ( 1 ) + i = 1 j w m i ( 1 ) x i
O m ( 1 ) = f m ( 1 ) ( φ m ( 1 ) )
φ m ( h ) = b m ( h ) + p = 1 h 1 k = 1 n O m k ( p ) w m k ( p ) + i = 1 j w m i ( h ) x i
O m ( h ) = f m ( h ) ( φ m ( h ) )
n 2 ( λ ) 1 =  B 1 λ 2 λ 2 C 1 + B 2 λ 2 λ 2 C 2 + B 3 λ 2 λ 2 C 3
T c = ( 1 r 1 ) ( 1 r 2 ) exp ( τ ) 1 r 1 r 2 exp ( 2 τ )
T d m = P t s P t t ,   R d m = P r s P r r
P r i = δ m R c i ( 1 R d d i s ) + R d i m α + ρ T d d i ( T d i + m T c i ) m α [ 1 ( m α + R d d i s ) ] [ 1 ( ρ T d d i ) 2 ] P i n ,   i = r , s
P t j = δ m T c j ( 1 R d d j s ) + T d j m α + ρ T d d j ( R d j + m R c j ) m α [ 1 ( m α + R d d j s ) ] [ 1 ( ρ T d d i ) 2 ] P i n ,   j = t , s
ρ = s 1 ( m α + R d d s )
R d d s = R d s ( 1 Δ r c s ) + Δ r c s
T d d s = T d s ( 1 Δ r c s )
R E = | X e s t i m a t e d X e x p e c t e d | X e x p e c t e d × 100 %
M R E = 1 n i = 1 n | X i , e s t i m a t e d X i , e x p e c t e d | X i , e x p e c t e d × 100 %
X a v g = 1 n i = 1 n X i , e x p e c t e d
R 2 = 1 i = 1 n ( X i , e s t i m a t e d X i , e x p e c t e d ) 2 i = 1 n ( X i , e s t i m a t e d X a v g ) 2
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.