OSA's Digital Library

Biomedical Optics Express

Biomedical Optics Express

  • Editor: Joseph A. Izatt
  • Vol. 5, Iss. 3 — Mar. 1, 2014
  • pp: 975–989
« Show journal navigation

Reliable recovery of the optical properties of multi-layer turbid media by iteratively using a layered diffusion model at multiple source-detector separations

Yu-Kai Liao and Sheng-Hao Tseng  »View Author Affiliations


Biomedical Optics Express, Vol. 5, Issue 3, pp. 975-989 (2014)
http://dx.doi.org/10.1364/BOE.5.000975


View Full Text Article

Acrobat PDF (1815 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

Accurately determining the optical properties of multi-layer turbid media using a layered diffusion model is often a difficult task and could be an ill-posed problem. In this study, an iterative algorithm was proposed for solving such problems. This algorithm employed a layered diffusion model to calculate the optical properties of a layered sample at several source-detector separations (SDSs). The optical properties determined at various SDSs were mutually referenced to complete one round of iteration and the optical properties were gradually revised in further iterations until a set of stable optical properties was obtained. We evaluated the performance of the proposed method using frequency domain Monte Carlo simulations and found that the method could robustly recover the layered sample properties with various layer thickness and optical property settings. It is expected that this algorithm can work with photon transport models in frequency and time domain for various applications, such as determination of subcutaneous fat or muscle optical properties and monitoring the hemodynamics of muscle.

© 2014 Optical Society of America

1. Introduction

Diffuse reflectance spectroscopy (DRS) is a noninvasive technique with good temporal resolution and rich spectral information. DRS has been applied to study the absorption and scattering properties of biological tissues, which can be further translated to clinically useful tissue functional information, such hemoglobin, water, and lipid concentrations [1

1. A. Cerussi, D. Hsiang, N. Shah, R. Mehta, A. Durkin, J. Butler, and B. J. Tromberg, “Predicting response to breast cancer neoadjuvant chemotherapy using diffuse optical spectroscopy,” Proc. Natl. Acad. Sci. U.S.A. 104(10), 4014–4019 (2007). [CrossRef] [PubMed]

, 2

2. S. H. Tseng, P. Bargo, A. Durkin, and N. Kollias, “Chromophore concentrations, absorption and scattering properties of human skin in-vivo,” Opt. Express 17(17), 14599–14617 (2009). [CrossRef] [PubMed]

]. The probing volume of DRS depends on the average photon travel lengths in tissue; thus, many studies manipulated the separation between source and detector to adjust the DRS probing volume. For instance, DRS with long and short SDSs (SDSs) was employed to study deep tissues, such as muscle and brain [3

3. A. Kienle and T. Glanzmann, “In vivo determination of the optical properties of muscle with time-resolved reflectance using a layered model,” Phys. Med. Biol. 44(11), 2689–2702 (1999). [CrossRef] [PubMed]

, 4

4. V. Toronov, A. Webb, J. H. Choi, M. Wolf, L. Safonova, U. Wolf, and E. Gratton, “Study of local cerebral hemodynamics by frequency-domain near-infrared spectroscopy and correlation with simultaneously acquired functional magnetic resonance imaging,” Opt. Express 9(8), 417–427 (2001). [CrossRef] [PubMed]

], and superficial tissues, such as skin [2

2. S. H. Tseng, P. Bargo, A. Durkin, and N. Kollias, “Chromophore concentrations, absorption and scattering properties of human skin in-vivo,” Opt. Express 17(17), 14599–14617 (2009). [CrossRef] [PubMed]

], respectively. In addition, DRS is a model based technique which involves solving an inverse problem by utilizing a nonlinear regression algorithm to fit the measurement data to a photon transport model to obtain the optical properties of samples. The adequacy of the photon transport model affects the accuracy of the recovered optical properties. For example, standard diffusion equation derived from radiative transport theory can be used to model the transport of photons that are in the diffusion regime [5

5. R. C. Haskell, L. O. Svaasand, T. T. Tsay, T. C. Feng, M. S. McAdams, and B. J. Tromberg, “Boundary Conditions for the Diffusion Equation in Radiative Transfer,” J. Opt. Soc. Am. A 11(10), 2727–2741 (1994). [CrossRef] [PubMed]

]. To reach this regime, photons have to experience sufficient scattering events and thus have long travel lengths (typically greater than 10 transport mean free paths) [6

6. S. H. Tseng, C. Hayakawa, J. Spanier, and A. J. Durkin, “Investigation of a probe design for facilitating the uses of the standard photon diffusion equation at short source-detector separations: Monte Carlo simulations,” J. Biomed. Opt. 14(5), 054043 (2009). [CrossRef] [PubMed]

]. In contrary, the standard diffusion model was shown to be ineffective in computing diffuse reflectance at SDSs shorter than 5mm [7

7. S. H. Tseng, C. Hayakawa, B. J. Tromberg, J. Spanier, and A. J. Durkin, “Quantitative spectroscopy of superficial turbid media,” Opt. Lett. 30(23), 3165–3167 (2005). [CrossRef] [PubMed]

]. Therefore, DRS works in conjunction with a standard diffusion equation at SDSs in the range from 10 to 30 mm has been successfully demonstrated for quantifying the optical properties of deep tissues [3

3. A. Kienle and T. Glanzmann, “In vivo determination of the optical properties of muscle with time-resolved reflectance using a layered model,” Phys. Med. Biol. 44(11), 2689–2702 (1999). [CrossRef] [PubMed]

, 8

8. B. J. Tromberg, N. Shah, R. Lanning, A. Cerussi, J. Espinoza, T. Pham, L. Svaasand, and J. Butler, “Non-invasive in vivo characterization of breast tumors using photon migration spectroscopy,” Neoplasia 2(1/2), 26–40 (2000). [CrossRef] [PubMed]

].

To investigate deep tissues using DRS, most studies employed a standard diffusion equation that assumed a semi-infinite sample geometry. Such semi-infinite assumption in modeling has been reported to introduce large errors for the determination of the optical properties of layered samples. Farrell et al. studied the influence of the top layer on the bottom layer optical properties recovery, and they reported that the semi-infinite assumption could lead to a recovery error higher than 100% [9

9. 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] [PubMed]

]. Because most biological tissues have a layered structure, diffusion models that take layered sample structure into account have been proposed. Martelli et al. used the eigenfunction method to solve the time-dependent diffusion equation for two- and three-layered samples [10

10. F. Martelli, A. Sassaroli, S. Del Bianco, Y. Yamada, and G. Zaccanti, “Solution of the time-dependent diffusion equation for layered diffusive media by the eigenfunction method,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 67(5), 056623 (2003). [CrossRef] [PubMed]

]. Liemert and Kienle found an analytical solution of the n-layered diffusion model in the Fourier domain [11

11. A. Liemert and A. Kienle, “Light diffusion in N-layered turbid media: frequency and time domains,” J. Biomed. Opt. 15(2), 025002 (2010). [CrossRef] [PubMed]

]. Although these multi-layered diffusion models can precisely describe photon transport in the diffusion regime, several studies indicated that using these models to determine the properties of a layered sample from measurements was still a difficult task [3

3. A. Kienle and T. Glanzmann, “In vivo determination of the optical properties of muscle with time-resolved reflectance using a layered model,” Phys. Med. Biol. 44(11), 2689–2702 (1999). [CrossRef] [PubMed]

, 12

12. A. Kienle, M. S. Patterson, N. Dognitz, R. Bays, G. Wagnieres, and H. van den Bergh, “Noninvasive determination of the optical properties of two-layered turbid media,” Appl. Opt. 37(4), 779–791 (1998). [CrossRef] [PubMed]

, 13

13. T. H. Pham, T. Spott, L. O. Svaasand, and B. J. Tromberg, “Quantifying the properties of two-layer turbid media with frequency-domain diffuse reflectance,” Appl. Opt. 39(25), 4733–4745 (2000). [CrossRef] [PubMed]

]. Pham et al. and Kienle et al. used layered diffusion models with frequency-domain and time resolved DRS system, respectively, to determine five parameters (absorption and reduced scattering coefficients of the first and second layers, and the thickness of the first layer) of two-layered samples. They both reported that simultaneous recovery of five parameters was possible only when the initial values provided for the inverse problem were very close to the true values; otherwise the errors of the recovered values would be unacceptably large [12

12. A. Kienle, M. S. Patterson, N. Dognitz, R. Bays, G. Wagnieres, and H. van den Bergh, “Noninvasive determination of the optical properties of two-layered turbid media,” Appl. Opt. 37(4), 779–791 (1998). [CrossRef] [PubMed]

, 13

13. T. H. Pham, T. Spott, L. O. Svaasand, and B. J. Tromberg, “Quantifying the properties of two-layer turbid media with frequency-domain diffuse reflectance,” Appl. Opt. 39(25), 4733–4745 (2000). [CrossRef] [PubMed]

]. As stated by Pham, the major reason for this phenomenon is that the inverse problem for the five-parameter fitting does not have a unique solution [13

13. T. H. Pham, T. Spott, L. O. Svaasand, and B. J. Tromberg, “Quantifying the properties of two-layer turbid media with frequency-domain diffuse reflectance,” Appl. Opt. 39(25), 4733–4745 (2000). [CrossRef] [PubMed]

]. The situation may become worse when the sample under investigation has a layer number larger than two.

