Method to improve the depth sensitivity of diffuse reflectance measurements to absorption changes in optically turbid medium.

We have studied the spatial distributions of the sensitivity of time-resolved near-infrared diffuse reflectance measurement. Sensitivity factors representing a change of parameters of a measured optical signal induced by absorption perturbation in a certain voxel of the medium were simulated using the diffusion equation solution. The parameters were statistical moments of measured distributions of time of flight of photons (DTOFs) i.e., the total number of photons, mean time of flight and variance. The distributions of the sensitivity of statistical moments of DTOFs to a change in absorption were generated for various source-detector separations and various optical properties of the medium. Furthermore, differential sensitivity distributions for two different source-detector separations were calculated. A measurement geometry, in which two detection spots, separated by 5 mm, in combination with two sources was proposed. For this setup differences between the signals obtained for both detectors were calculated independently for both sources and afterward summed up for both source positions. Obtained differences in moments of DTOFs assessed at two source-detector separations and summed up for different positioning of the sources allowed to shape up the sensitivity profiles. Calculated sensitivity profiles show that positive sensitivities of the mean time of flight of photons and variance of the DTOF can be obtained. These positive sensitivity areas are located just between both detection spots and cover the compartment located deeply in the medium. The sensitivity in superficial compartments of the medium is negative and much smaller in amplitude. The proposed technique can be used for improved discrimination of optical signals related to the intracerebral change in absorption which remains a serious obstacle in the application of the NIRS technique in the assessment of brain oxygenation or perfusion.