2. Theoretical models

2.1 Layered diffusion model

Liemert and Kienle derived the solution of diffusion equation for N-layered turbid medium in frequency domain [11

11. A. Liemert and A. Kienle, “Light diffusion in N-layered turbid media: frequency and time domains,” J. Biomed. Opt. 15(2), 025002 (2010). [CrossRef] [PubMed]

]. The two-layered and three-layered photon diffusion models used in this study were adopted from their work and are summarized in the following. In a layered turbid medium system, the time dependent diffusion equation can be written as:
[1cit+μai[Di(r)]]Φi(r,t)=Si(r,t),
(1)
where D = 1/(3μa + 3μs') and Φ are the diffusion constant and the fluence rate, respectively. S is the source term, c is the speed of light in the medium, and i is the number of layer. The pencil beam light source in frequency domain can be approximated as a point source beneath the surface and is expressed as S1 = δ(x,y,z-z0), where z0 = 1/(μa + μs') is the location of the point source and ω = 2πf with f representing the source modulation frequency. By applying the extrapolated boundary condition and assuming the fluence and the flux are continuous at the boundary, the fluence rate of the diffusion equation system can be solved in the Fourier domain using the 2-D Fourier transform:

Φi(z,s1,s2,t)=Φi(x,y,z,t)exp[j(s1x+s2y)]dxdy.
(2)

The fluence rate at the air-sample interface for a layered sample whose bottom layer is infinitely thick has the following form in the Fourier domain:
Φ1(z,s)=sinh[α1(zb1+z0)]D1α1Z(z,s)N(z,s)sinh[α1(z0z)]D1α1,
(3)
,where s2=s12+s22, and zb1=2D1(1+Reff1)/(1+Reff1) with Reff1 representing the fraction of photons that is internally reflected at the air-sample boundary. For a two-layered sample, Z(z,s) and N(z,s) have the following form:
Z(z,s)=D1α1cosh[α1(l1z)]+(n22/n12)D2α2sinh[α1(l1z)],
(4)
N(z,s)=D1α1cosh[α1(l1+zb1)]+(n22/n12)D2α2sinh[α1(l1+zb1)].
(5)
For a three-layered sample, Z(z,s) and N(z,s) are:
Z(z,s)=D1α1[D2α2cosh(α2l2)+(n32/n22)D3α3sinh(α2l2)]cosh[α1(l1z)]+(n22/n12)D2α2[D2α2sinh(α2l2)+(n32/n22)D3α3cosh(α2l2)]×sinh[α1(l1z)],
(6)
N(z,s)=D1α1[D2α2cosh(α2l2)+(n32/n22)D3α3sinh(α2l2)]cosh[α1(l1+zb1)]+(n22/n12)D2α2[D2α2sinh(α2l2)+(n32/n22)D3α3cosh(α2l2)]×sinh[α1(l1+zb1)].
(7)
li and ni are the thickness and refractive index of the layer, respectively and αi2=(Dis2+μai+jω/ci)/Di. By performing inverse Fourier transform to Eq. (3) numerically, we can obtain the fluence rate of the diffuse reflectance. The spatially resolved reflectance can be calculated as the integral of the radiance L1 at the boundary, where L1=Φ1+3D1(Φ1/z)cosθ, over the backward hemisphere:
R(ρ)=2π[1Rfres(θ)]cosθ(L1/4π)dΩ.
(8)
Here ρ=x2+y2, and Rfres(θ) is the Fresnel reflection coefficient for a photon with an incident angle θ relative to the normal to the boundary.

2.2 Layered Monte Carlo model

3. Results and discussion

In this section, we will first explore the outcome of directly applying a layered diffusion model for solving the five parameters of two-layered samples in 3.1. To reliably solve the two-layer problem, an iterative algorithm was developed and will be described in detail in 3.2. Furthermore, the iterative algorithm was extended to solve the three-layer problem and will be discussed in 3.3.

3.1 Direct optical properties recovery of two-layered samples with the layered diffusion model

To test the layered diffusion model’s capability, we fit the Monte Carlo simulated frequency domain diffuse reflectance, including amplitude demodulation and phase delay, of two-layered samples to the diffusion model to simultaneously recover the five parameters of samples. In the Monte Carlo simulations, the dermis-fat structure was employed with the top dermis layer thickness of 1 mm. The five parameters of two-layered samples recovered using the layered diffusion model from 16 independent Monte Carlo simulations of 16 SDSs ranging from 5 to 20 mm are shown in Fig. 1
Fig. 1 Recovered (a) absorption coefficients and (b) reduced scattering coefficients of the top layer (squares) and bottom layer (triangles) at various source-detector separations. (c) Recovered top layer thickness L1 (asterisks) at various source-detector separations. Solid lines represent benchmark values.
. In this study, the initial value of the least-squares routine (MATLAB “lsqcurvefit” function with the Trust-Region-Reflective algorithm) at all SDSs was set to the benchmark values of the sample unless otherwise noted. It can be seen in Fig. 1(a) and 1(b) that, in general, the recovery percent errors of the bottom layer optical properties (maximal recovery percent errors are 20% and 32% for μa2 and μs2', respectively) are smaller than those of the top layer optical properties (maximal recovery percent errors are 71% and 49% for μa1 and μs1', respectively). The recovery error of the top layer thickness can be as high as 177%. It is found that the bottom layer optical properties and the top layer reduced scattering coefficient could be more accurately recovered at large SDSs than at short SDSs. Among the five parameters, the top layer absorption coefficient and the top layer thickness both have higher recovery errors than other parameters in general.