Introduction
In last decade the time-resolved measurement of diffuse reflectance became to be serious alternative to the continuous wave spectroscopy for the assessment of tissue oxygenation and perfusion [1]. This trend is related to the reducing costs of components as well as appearance of new technologies allowing for estimation of distributions of times of flight of photons. Moreover, it was presented in multiple studies that in some application continuous wave near infrared spectroscopy (cwNIRS) technique is insufficient because the results are strongly contaminated by an influence of extracerebral tissues [2][3][4][5][6][7]. The distance-resolved NIRS technique may allow to reduce such extracerebral contamination [8,9], however the time-resolved near infrared spectroscopy (trNIRS) technique still remains the most advanced approach. The trNIRS allows to gather more information about the tissue of interest, in comparison to cwNIRS technique [10][11][12]. It was reported that trNIRS allows to estimate optical properties of the tissue (absorption and reduced scattering coefficients µ a and µ s ') considered as a homogenous medium [13]. Furthermore, changes in absorption caused by the change in concentration of chromophores at can be assessed as a function of depth in the tissue [11,14]. In measurements of the cerebral perfusion or oxygenation the technique of near infrared spectroscopy is typically used in diffuse reflectance mode. It was shown in series of studies that the time-resolved technique allows to solve the problem of the contamination of the brain cortex oxygenation by components related to the extracerebral tissues. Assessment of changes in absorption coefficient as a function of depth in the tissue is possible when the DTOFs are analyzed using time-windows [11,15,16], statistical moments [17] or Mellin-Laplace moments [18].
Recently, it was shown that the estimation of difference in the statistical moments of DTOFs (mean time of flight and variance) obtained at two source-detector separations allows to estimate optical properties of a homogeneous medium [19]. It allows also to improve sensitivity of the time-resolved method to changes in absorption appearing deeply in the medium [20].
In the present paper we will show the results of theoretical studies on spatial distributions of factors describing sensitivity of the statistical moments of DTOFs to changes in absorption coefficient. We will apply a three-dimensional model of the medium with assumed small changes in absorption coefficient in a defined voxel of the structure. Furthermore we will apply diffusion theory in order to simulate DTOFs of photons penetrating the medium in reflectance geometry. Generation of DTOFs for a set of defined optical properties and for the medium with the absorption perturbation in the selected voxel allowed to estimate changes in statistical moments related to presence of this small absorbing inclusion. We will show how the selection and positioning of source-detector pairs located on the surface of the medium as well as combined analysis of statistical moments or their spatial derivatives leads to shaping up the sensitivity profiles of moments of DTOFs and optimize the depth-selectivity of the time-resolved diffuse reflectance.

Methods
The method of simulations which was used in the present study was previously described by Kacprzak et al. [21]. Sensitivity factors were computed using a diffusive model of light propagation in the tissue. The sensitivity factors describe relation between the calculated statistical moments of DTOFs for the homogenous medium and the ones calculated for the medium with the absorption coefficient perturbation appearing in different voxels [11,17,22].
Assuming that the source-detector separation is much larger than the mean freepath of photons the light propagation in the turbid medium can be described by the time dependent diffusion equation [23,24]. The photon fluence rate Φ(r,t) can be expressed as a function of time t, position in the medium described by vector r and the optical properties of the medium: where S(r,t) describes the photon source, c is the velocity of the light in the medium, D is the diffusion coefficient D = (3µ s ') −1 , µ a is the absorption coefficient of the medium and µ s ' is the reduced scattering coefficient [25,26]. For semi-infinite, homogenous medium and an infinitely short light pulse, the diffuse reflectance can be derived, which is number of photons reemitted from the medium per area unit per time unit. For the homogeneous medium R h (ρ,t) depends on time t and distance ρ between the source and the detector [23]: For a small change in absorption coefficient ∆µ a appearing in a small voxel of volume dV S which location inside the investigated medium is defined by vector r, the time-dependent change of diffuse reflectance ∆R(r,ρ,t) at source-detector separation ρ is: where E describes the probability per unit area that the photon emitted from location r will reach the detection spot on the surface of the medium located at the distance ρ from the source position (see definition in [27] Eq.13). Next, the distributions of times of flight of photons R diffusely reflected from the medium in which an inclusion is located at r can be calculated: Normalized statistical moments of order k of the diffuse reflectance distribution R(r,ρ,t) can be defined by: As a consequence, the changes in statistical moments of the DTOFs described by R, which are caused by a local change in absorption ∆µ a of a sampling volume dV S located at r can be calculated using following formulas [17]: where ∆A is the change in attenuation, N tot is the total number of photons, t = m 1 is the mean time of flight of photons and V = m 2 − m 2 1 is the variance of the DTOF. The quantities with index h refer to the statistical moments of the diffuse reflectance obtained for homogeneous medium R h (ρ, t).
The spatial distributions of sensitivity factors can be derived from these changes in the statistical moments obtained for each small volume dV S (voxel) indexed by i. These sensitivity factors (SFs) are defined as follows: and are called mean partial pathlength (MPP [cm]), mean time of flight sensitivity factor (MTSF [ns*cm]) and variance sensitivity factor (VSF [ns 2 *cm]). The simulated SFs can be summed up over all voxels of the medium and these sums fulfill formulas for sensitivity of the measurement of statistical moments in the homogeneous medium [11,28]: Simulations of 3D spatial distributions of the SFs were carried out for the following optical properties of the semi infinite medium: µ a =0.1cm −1 , µ s '=10 cm −1 and refractive index n = 1.4.
These optical properties are typically used in simulations related to living tissues [29]. Integrals of SFs were obtained in 3 projections by summing up the SF values along the X, Y and Z axis of the medium. The XZ surface is formed by the axes of source and detector whereas the Z axis directs towards depth of the model. The solutions of the diffusion equations which were used in this study were based on the zero boundary conditions which do not allow to account of the refractive index mismatch between the studied medium and the external medium. The extrapolated boundary conditions can be used to consider reflections of light at the boundaries [30]. However, the quantities derived in this modeling (sensitivity factors) represent ratios or differences between the quantities obtained from the homogeneous and perturbated medium. Thus, the differences between sensitivity factors obtained using solution of the diffusion equation with the zero boundary conditions and with the extrapolated boundary conditions are small.
Using the simulation procedure described above, the SFs were obtained for different sourcedetector separations and for different optical properties of the medium. Furthermore, the SFs for differences (or sums) of moments estimated at different source-detector separations can be obtained by voxel-wise subtracting (or summing up) SFs.

Results
The spatial distributions of sensitivity factors (MPP, MTSF and VSF) were first derived for source-detector separation of ρ=3cm i.e. for the source located at x = 0, y = 0,z = 0 and the detection spot located at x = 3cm, y = 0, z = 0 (geometry presented in Fig. 1). The results of simulations were presented in Fig. 2. Projections of SFs distributions on three planes (XZ, YZ and XZ) were presented. These results confirm the previous findings -the sensitivity to changes in absorption appearing in the deeper compartment of the medium increases with an increase in the order of the statistical moment [17,31]. A symmetrical character of the SFs distributions is clearly visible which is typically observed in the homogeneous medium. The negative sensitivity of the mean time of flight of photons and the variance of the DTOF is noted for changes in the absorption coefficient appearing superficially in the medium. This effect was also noted in previous simulations and experiments [21].
In Fig. 3. the distributions of SFs in XZ projection were presented for different source-detector (SD) separations. It is clear that the increase in SD separation leads to improved sensitivity of i.e. for the source located at x = 0, y = 0,z = 0 and the detection spot located at x = 3cm, y = 0, z = 0 (geometry presented in Fig. 1). the moments to changes in absorption appearing in deeper compartments of the medium. This effect was well documented in the previous theoretical and experimental studies [16,17]. Next, the sensitivities of differences between the statistical moments ∆SF acquired at two source-detector separations were analyzed. In Fig. 4, the distributions of SFs in XZ projection for two SD separations of ρ 1 =3 cm and ρ 2 =2.5 cm were presented (two detection spots located at x = 2.5 cm, y = 0 z = 0 and x = 3 cm, y = 0, z = 0). In the lowest row of the Fig. 4 the differences in the sensitivity distributions ∆SF obtained for these two SD separations were shown. The ∆SF distributions are non-symmetrical with large negative and positive sensitivity areas located close to both detection spots. Influence of the separation between these detection spots on the ∆SF distributions can be observed in Fig. 5. A smaller size of the ∆SF sensitivity area is noted when the distance between detection spots is lower. In the next step we considered potential measurement scenario in which two sources on the outside combined with two detectors in the middle were applied. Geometry of the measurement scenario was presented in Fig. 6.
The sums of the distributions of sensitivity differences ∆SF obtained for both sources Σ s1s2 ∆SF were derived and presented in Fig. 7. The position of the sources in respect to the detection spots were shown schematically in the left column of the Fig. 7. When detection spots for both considered sources are in line (middle row in Fig. 7), the sensitivity distributions were confined to small area beneath the positions of both detection spots. Interestingly, in this case in the MTSF Fig. 6. The geometry of measurement scenario in which two sources were used and for each of them two detection spots were applied. and the VSF the large positive sensitivity area is located deeper in the medium with only small components of the negative sensitivity located superficially in the medium. This case is of potential usefulness because for its technical realization only two source positions and two detection spots are necessary. Therefore, we considered influence of the distance between detection spots on the Σ s1s2 ∆SF distributions. The results of simulations obtained for ρ 1 =3 cm and three ρ 2 (2 cm, 2.5 cm and 2.75 cm) were presented in Fig. 8.
With shortening the distance between the detection spots the sensitivity area moves towards deeper compartments of the medium. Unfortunately, this effect is combined with appearance of non-homogeneities of the sensitivity in the superficially located compartments of the medium.
Influence of the reduced scattering coefficient on the distributions of sums of the sensitivities Σ s1s2 ∆SF was considered in further simulations. The results are presented in Fig. 9. It was observed that for smaller µ s ' the sensitivity area moves into deeper compartments of the medium, whereas for larger scattering a significant non-homogeneity of the sensitivity area in the superficial compartments was noted.   9. Distributions of sums of the sensitivities Σ s1s2 ∆SF obtained for two sources and two detectors (at the same location for both sources) presented in XZ projection for ρ 1 =3cm and ρ 2 = 2.5 cm, obtained for three different values of reduced scattering coefficient µ s '=20 cm −1 , 10 cm −1 and 5 cm −1 .

Discussion and conclusions
We have proposed a novel method allowing to improve depth-selective assessment of changes in absorption appearing in the optically turbid medium. The method is based on the application of time-resolved near infrared spectroscopy technique in diffuse reflectance geometry and analysis of the DTOFs acquired for two closely located detection spots and two sources located in equal distances from these detection spots (the geometry was presented in Fig. 6). It was observed that the sensitivity distributions obtained in such data acquisition scenario forms a spatially confined area located beneath the couple of detection spots. Moreover, the analysis of the spatial distributions of sensitivity of higher order moments (MTSF and VSF) shows that when the changes in absorption appear deeply in the medium the sensitivity is positive whereas when the absorption appears superficially the sensitivity is negative. Furthermore, we have studied the influence of the interoptode distance and the scattering properties of the medium on the spatial distribution of sensitivity. As expected the high sensitivity volume was shifted into deeper located compartments of the medium when the higher order statistical moments of DTOFs (especially variance) were used. It was also shown that the position of this volume of high sensitivity depends on the reduced scattering coefficient of the medium -decrease in µ s ' causes that this volume was shifted towards deeper compartments of the medium.
The results obtained in this study provided more detailed explanation of the improved sensitivity of the moments based on subtraction technique described by Milej et al. to the changes in absorption appearing in deeper compartments of the medium [19,20]. Whereas in the previous studies only single source and two detection spots were used, we propose to apply the set of two sources and two detectors to improve spatial confinement of the measurement sensitivity. Thus, the proposed method may lead to improved spatial resolution of the time-resolved NIRS measurement.
A trade-off between the detection spots distance and the signal to noise of the measurement needs to be considered [14]. It was observed that the shorter the distance the higher the sensitivity in deeper compartments of the medium can be obtained. However, shortening of the distance between detection spots leads also to decrease in amplitude of the difference in the derived statistical moments, thus, the signal to noise ratio is reduced. We have to note that the presented concept based on two detection spots and two source positions is completely symmetricalthe sources can be replaced with detectors and vice versa depending on measurement setup requirements. Moreover, the proposed method of the improvement in selectivity of the diffuse reflectance measurement works well for the mean time of flight of photons. Thus, this method can be used for data acquired using frequency-domain systems in which the phase change is directly related to the mean time of light of photons [32 -36].
As expected from the previous studies the sensitivity area of the attenuation (which reflects measurement using continuous wave NIRS technique) is mostly located in superficial compartments of the medium [17,31]. Unfortunately, the proposed technique does not improve significantly the sensitivity of attenuation to changes in absorption appearing in deeper compartments of the medium.
Concluding, the proposed technique provides improved discrimination of optical signals related to the change in absorption appearing in the deeper compartment of the optically turbid medium. The proposed method may be of use in the assessment of cerebral tissue oxygenation providing better sensitivity of the measurement to changes in oxygenation appearing in brain cortex.

Disclosures
The authors declare that there are no conflicts of interest related to this article.