The fact that the optical properties of a two-layered sample cannot be reliably recovered using the frequency domain two-layered diffusion model at a single SDS can be further understood by Fig. 2
Fig. 2 χ2 computed at various top layer μa1 and μs1' combinations for (a) top layer thickness = 1 mm, and (b) top layer thickness = 1.5 mm, at source detector separation of 6 mm.
. To generate Fig. 2(a), we first ran a Monte Carlo simulation to calculate the frequency domain reflectance of a sample with dermis-fat structure (Y) using the following parameters: μa1 = 0.05 mm−1, μs1' = 1.5 mm−1, μa2 = 0.002 mm−1, μs2' = 1.0 mm−1, n = 1.43, top layer thickness = 1 mm, modulation frequency = 0 to 500 MHz with 1 MHz step, and the SDS = 6 mm. Next, by keeping the parameters other than μa1 and μs1' the same as those used in the Monte Carlo simulation, we utilized the two-layered diffusion model to calculate the frequency domain reflectance at various μa1 and μs1' combinations (Yj). Here the subscript j is the number of the μa1 and μs1' combinations, and it was set in the range from 1 to 2500 for plotting Fig. 2. Finally, we compute χ2 = (Y –Yj)2 at various μa1 and μs1' combinations and plot Fig. 2(a). In computing χ2 displayed in Fig. 2 and 3
Fig. 3 χ2 computed at various top layer μa1 and μs1' combinations for (a) top layer thickness = 1 mm, and (b) top layer thickness = 1.5 mm, at source detector separation of 13 mm.
, the amplitude magnitudes were magnified by the ratio of average phase to average amplitude, and we summed up the results at all modulation frequencies. It can be clearly seen in Fig. 2(a) that as there are only two variables to be fit (μa1 and μs1'), a least-squares fitting method can lead to a single solution for this two-layered sample problem since there is only one global minimum. It should be noted that the trough of the curve in Fig. 2(a) corresponds to (μa1 = 0.0449 mm−1, μs1' = 1.724 mm−1) which deviates from the benchmark top layer optical properties set in the Monte Carlo simulation. This suggests that there exists system differences between Monte Carlo method and the diffusion equation as indicated by Kienle et al. [12

12. A. Kienle, M. S. Patterson, N. Dognitz, R. Bays, G. Wagnieres, and H. van den Bergh, “Noninvasive determination of the optical properties of two-layered turbid media,” Appl. Opt. 37(4), 779–791 (1998). [CrossRef] [PubMed]

] and the diffusion equation has non-negligible errors even in the diffusion regime. Following similar procedure in creating Fig. 2(a), we generated Fig. 2(b) by changing top layer thickness to 1.5 mm and keeping other parameters in the diffusion model the same as those used for generating Fig. 2(a). Note that we employed the same benchmark Monte Carlo simulation data (Y) for generating Fig. 2(a) and 2(b). From Fig. 2(b), we can find that the value of (Yi -Yi)2 at the global minimum is 0.4272 which is smaller than that of Fig. 2(a) (0.4558). This result implies that as the free parameters increase to 3 (μa1, μs1', and top layer thickness), a least-squares fitting routine may find a solution that greatly deviates from the benchmark values for this two-layer problem. In addition, we found that increasing the number of free parameters to 4 or 5 would lead to multiple local minima in the χ2 space for this two-layered sample problem (data not shown) and thus a least-squares fitting routine employed for solving this problem would find different answers for different initial value settings.

We also investigated the effect of simultaneously employing multiple SDSs on the accuracy of the recovered optical properties of a dermis-fat structured sample. We tested various combinations of SDSs and part of the results are listed in Table 1

Table 1. Optical properties of top (μa1, μs1') and bottom (μa2, μs2') layers as well as top layer thickness (L1) recovered at various source-detector separation combinations for a dermis-fat structured sample. Recovery errors are listed below the values.

table-icon
View This Table
| View All Tables
. In Table 1, data recovery results listed were obtained with 2, 3, and 4 source-detector pairs. All three combinations include the relatively short (6 mm) and long (20 mm) SDSs, and two of them also include intermediate SDSs. It can be seen that μs1', μa2 and μs2' have recovery errors less than 10%, while the recovery errors for μa1 and top layer thickness are greater than 25.2% and 192.6%, respectively. In general, employing multiple source-detector pairs produced better recovery results than the single source-detector pair counterpart. For reference, the errors of (μa1, μs1', μa2, μs2', L1) recovered at 6 mm SDS were −81.0%, 46.5%, 5.0%, −39.4%, and 33.3%, respectively. However, the benefit of increasing the source-detector pair number to more than 2 is found only for μs1' and μs2' but is not evident for other parameters especially the top layer thickness. Our results indicate that the two-layered sample problem cannot be reliably solved by simultaneously employing multiple source-detector pairs.

3.2 Optical properties recovery of a two-layered sample by iteratively using a layered diffusion model at multiple SDSs

3.2.1 Iterative algorithm and selection of SDSs

As investigating the two-layered sample problem, we discovered an interesting fact that the top layer thickness could be better recovered at some SDSs than others. This fact can be seen in Fig. 1(c) where the top layer thickness can be more accurately recovered at intermediate SDSs roughly ranging from 8 to 13 mm than at other SDSs. Besides, in this intermediate range, top layer absorption and scattering properties can sometimes be properly recovered. We generated plots in a similar way to those demonstrated in Fig. 2 to understand this phenomenon. The reference data was generated from our Monte Carlo model to simulate the dermis-fat structure with top layer thickness = 1 mm and SDS = 13 mm. We then calculated the χ2 between the reference data and the reflectance determined from the two-layered diffusion model at various (μa1, μs1') combinations to generate Fig. 3. The top layer thickness in the diffusion model was set at 1mm for Fig. 3(a) and 1.5 mm for Fig. 3(b). It can be seen in Fig. 3 that χ2 is smaller at 1 mm than at 1.5 mm which is in contrary to the trend demonstrated in Fig. 2. In fact, as the bottom layer optical properties were fixed at the benchmark values and there were three floating variables, the minimum of the χ2 for SDS = 13 mm was found at L1 = 1.028 mm, μa1 = 0.057 mm−1, and μs1' = 1.89 mm−1. This implied that top layer properties could be properly estimated at intermediate SDSs as the bottom layer optical properties were known.

Furthermore, it can be noticed in Fig. 2(a) that when the top layer thickness is close to 1 mm and the SDS = 6 mm, the best top layer optical properties that minimize χ2 are near μa1 = 0.045 mm−1, μs1' = 1.72 mm−1, which are closer to the benchmark top layer optical properties than those determined at SDS = 13 mm (μa1 = 0.057 mm−1, μs1' = 1.89 mm−1). This suggests that if the top layer thickness is fixed at a value that is close to the true value, 6 mm SDS measurement data could lead to a better top layer optical properties recovery than the 13 mm counterpart. It can thus be inferred that measurements conducted at short and intermediate SDSs are more sensitive to the variations of top layer optical properties and top layer thickness, respectively, than those conducted at other SDSs.

The algorithm described above was used to recover the properties of a dermis-fat sample with 1 mm top layer thickness, and the results are shown in Fig. 5
Fig. 5 Recovered (a) absorption coefficients and (b) reduced scattering coefficients of the top layer (squares) and bottom layer (triangles), as well as (c) recovered top layer thickness L1 (asterisks) versus iteration number for a dermis-fat structure. Solid lines represent benchmark values.
. Three distances 6, 13, and 20 mm were chosen to cover the short, intermediate, and long SDS regimes, respectively. We can see that the bottom layer optical properties were quite stable throughout the iterative calculation process, while the other three parameters gradually approached the values that were near the benchmark values. The iteration algorithm was terminated under the following conditions: 1. One of the determined top layer properties of this iteration had less than 0.1% variation from those of the last iteration. (convergence termination) 2. One of the determined top layer properties showed a reverse trend in the iteration and the variation was larger than 1%. (oscillation/divergence termination) Termination condition 1 usually implies successful iteration process is achieved while termination condition 2 suggests that the calculated values are possibly oscillating or starting to diverge. When the variation percentage setting is judiciously adjusted in the step 3 of an iteration described above, termination condition 2 is satisfied most likely due to oscillating results rather than diverging results. To avoid the algorithm from premature termination, the algorithm was set to have a minimal iteration number of 5. The iteration process for obtaining the results demonstrated in Fig. 5 was terminated because the determined top layer thickness satisfied the termination condition 2. It is worth noting that the recovered values of the top layer properties in the first round of iteration were not close to the benchmark values which suggested that our algorithm was robust. In addition, we also found that the algorithm was not sensitive to the initial value settings.

The final recovery results of Fig. 5 are listed in the first row of Table 2

Table 2. Optical properties of top (μa1, μs1') and bottom (μa2, μs2') layers as well as top layer thickness (L1) recovered at various source-detector separation combinations for a dermis-fat structured sample. Recovery errors are listed below the values.

table-icon
View This Table
| View All Tables
. To understand the effect of measurement noise on the recovered results, we superimposed additional 5% amplitude and 0.6° phase noises to the Monte Carlo data to simulate the measurement noise of real frequency domain reflectance systems [13

13. T. H. Pham, T. Spott, L. O. Svaasand, and B. J. Tromberg, “Quantifying the properties of two-layer turbid media with frequency-domain diffuse reflectance,” Appl. Opt. 39(25), 4733–4745 (2000). [CrossRef] [PubMed]

]. The percent errors of the recovered μa1, μs1', μa2, μs2', and L1 were −11.5%, 15.7%, 5.0%, −4.0%, and 11.2%, respectively, which are comparable to the values listed in the second row of Table 2. This result indicates that the iterative algorithm is still robust under typical measurement noise.

Moreover, we investigated the influence of the choice of SDSs on the recovered results. The sample setup here was again a dermis-fat structure with 1 mm top layer thickness. In making the combination of the SDSs, we kept the longest SDS fixed at 20 mm and adjusted the shortest SDS. The intermediate SDS was approximately set at a position which was close to the average of the longest and the shortest SDSs. The results are listed in Table 2 and it can be seen that all recovered values are similar as the shortest SDS varies from 6 to 8 mm, and the (8, 14, 20) mm combination produced the best top layer results among them. As the shortest SDS increases to 12 and 14 mm, the errors of the recovered top layer properties grow accordingly. This means that the shortest SDS has to be short enough to enhance the accuracy of the top layer properties recovery. It is worth noting that the bottom layer optical properties were properly recovered for all combinations with 20 mm SDS. In the last set of data shown in Table 2, we reduced the longest and intermediate SDSs. The recovery errors for the bottom layer properties of this data set were larger than those of the other data sets. Larger bottom layer optical properties recovery errors as well as a not optimized intermediate SDS setting led to the great deviation of the recovered top layer properties from the benchmark values in the last data set.

3.2.2 Influence of the layer thickness

3.2.3 Hemodynamics

Several researchers studied the muscle optical properties as well as muscle hemodynamics using a semi-infinite diffusion model and reported compromised accuracy of the recovered muscle optical properties due to the presence of dermis or fat [3

3. A. Kienle and T. Glanzmann, “In vivo determination of the optical properties of muscle with time-resolved reflectance using a layered model,” Phys. Med. Biol. 44(11), 2689–2702 (1999). [CrossRef] [PubMed]

, 15

15. G. Alexandrakis, D. R. Busch, G. W. Faris, and M. S. Patterson, “Determination of the optical properties of two-layer turbid media by use of a frequency-domain hybrid monte carlo diffusion model,” Appl. Opt. 40(22), 3810–3821 (2001). [CrossRef] [PubMed]

]. This problem can be approximated to a dermis-muscle structure problem and solved with a two-layered diffusion model. We investigated that whether the bottom layer absorption coefficient recovered using the iterative algorithm can faithfully reflect the variation of the muscle layer absorption property. To this end, in the Monte Carlo simulations, the bottom layer absorption coefficient was varied from 0.05 to 0.07 mm−1 to simulate muscle hemodynamics and the other parameters were kept fixed. The layer thickness of dermis was 1 mm. The recovery results are shown in Fig. 9
Fig. 9 Recovered (a) absorption coefficients and (b) reduced scattering coefficients of the top layer (squares) and bottom layer (triangles), as well as (c) recovered top layer thickness L1 (asterisks) versus top layer thickness for a dermis-muscle structure. Sample optical properties recovered using a semi-infinite diffusion model are depicted as empty circles. Solid lines represent benchmark values.
. In Fig. 9(a), the recovered muscle absorption coefficients are systematically lower than the benchmark values by 0.0039 mm−1 in average (6% error in average) and the slope of the regression line of the two sets of data is 1.069. Our results suggest that muscle hemodynamics could be tracked using our method with a modest error. The other four parameters plotted in Fig. 9(a) to 9(c) do not show clear dependence with the bottom layer absorption coefficient.

In addition, using the same set of two-layered Monte Carlo data, we calculated the optical properties of sample using a frequency domain semi-infinite diffusion model at SDS of 20 mm, and the results are shown as empty circles in Fig. 9(a) and 9(b). It can be seen that although the recovered absorption coefficients recovered with a semi-infinite diffusion model vary in accordance with the benchmark values, they remarkably deviate from the benchmark values by 0.023 mm−1.On the other hand, the reduced scattering coefficients recovered using the semi-infinite model do not show strong dependence on the variation of bottom layer absorption, but they are less accurate than those recovered using the iterative algorithm.

3.3 Optical properties recovery of a three-layered sample by iteratively using a layered diffusion model at several SDSs

3.3.1 Iterative algorithm and selection of SDSs

The skin of trunk and limbs usually sits on a fat-muscle structure which forms a three-layered sample model. The thickness of fat layer is highly body-site, gender, and subject dependent, and can vary between 1 to 70 mm [17

17. S. Leahy, C. Toomey, K. McCreesh, C. O’Neill, and P. Jakeman, “Ultrasound measurement of subcutaneous adipose tissue thickness accurately predicts total and segmental body fat of young adults,” Ultrasound Med. Biol. 38(1), 28–34 (2012). [CrossRef] [PubMed]

]. If the fat layer is very thin, the three layer problem can be approximated as a dermis-muscle two-layer problem. On the other hand, if the fat layer is extremely thick, the three-layer problem can be treated as a dermis-fat two-layer problem. When the fat layer has a moderate thickness, the three-layer problem is practical that one needs to consider when studying muscle hemodynamics using DRS techniques.

The SDS combination of (6, 10, 13, 17, 20) mm does not work for all kinds of three-layered sample structures such as the structure with thick middle layer. For such cases, larger intermediate SDSs may be required to obtain reasonable results for the middle layer. For example, we ran another set of Monte Carlo simulation where the middle fat layer thickness was increased to 10 mm and kept other parameter values untouched. By using the SDS combination of (6, 10, 13, 17, 20) mm, the recovered (μa3, μs3') had percent errors of −62.7% and 93.9%, respectively. When the SDS combination was adjusted to (6, 10, 20, 25, 35) mm, the errors of recovered (μa3, μs3') were reduced to −36.0% and −34.1%, respectively. Therefore, picking a proper combination of SDS based on the structure of the site under investigation is the prerequisite for obtaining accurate properties of three-layer samples.

3.3.2 Hemodynamics

As mentioned earlier that the three-layered diffusion model has its significance in studying muscle hemodynamics. Thus, it is crucial to understand that whether the iterative algorithm with a three-layered diffusion model can be used to track the hemodynamics of muscle tissue in a three-layered tissue structure. Monte Carlo simulations were carried out to generate the experimental data. In the Monte Carlo simulations, the thicknesses of dermis and fat layers were 1 and 5 mm, respectively, and the absorption coefficient of muscle layer was varied from 0.05 to 0.07 mm−1. The recovered muscle layer absorption and reduced scattering coefficients using the SDS combination of (6, 10, 13, 17, 20) mm are shown in Fig. 10
Fig. 10 Recovered (a) absorption coefficients and (b) reduced scattering coefficients of the muscle layer (triangles) at various muscle layer absorption for a dermis-fat-muscle structure. Crosses are the results obtained with fixed top and middle layer thicknesses. Sample optical properties recovered using a semi-infinite diffusion model are depicted as empty circles. Solid lines represent benchmark values.
as filled triangles. In Fig. 10(a), the recovered muscle absorption coefficients increase with the benchmark values and the linear regression slope is 0.538. Also depicted in Fig. 10 as empty circles are the recovered muscle optical properties using the semi-infinite diffusion model at 20 mm SDS. It can be seen that the muscle optical properties determined using the iterative algorithm are closer to the benchmark values than those determined using the semi-infinite model. The absorption coefficients calculated using the semi-infinite model barely follow the variation trend of the benchmark values and have a linear regression slope of 0.145. Our simulation results suggest that the semi-infinite model is not suitable for tracking muscle hemodynamics in a three-layered tissue model.

4. Conclusion

We observed that the measurements conducted at short and long SDSs were more suitable for recovering the optical properties of top and bottom layers of a layered sample, respectively. Based on our observation, we developed an algorithm in an attempt to accurately determine the properties of layered samples. The iteratively algorithm revealed in this study is capable of solving the optical properties and layer thicknesses of layered samples. We tested this algorithm using Monte Carlo simulations and found that it was robust and could reliably recover the optical properties and layer thicknesses of s two-layered or three-layered sample structures. For a two-layered dermis-fat or dermis-muscle structure with maximum top layer (dermis) thickness of 3 mm, we found that (6, 13, 20) mm was a good SDS combination for the iterative algorithm. For all the two-layered skin cases we tested, the recovered optical properties and layer thickness had percent errors less than 25%. On the other hand, for the three-layered skin structure, it is difficult to determine the best SDS combination since the thickness of the middle fat layer is highly site and subject dependent and has a great variation. For example, as the top and middle layer thicknesses were 1 and 5 mm, respectively, we found that the reflectance data measured at a SDS combination of (6, 10, 13, 17, 20) mm was sufficient for accurately deriving sample optical properties and thicknesses. However, when the middle layer thickness was increased to 10mm, the longest SDS had to be increased to at least 35 mm to obtain reasonable bottom layer optical properties.

In this study, it typical required 1000-3000s to obtain the final results on an Intel i7-4930K based computer. We found that the computation time could be reduced to less than 60s by only employing 11 frequencies in the 0-500 MHz range; however, the recovery accuracy was slightly compromised. Due to the speed limitation, this algorithm may not be suitable for clinical settings that require real-time information of tissues.

We employed a frequency domain layered diffusion model to work with the iterative algorithm in this study; nevertheless, it is reasonable to expect that other photon transport models in frequency or time domain can also work with the algorithm. For example, photon diffusion models are generally known to be incapable of describing the photon transport in highly absorbing materials or at short SDSs, while Monte Carlo methods do not have such limitations. Employing Monte Carlo models together with the iterative algorithm could help determine the optical properties of highly absorbing layered samples or the thickness of samples with very thin layers. Furthermore, the algorithm could be extended to study samples that have layer number larger than three such as human or animal head models.

Acknowledgments

This research was supported by the National Science Council of Taiwan under Grant No. NSC 102-2221-E-006-248-MY2.

References and links

1.

A. Cerussi, D. Hsiang, N. Shah, R. Mehta, A. Durkin, J. Butler, and B. J. Tromberg, “Predicting response to breast cancer neoadjuvant chemotherapy using diffuse optical spectroscopy,” Proc. Natl. Acad. Sci. U.S.A. 104(10), 4014–4019 (2007). [CrossRef] [PubMed]

2.

S. H. Tseng, P. Bargo, A. Durkin, and N. Kollias, “Chromophore concentrations, absorption and scattering properties of human skin in-vivo,” Opt. Express 17(17), 14599–14617 (2009). [CrossRef] [PubMed]

3.

A. Kienle and T. Glanzmann, “In vivo determination of the optical properties of muscle with time-resolved reflectance using a layered model,” Phys. Med. Biol. 44(11), 2689–2702 (1999). [CrossRef] [PubMed]

4.

V. Toronov, A. Webb, J. H. Choi, M. Wolf, L. Safonova, U. Wolf, and E. Gratton, “Study of local cerebral hemodynamics by frequency-domain near-infrared spectroscopy and correlation with simultaneously acquired functional magnetic resonance imaging,” Opt. Express 9(8), 417–427 (2001). [CrossRef] [PubMed]

5.

R. C. Haskell, L. O. Svaasand, T. T. Tsay, T. C. Feng, M. S. McAdams, and B. J. Tromberg, “Boundary Conditions for the Diffusion Equation in Radiative Transfer,” J. Opt. Soc. Am. A 11(10), 2727–2741 (1994). [CrossRef] [PubMed]

6.

S. H. Tseng, C. Hayakawa, J. Spanier, and A. J. Durkin, “Investigation of a probe design for facilitating the uses of the standard photon diffusion equation at short source-detector separations: Monte Carlo simulations,” J. Biomed. Opt. 14(5), 054043 (2009). [CrossRef] [PubMed]

7.

S. H. Tseng, C. Hayakawa, B. J. Tromberg, J. Spanier, and A. J. Durkin, “Quantitative spectroscopy of superficial turbid media,” Opt. Lett. 30(23), 3165–3167 (2005). [CrossRef] [PubMed]

8.

B. J. Tromberg, N. Shah, R. Lanning, A. Cerussi, J. Espinoza, T. Pham, L. Svaasand, and J. Butler, “Non-invasive in vivo characterization of breast tumors using photon migration spectroscopy,” Neoplasia 2(1/2), 26–40 (2000). [CrossRef] [PubMed]

9.

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] [PubMed]

10.

F. Martelli, A. Sassaroli, S. Del Bianco, Y. Yamada, and G. Zaccanti, “Solution of the time-dependent diffusion equation for layered diffusive media by the eigenfunction method,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 67(5), 056623 (2003). [CrossRef] [PubMed]

11.

A. Liemert and A. Kienle, “Light diffusion in N-layered turbid media: frequency and time domains,” J. Biomed. Opt. 15(2), 025002 (2010). [CrossRef] [PubMed]

12.

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

13.

T. H. Pham, T. Spott, L. O. Svaasand, and B. J. Tromberg, “Quantifying the properties of two-layer turbid media with frequency-domain diffuse reflectance,” Appl. Opt. 39(25), 4733–4745 (2000). [CrossRef] [PubMed]

14.

L. H. Wang, S. L. Jacques, and L. Q. Zheng, “Mcml - Monte-Carlo Modeling of Light Transport in Multilayered Tissues,” Comput. Meth. Prog. Bio. 47(2), 131–146 (1995). [CrossRef]

15.

G. Alexandrakis, D. R. Busch, G. W. Faris, and M. S. Patterson, “Determination of the optical properties of two-layer turbid media by use of a frequency-domain hybrid monte carlo diffusion model,” Appl. Opt. 40(22), 3810–3821 (2001). [CrossRef] [PubMed]

16.

A. Liemert and A. Kienle, “Bioluminescent light diffusion in a four-layered turbid medium,” Med. Laser Appl. 25(3), 161–165 (2010). [CrossRef]

17.

S. Leahy, C. Toomey, K. McCreesh, C. O’Neill, and P. Jakeman, “Ultrasound measurement of subcutaneous adipose tissue thickness accurately predicts total and segmental body fat of young adults,” Ultrasound Med. Biol. 38(1), 28–34 (2012). [CrossRef] [PubMed]

OCIS Codes
(170.5270) Medical optics and biotechnology : Photon density waves
(170.5280) Medical optics and biotechnology : Photon migration
(290.1990) Scattering : Diffusion

ToC Category:
Optics of Tissue and Turbid Media

History
Original Manuscript: January 1, 2014
Revised Manuscript: February 9, 2014
Manuscript Accepted: February 21, 2014
Published: February 27, 2014

Citation
Yu-Kai Liao and Sheng-Hao Tseng, "Reliable recovery of the optical properties of multi-layer turbid media by iteratively using a layered diffusion model at multiple source-detector separations," Biomed. Opt. Express 5, 975-989 (2014)
http://www.opticsinfobase.org/boe/abstract.cfm?URI=boe-5-3-975


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. A. Cerussi, D. Hsiang, N. Shah, R. Mehta, A. Durkin, J. Butler, and B. J. Tromberg, “Predicting response to breast cancer neoadjuvant chemotherapy using diffuse optical spectroscopy,” Proc. Natl. Acad. Sci. U.S.A.104(10), 4014–4019 (2007). [CrossRef] [PubMed]
  2. S. H. Tseng, P. Bargo, A. Durkin, and N. Kollias, “Chromophore concentrations, absorption and scattering properties of human skin in-vivo,” Opt. Express17(17), 14599–14617 (2009). [CrossRef] [PubMed]
  3. A. Kienle and T. Glanzmann, “In vivo determination of the optical properties of muscle with time-resolved reflectance using a layered model,” Phys. Med. Biol.44(11), 2689–2702 (1999). [CrossRef] [PubMed]
  4. V. Toronov, A. Webb, J. H. Choi, M. Wolf, L. Safonova, U. Wolf, and E. Gratton, “Study of local cerebral hemodynamics by frequency-domain near-infrared spectroscopy and correlation with simultaneously acquired functional magnetic resonance imaging,” Opt. Express9(8), 417–427 (2001). [CrossRef] [PubMed]
  5. R. C. Haskell, L. O. Svaasand, T. T. Tsay, T. C. Feng, M. S. McAdams, and B. J. Tromberg, “Boundary Conditions for the Diffusion Equation in Radiative Transfer,” J. Opt. Soc. Am. A11(10), 2727–2741 (1994). [CrossRef] [PubMed]
  6. S. H. Tseng, C. Hayakawa, J. Spanier, and A. J. Durkin, “Investigation of a probe design for facilitating the uses of the standard photon diffusion equation at short source-detector separations: Monte Carlo simulations,” J. Biomed. Opt.14(5), 054043 (2009). [CrossRef] [PubMed]
  7. S. H. Tseng, C. Hayakawa, B. J. Tromberg, J. Spanier, and A. J. Durkin, “Quantitative spectroscopy of superficial turbid media,” Opt. Lett.30(23), 3165–3167 (2005). [CrossRef] [PubMed]
  8. B. J. Tromberg, N. Shah, R. Lanning, A. Cerussi, J. Espinoza, T. Pham, L. Svaasand, and J. Butler, “Non-invasive in vivo characterization of breast tumors using photon migration spectroscopy,” Neoplasia2(1/2), 26–40 (2000). [CrossRef] [PubMed]
  9. 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] [PubMed]
  10. F. Martelli, A. Sassaroli, S. Del Bianco, Y. Yamada, and G. Zaccanti, “Solution of the time-dependent diffusion equation for layered diffusive media by the eigenfunction method,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys.67(5), 056623 (2003). [CrossRef] [PubMed]
  11. A. Liemert and A. Kienle, “Light diffusion in N-layered turbid media: frequency and time domains,” J. Biomed. Opt.15(2), 025002 (2010). [CrossRef] [PubMed]
  12. A. Kienle, M. S. Patterson, N. Dognitz, R. Bays, G. Wagnieres, and H. van den Bergh, “Noninvasive determination of the optical properties of two-layered turbid media,” Appl. Opt.37(4), 779–791 (1998). [CrossRef] [PubMed]
  13. T. H. Pham, T. Spott, L. O. Svaasand, and B. J. Tromberg, “Quantifying the properties of two-layer turbid media with frequency-domain diffuse reflectance,” Appl. Opt.39(25), 4733–4745 (2000). [CrossRef] [PubMed]
  14. L. H. Wang, S. L. Jacques, and L. Q. Zheng, “Mcml - Monte-Carlo Modeling of Light Transport in Multilayered Tissues,” Comput. Meth. Prog. Bio.47(2), 131–146 (1995). [CrossRef]
  15. G. Alexandrakis, D. R. Busch, G. W. Faris, and M. S. Patterson, “Determination of the optical properties of two-layer turbid media by use of a frequency-domain hybrid monte carlo diffusion model,” Appl. Opt.40(22), 3810–3821 (2001). [CrossRef] [PubMed]
  16. A. Liemert and A. Kienle, “Bioluminescent light diffusion in a four-layered turbid medium,” Med. Laser Appl.25(3), 161–165 (2010). [CrossRef]
  17. S. Leahy, C. Toomey, K. McCreesh, C. O’Neill, and P. Jakeman, “Ultrasound measurement of subcutaneous adipose tissue thickness accurately predicts total and segmental body fat of young adults,” Ultrasound Med. Biol.38(1), 28–34 (2012). [CrossRef] [PubMed]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.


« Previous Article

OSA is a member of CrossRef.

CrossCheck Deposited