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

High-density diffuse correlation tomography with enhanced depth localization and minimal surface artefacts

Open Access Open Access

Abstract

A spatially weighted filter applied to both the measurement and the Jacobian is proposed for high-density diffuse correlation tomography (DCT) to remove unwanted extracerebral interferences and artefacts along with better depth localization in the reconstructed blood flow images. High-density DCT is implemented by appropriate modification of recently introduced Multi-speckle Diffuse Correlation Spectroscopy (M-DCS) system. Additionally, we have used autocorrelation measurements at multiple delay-times in an iterative manner to improve the reconstruction results. The proposed scheme has been validated by simulations, phantom experiments and in-vivo human experiments.

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

1. Introduction

Diffuse Correlation Spectroscopy (DCS) is one of the emerging technologies to non-invasively measure blood flow based on principles of near-infrared (NIR) spectroscopy [1,2]. Diffuse intensity speckles emanating from the tissue is measured using sensitive detectors like avalanche photodiodes (APD) [3,4] or photomultiplier tubes (PMT) [2]. Temporal decorrelation of these speckles measured at multiple source-detector (SD) separations is used to quantify the blood flow. The tomographic counterpart of DCS termed as Diffuse Correlation Tomography (DCT) employs the intensity speckles obtained at multiple SD pairs along with an inversion algorithm to reconstruct three-dimensional distribution of blood flow [58].

We have recently introduced high-density DCS termed as multi-speckle Diffuse Correlation Spectroscopy (M-DCS) where a low frame rate camera was used to capture fast decorrelation of speckles due to blood flow [9,10]. In this paper, we have extended M-DCS to a high-density multi-speckle diffuse correlation tomography imaging system termed as M-DCT to reconstruct three-dimensional distribution of blood flow. One of the major advantages of using the proposed M-DCT system is that we can obtain simultaneous high-density detection of autocorrelation measurements by employing camera with low frame rate. The proposed system not only reduces the cost by several folds, but also avoids multiplexing of fiber-based detection schemes which are usually employed in DCT measurements [11].

One of the drawbacks of DCT or M-DCT is that the autocorrelation measurements are sensitive to both optical properties and blood flow in the outer layers such as scalp and skull [12,13]. This sensitivity of the measurements to extracerebral tissue properties results in poor depth localization and surface artefacts [14] in the tomographic reconstructions. Obviously, the same problem exists in Diffuse Optical Tomography (DOT) as well. One of the widely accepted and easily implementable solution to the above problem is utilizing the data from shorter SD separation for a correction either in the measurement or in the estimated flow [1519]. However, there is a lack of understanding of an optimal shorter SD separation from which these corrections can be made to minimize the influence of extracerebral blood flow changes. Additionally, the correction is not applied to the “inversion algorithm” or more specifically to the sensitivity matrices while attempting a tomography as in DCT or DOT. In order to address this problem, we propose a recipe to find the optimal shorter SD separation to design a filter for both measurement and sensitivity matrices.

Another challenge faced in DOT/DCT is that the tomographic reconstructions are always biased towards shallower regions owing to the higher sensitivity at lower SD separation, leading to poor depth localization [20,21]. To overcome this, exponentially weighted spatial regularizations are often employed [22,23]. We address this by applying the spatially dependent weight to sensitivity matrix which will account for both the biasness of data towards shorter SD separations and poor signal to noise ratio (SNR) at longer SD separations.

In most of the previously reported DCT methods, the intensity autocorrelation measurement at the characteristic decay time (τc) is always used for reconstruction [6,7], exceptions being the work reported in [24,25] where an integrated autocorrelation is used. Although this approach works reasonably well in case of simple single-layer models (or single exponent models), it fails in case of a multi-layer model. It is reported in Ref. [8] that choosing appropriate correlation delay-time for the reconstruction in tomography is essential and choosing inappropriate correlation delay-time would lead to significant errors in reconstructions. In this paper, we incorporate iterative multi delay-time scheme, wherein we iteratively reconstruct the blood flow at multiple delay-times rather than a single delay-time as done in conventional DCT systems.

To address the above-mentioned challenges associated with DCT, we propose a spatially weighted filter along with iterative multi delay-time reconstruction scheme to remove surface artefacts and achieve better depth localization in tomographic reconstruction of blood flow. We have validated the system using extensive simulation studies, experimental studies using tissue mimicking phantoms and in-vivo human studies, by measuring blood flow changes in the human brain associated with functional activation tasks.

2. Method

2.1 Brief introduction to DCT

The propagation of field autocorrelation ${G_1}({r,\tau } )$ in tissue obeys Correlation Diffusion Equation (CDE) [1,2], which is

$$\left( {\nabla \cdot ({D(r )\nabla } )- {\mu_a}(r )- \frac{1}{3}\mu_s^{\prime}k_0^2\varDelta {r^2}(\tau )} \right){G_1}({r,\tau } )={-} Q(r ). $$

Here, $\langle \varDelta {r^2}(\tau )\rangle $ is the mean-square displacement of the scattering particles which depends on both spatial variable r and delay-time $\tau $, ${k_0}\; = \;2\pi n/\lambda $ is the wave number of the continuous (CW) light diffusing through the medium, n and λ are the refractive index of media and wavelength of the incident light respectively. D, ${\mu _a}\;$ and $\mu _s^{\prime}$ are the tissue optical properties such as diffusion coefficient, absorption and scattering coefficient respectively and Q denotes the source. Assuming a Brownian motion model for blood flow, mean square displacement is given by $\langle \varDelta {r^2}(\tau )\rangle = 6{D_B}\tau $, where ${D_B}$ is the particle diffusion coefficient which is proportional to blood flow.

Using the CDE as forward model, the inverse model for DCT can be defined using the first-Born approximation as

$$g_1^\delta ({{r_s},{r_d},\tau } )\equiv \frac{{G_1^\delta ({{r_s},{r_d},\tau } )}}{{G_1^0({{r_s},{r_d},0} )}} ={-} \frac{{2u_s^{\prime}k_0^2}}{{G_1^0({{r_s},{r_d},0} )}}\smallint G({{r_d},r,\tau } )G({r,{r_s},\tau } )D_B^\delta (r )\tau dr, $$
where ${g_1}({{r_s},{r_d},\tau } )$ is the normalized field autocorrelation and G is the Green’s function solution to Eq. (1). The variables ${r_d}$ and ${r_s}$ denote the spatial coordinates of the detector and source respectively. The above equation can be discretized to form a system of linear equation as $Jx = b$, where J is the sensitivity/ Jacobian matrix, b is the change in normalized field autocorrelation measurement $g_1^\delta $ and x is the change in blood flow $D_B^\delta $. Owing to the illposedness, the above system is usually solved using Tikhonov regularized least-square minimization.

2.2 Spatially weighted filter design

The autocorrelation measurement at the boundary is sensitive to both the deep tissue blood flow and the flow from extracerebral layers like scalp [2628]. In order to address this, we design a filter by studying the perturbation in the spectrum of measurement due to perturbation in the flow at individual layer of the object representing scalp, skull and brain. The designed filter is then applied to both measurement and the Jacobian to filter out the surface artefacts and interferences in the reconstructed images. We also address the biased sensitivity of tomography to shallower regions by rescaling the Jacobians with a spatially weighted skewed Gaussian function.

We define the Fourier transform of the perturbation in measurement as

$$\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )\, = \;{\widetilde {{g_1}}^{\;\delta ,sb}} - \;{\widetilde {{g_1}}^{\delta ,\;s}}, $$
where ${\widetilde {{g_1}}^{\delta ,sb}}$ corresponds to perturbation in both scalp and brain while ${\widetilde {{g_1}}^{\delta ,s}}$ corresponds to perturbation only in scalp. In order to look at the dependence of $\varDelta \widetilde {{g_1}}$ on SD separation, we define the sensitivity S as
$$S = \frac{\partial }{{\partial r}}[{{}_\omega^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |} ].\;$$

The nature of variation of S with respect to SD separation will help us to determine the parameters of the filter as explained below.$S$ is a monotonically decreasing function of SD for a single layer. However, for multiple layers, S shows isolated abrupt changes corresponding to different layers. We define $S{D_{opt}}$ to be the SD separation at which S shows the first abrupt change as it increases. Clearly $S{D_{opt}}$ lies near to the shorter SD separations which corresponds to perturbations in scalp.

To design the filter, we consider the Fourier transform of the function ${e^{ - k\surd \tau }}$ [29] (denoted as $Q(\omega ))$, which is similar to the Green’s function solution to CDE in terms of functional dependence on $\tau $. We perform a non-linear least-square fitting of ${\widetilde {{g_1}}^\delta }({S{D_{opt}},\omega } )$ to $Q(\omega )$ which is again denoted as Q itself. The filter is then defined as

$$\tilde F(\omega )= 1 - \frac{{Q(\omega )}}{{\max ({Q(\omega )} )}} + \varepsilon , $$
where $\varepsilon $ is a constant shift. $\tilde F(\omega )$ essentially acts as a high pass filter. To compute $Q(\omega )$, we make the following assumption based on finite weighted sum of exponentials
$$Q(\omega )\equiv FT\left( {{e^{ - k\sqrt \tau }}} \right) \cong FT\left( {\mathop \sum \limits_{i = 1}^N {a_i}{e^{{b_i}\tau }}} \right), $$
where FT denotes the Fourier transform. The rationale of this approximation is that FT of simple exponentials can be easily computed (analytically or numerically).

We apply the filter to both measurement and the Jacobian by modifying the Born approximation given in Eq. (2). The filtered measurement ${\widetilde {{g_1}}^{\;\delta ,F}}({{r_s},{r_d},\omega } )$ is given by

$${\widetilde {{g_1}}^{\;\delta ,F}}({{r_s},{r_d},\omega } )= \;{\widetilde {{g_1}}^{\;\delta }}({{r_s},{r_d},\omega } )\times Q(\omega ). $$

Similarly, we define the filtered Green’s function

$${\tilde G^{\;F}}({{r_1},{r_2},\omega } )= \;{\tilde G}({{r_1},{r_2},\omega } )\times Q(\omega )$$
for any two spatial points ${r_1}$ and ${r_2}$.

In order to address the biased sensitivity of the reconstruction towards shallower region, we introduce a weight function W defined as

$$W(z )= 2CDF\left( {\frac{1}{{\sqrt {2\pi {\sigma^2}} }}{e^{ - \frac{{\alpha {{({z - \mu } )}^2}}}{{2{\sigma^2}}}}}} \right)\frac{1}{{\sqrt {2\pi {\sigma ^2}} }}{e^{ - \frac{{{{({z - \mu } )}^2}}}{{2{\sigma ^2}}}}}, $$
where CDF denotes the cumulative distribution function and z denotes the z co-ordinate representing the depth of the sample. The rationale of choosing a skewed Gaussian function as weight is that a reduced weight is needed at lower SD separation to account for biased sensitivity. At the same time, for longer SD separations the SNR of the measurement will be very low and hence again necessitates reduced weight to avoid noise triggered artefacts in the reconstruction. The depth of the brain is determined using conventional DCS as explained in section 2.5.3. We plot the SNR of the measurement data against SD separations to find the SD upto which SNR is greater than dark noise of the camera. Based on these data, we fit the parameters mean $\mu $, variance ${\sigma ^2}$ and skewness $\alpha $ to the function $W(z )$. With the weighted filter being applied to Eq. (2), the modified Born approximation is given as
$$g_1^{\delta ,F}({{r_s},{r_d},\tau } )={-} \frac{{2u_s^{\prime}k_0^2}}{{G_1^0({{r_s},{r_d},0} )}}\smallint \tau {G^F}({{r_d},r,\tau } ){G^F}({r,{r_s},\tau } )W(r )D_B^\delta dr, $$
where $g_1^{\delta ,F}$ is the inverse FT of ${\widetilde {{g_1}}^{\;\delta ,F}}$ and similarly for other functions.

2.3 Multi delay-time iterative scheme

The selection of optimal delay-time for reconstruction of blood flow is very important as improper choice of τ can lead to misinterpreted results. In many literatures, it is common to choose optimal delay-time as characteristic decay time, ${\tau _c}$. Though it works well for single layer model, it fails to hold true for multi-layer models. The discretization of Eq. (10) results in ${J_\tau }x = {b_\tau }\;$ which is solved by minimizing the functional $\left\|J_{\tau_i} x_i-b_{\tau_i}\right\|{ }_2^2+\lambda\left\|x_i-x_{i-1}\right\|$ for every τ in an iterative manner. By this method, a prior knowledge on ${\tau _c}$ is not required as it is required for single τ based tomographic reconstruction. As the measurement and the Jacobians are already filtered, multiple τ scheme can be safely employed without concerns about the bias towards lower SD separation and noises in the measurements.

2.4 High-density multi speckle diffuse correlation tomography (M-DCT)

We have recently introduced a high-density DCS system called M-DCS capable of employing a low frame rate camera to capture the intensity autocorrelation and thereby quantify the blood flow. We extend the M-DCS to M-DCT using the filtered Born approximation as given in Eq. (10) in section 2.2. The details of the M-DCS theory and implementation for in-vivo human studies can be found in Ref. [9,10,30].

2.5 Validation methods

2.5.1 Simulation studies

We have implemented both forward and inverse algorithm using Finite Element Method (FEM). The forward data is simulated by solving the CDE given in Eq. (1) implemented in MATLAB platform. The simulation studies are carried out with a Windows 10, 8 GB RAM 10th generation Intel system. The mesh size used for simulation is 4cm×4cm×3cm with 52111 nodes and 288000 tetrahedral elements. The proposed method takes around 55 minutes for solving the CDE in the high-density mesh and reconstruction of the blood flow. The background absorption coefficient ${\mu _a}$, reduced scattering coefficient $\mu _s^{\prime}$ and particle diffusion coefficient ${D_B}$ used are 0.19 cm-1, 6.6cm-1 and 1 × 10−8 cm2/s respectively for scalp, 0.136 cm-1, 8.6 cm-1 and 0 cm2/s respectively for skull and 0.186 cm-1, 11.10 cm-1 and 6 × 10−8 cm2/s respectively for brain [27,31,32]. A spherical inhomogeneity in ${D_B}$ having diameter 0.3cm is incorporated in the scalp where it is perturbed by 20% from baseline and another spherical perturbation of 0.8cm diameter is given in brain where ${D_B}$ is increased by 50% as shown in Fig. 1(a). The source and detector positions are depicted as in Fig. 1(b). The reconstructions are done using Tikhonov least square minimization method.

 figure: Fig. 1.

Fig. 1. Mesh to simulate 3-layer head model: (a) scalp, skull and brain with two inhomogeneities and (b) positions of sources and detectors indicated by red stars and green dots respectively

Download Full Size | PDF

2.5.2 Multi-layer phantom experiments

The experimental setup for the study is shown in Fig. 2. A 785nm collimated laser source with beam shaping optics is used to form a point source of diameter 0.9 mm. A Galvo mirror is used to scan the surface of the phantom. A tissue mimicking phantom is prepared using Polydimethylsiloxane, Intralipid and glycerol [3335] to mimic a 3-layer human brain with scalp, skull and the brain having thickness of 0.3cm, 0.4cm and 2cm respectively. We perform the phantom experimental study and reconstruct the flow in 3-layer phantom where perturbation in flow is introduced using syringe pump to the entire top layer and inside a pipe of diameter 0.4cm kept at 0.9cm depth from the surface in the third layer. The ${\mu _a}$ and $\mu _s^{\prime}$ of the top layer of the phantom is estimated to be 0.2 cm-1 and 6 cm-1, the middle layer is estimated to be 0.12 cm-1 and 8.1 cm-1 and the third layer is estimated to be 0.19 cm-1 and 10.5cm-1. The light source is focused on the sample and the scattered intensity is measured using a sCMOS camera (Photometrics Prime BSI) with an objective lens (Navitar, Zoom 7000), in reflection geometry. The Field of View (FOV) captured is 4 cm ×4 cm with 900 × 900 pixels in the FOV. We employ M-DCS to measure the normalized field autocorrelation ${g_1}$ from multi exposure speckle contrast data at several sources. The filter is then applied to ${g_1}$ as explained in section 2.2.

The total number of source detector pairs used are 13689. A total of 10 exposures between 100us to 10ms and 100 images for each exposure are acquired which involves an acquisition time of 25s per source. The spatially weighted filter is applied to Jacobians corresponding to every SD pairs. The regularization parameter λ for reconstruction is chosen to be 0.9 by trial-and-error.

 figure: Fig. 2.

Fig. 2. Experimental setup for flow reconstruction using two different flow rates at different layers in a 3-layer phantom.

Download Full Size | PDF

2.5.3 In-vivo human experiments

Healthy human subjects (n = 7) including both male and female aged between 20 and 30 years have voluntarily participated for the study. The tomographic scans were carried out in the subjects in a supine position. A laser source (785 nm, 3.8mW, 0.95mm diameter) was focused in the forehead to access the prefrontal cortex at FP1 area. A Galvo mirror was used to scan the laser at nine source points. The camera (Photometrics Prime BSI) was focused to a FOV of 4cm×4cm and was used to capture 100 images at multi exposure time ranging from 100us to 10ms to measure ${g_1}$ using M-DCS. The total number of nodes and elements in the rectangular FEM mesh are 52111 and 288000 forming a mesh of size 4cm×4cm×3cm. For visualization purpose, a mesh based on human-atlas [36] is used for FEM solution of the forward and inverse problem to reconstruct $D_B^\delta $ for one of the subjects. In the prefrontal region of human atlas, a region of 4cm×4cm×3cm consisting of 194967 nodes and 262368 tetrahedral elements is used for the tomographic reconstruction.

Protocol: We have performed a controlled experimental protocol to perturb flow in extracerebral layer and deep tissue independently. While a functional activation task perturbs blood flow in the prefrontal cortex, a heating pad will help to control the blood flow in the scalp. The blood flow in the brain was perturbed by using a number processing task, wherein the subject was given a random number between 5000 and 6000 and was asked to subtract seven from the random number serially. The prefrontal cortex becomes activated during this cognitive task [9,37].

Seven subjects with a total of 21 readings were taken for the analysis. Three sets of studies were performed for each trial: (a) a serial 7 subtraction task which perturbs the blood flow in prefrontal cortex of the brain (b) serial 7 subtraction with a heating pad placed on the subject’s forehead to increase the scalp blood flow along with perturbation in prefrontal cortex and (c) the heating pad alone is kept to perturb only the surface blood flow and the subject is asked not to perform any task which is used for validation. The temperature at the surface was recorded and it was maintained at a constant surface temperature of 38.3℃. The order of the tasks was randomized for 3 different trials of each subject. A 4 minutes baseline data was taken before each task along with a 1-minute break post task. This is done to make sure that the subject comes back to the initial state after each task and the brain activity during the three tasks remain independent of each other. The timeline of the protocol along with experimental setup and SD positions are shown in Fig. 3.

In order to show the perturbation of surface blood flow by using the heating pad, we have also performed Laser Speckle Contrast Imaging (LSCI) to quantify the surface blood flow using a uniform laser source of 785nm. A speckle contrast data was obtained for 200 images at an exposure of 1ms. The heating pad was kept in the proximity of the FP1 region which does not fall under the FOV of the camera.

 figure: Fig. 3.

Fig. 3. Experimental setup for tomographic imaging of blood flow in prefrontal cortex: (a) protocol timeline, (b) M-DCT set-up for the experiment and (c) sources and detectors indicated by red and green dots respectively.

Download Full Size | PDF

For estimating the depth of the brain for different subjects, we measure ${g_1}$ by DCS system at various SD separations using APD, hardware autocorrelator and optical fibers. Hereafter, the above system is referred as conventional DCS system to differentiate it from the M-DCS system. The quantity ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$, the sensitivity S and the spatial variance of S for these SD separations are computed. The spatial variance of S calculated over a sliding window of 3 SDs, shows a maxima where the measured data has maximum perturbation in measurement from the baseline. We show that this SD separation with maximum perturbation can estimate the approximate depth of the brain. According to the photon diffusion theory, the depth at which NIR light penetrates biological tissues is roughly equal to half the SD separation distance [38]. The depth of the brain estimated by the above method is consistent with the depth of brain reported in Ref. [2628]. As mentioned in section 2.2, we design the proposed filter based on the depth estimated using the above method.

2.6 Validation parameters

The key parameters for the simulation, phantom and human studies ($\mu _s^{\prime}$, ${\mu _a}$, layer thickness, DB) for scalp, skull and brain are provided in the Table 1. All the parameters for simulation studies and $\mu _s^{\prime}$, ${\mu _a}$ and scalp thickness for human studies are taken from Ref. [27]. The parameters used for the human studies are provided for the subject for whom the results are displayed in result section and these properties vary among subjects.

In order to find the effect of the filter along with iterative multi delay-time scheme on the reconstruction of blood flow in extracerebral layer and deep tissue, the parameters relative scalp blood flow (rSBF) and relative cerebral blood blow (rCBF) have been used for simulation and human studies, where $rCBF = {({D_B^\delta /D_B^0} )_{brain}}$ and $rSBF = {({D_B^\delta /D_B^0} )_{scalp}}$.

Tables Icon

Table 1. Optical properties and thickness of scalp, skull and brain used for simulation, phantom and in-vivo human studies

To quantitatively compare the reconstruction obtained by the proposed method, we use the evaluation metrics such as mean square error (MSE) and structural similarity index (SSIM) [39,40] formulated as

$${\rm{MSE = }}{\left[ {\left( {\sum\limits_{i = 1}^N {|{{x_{true}}(i) - {x_{recon}}(i)} |} } \right)/N} \right]^{\rm{2}}}{\rm{ + }}\left( {\sum\nolimits_{i = 1}^N {|{{x_{true}}(i) - {{\overline x }_{recon}}} |} } \right)/N, $$
$${\rm{SSIM}} = \frac{{({2{{\bar x}_{true}}{{\bar x}_{recon}} + {c_1}} )({2{\sigma_{true,}}_{recon} + {c_2}} )}}{{({\bar x_{true}^2 + \bar x_{recon}^2 + {c_1}} )({\sigma_{true}^2 + \sigma_{recon}^2 + {c_2}} )}}. $$

Here, ${x_{true}}$ is the known original diffusion coefficient and ${x_{recon}}$ is the corresponding reconstructed diffusion coefficient. The constants c1 and c2 are the measurements of stability expressed as ${c_1} = {({0.01 \times L} )^2}$ and ${c_2} = {({0.03 \times L} )^2}$, where L is the dynamic range of the pixel values. ${\sigma _{true}}$ and ${\sigma _{recon}}$ are the local standard deviation of the original and the reconstructed diffusion coefficient respectively, while ${\sigma _{true,recon}}$ denotes the covariance between original and the reconstructed diffusion coefficients.

3. Results

3.1 Optimal parameter selection for filter design of high-density DCT system

The filter is designed based on frequency spectrum of field autocorrelation measurement from human head. In human, the extracerebral blood flow is prominent in the first layer that is scalp whereas the blood flow in skull counts to almost zero. We plot ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$ for different SDs as shown in Fig. 4(a) and the corresponding sensitivity S as shown in Fig. 4(b). We select 0.7cm as $S{D_{opt}}$ from our simulation because the slope of S becomes flat from 0.75cm in Fig. 4(b) which means the signal from 0.8cm starts getting affected by deep tissue perturbation whereas the signals up to SD = 0.7cm contain information mostly from the scalp.

 figure: Fig. 4.

Fig. 4. Choice of $S{D_{opt}}:$ (a) ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$ is plotted for each SD and (b) the sensitivity S plotted as function of mean SD.

Download Full Size | PDF

In addition, an appropriate constant shift $\varepsilon $ is decided such that the DC component of the filter gets non zero low weightage and the higher frequency components get more weightage. The rationale behind this selection is that the deeper layer perturbation affects dominantly the higher frequency components than lower frequency components. An error analysis (MSE) between the reconstruction of $D_B^\delta $ and background $D_B^\delta $ without extracerebral perturbation as shown in Fig. 5(b) has been performed for different values of $\varepsilon $ ranging from 1.3 to 2.1 and $\varepsilon $ = 1.5 gives minimum error in reconstruction. It should also be noted here that $\varepsilon $ ranging from 1.5 to 1.9 gives relatively less error. In this paper, we have used $\varepsilon $ = 1.5 in the design of filter using the proposed method. The Fig. 5(a) shows the filter $\tilde F(\omega )$ design as given in Eq. (5) using $\varepsilon $ = 1.5.

We consider three z values to find corresponding spatial weight function W(z). The first point is at z = 0, which is the surface at which the data is measured. The second point is at the depth of brain which is estimated at approximately z = 0.8 cm in our case. The third point denotes the point beyond which noise dominates the measurement and for Photometrics Prime BSI camera we found it at 2.2 cm. We fit the points 0 cm, 0.8cm and 2.2cm according to the maxima and minima to a skewed Gaussian curve as shown in Fig. 6.

 figure: Fig. 5.

Fig. 5. Filter design and choice of $\varepsilon $: (a) $\tilde F(\omega )$ with $S{D_{opt}} = 0.7cm$ and $\varepsilon = 1.5$ and (b) MSE for different values of $\varepsilon $ ranging from 1.3 to 2.1.

Download Full Size | PDF

 figure: Fig. 6.

Fig. 6. Spatial weight function W(z) designed based on depth of brain and SNR of the measurement data.

Download Full Size | PDF

After we apply the filter on the measurement data and Green’s function solution to CDE, the curves shift towards the left as shown in Fig. 7 which means they have moved towards the earlier delay-time values which has more information from the deep tissue than the extracerebral layers.

 figure: Fig. 7.

Fig. 7. Effect of filter: (a) plot of filtered and unfiltered $g_1^\delta $ and (b) plot of filtered and unfiltered $G$. The filtered curves shift towards earlier part of τ, which corresponds mostly to the deep tissue.

Download Full Size | PDF

Simulation study:

The reconstruction of $D_B^\delta $ for the three-layer model with two different spherical perturbations in ${D_B}$, one denoting the perturbation in the extracerebral layer and the other in the brain are shown for the following cases: (a) original data (b) single delay-time (τ = 2.5 × 10−5s, which is approximately τc at 2cm) (c) single delay-time (τ = 2.5 × 10−5s) with spatial weight W(z) (d) single delay-time (chosen at different τc at different SD separations) with W(z), (e) single delay-time with filter applied only on measurement data along with W(z), (f) multiple delay-times with W(z), (g) multiple delay-times with the proposed filter without W(z) and (h) multiple delay-times with the proposed filter and W(z). The above cases are chosen to show the study and analyse the effects of multiple delay times, spatial weight W(z) and the proposed filter.

 figure: Fig. 8.

Fig. 8. Reconstruction of $D_B^\delta $ in XY plane for different depths (Z axis) of 3-layer model a$\to $ original data, b$\to $ single delay-time (${\tau _c}$ at single SD), c$\to $ single delay-time (${\tau _c}$ at single SD) and W(z), d$\to $ single delay-time (different ${\tau _c}$ at different SDs) and W(z), e$\to $ single delay-time with filter on measurement only and W(z), f$\to $ multi delay-time and W(z), g$\to $ multi delay-time and filter without W(z), h$\to $ multi delay-time and filter with W(z).

Download Full Size | PDF

From the results shown in Fig. 8, we can see that the reconstruction using a single delay-time (row 2) are not well localised in depth and both the perturbations are prominent in all the layers. When spatial weight function W(z) is applied (row 3), both the perturbations are reduced in the top layer but they get equal weightage in the deep layer. While there is an improvement in the reconstruction, the effect of extracerebral perturbation can still be seen when different τc for different SDs (row 4) are used for reconstruction. In order to remove the extracerebral perturbation, we use the proposed filter on both measurement and Jacobian along with spatial weight function W(z) and iterative multi delay-time scheme. If the filter along with W(z) is applied only on the measurement data (row 5), the extracerebral perturbation reduces significantly but is not removed completely. The iterative multi delay-time scheme (row 6) improves the reconstruction results by several folds when compared to case c and d. When the proposed filter is applied to both the measurement data and Jacobian along with the iterative multi delay-time scheme (row 7), the extracerebral perturbation and artefacts are removed. The deep layer perturbation gets better spatial localisation when the filter, the iterative multi delay-time scheme and W(z) are applied together (row 8). The Fig. 9(a) and (b) show the rSBF and rCBF of the reconstruction of the extracerebral perturbation located at (1.5,1.5,2.85) cm for different cases at a depth of 0.2cm and 1cm respectively. The MSE between the reconstruction of $D_B^\delta $ and background $D_B^\delta $ without extracerebral perturbation for different cases (case b-h) is shown in Fig. 9(c).

 figure: Fig. 9.

Fig. 9. Simulation results: (a) rSBF comparison for reconstruction of the extracerebral perturbation at a depth of 0.2cm (b) rCBF comparison for reconstruction of the extracerebral perturbation at a depth of 1cm and (c) MSE for reconstruction of $D_B^\delta $

Download Full Size | PDF

We use the iterative multi delay-time scheme for reconstruction instead of choosing a single delay-time. The number of iterations used for the multiple delay-time scheme is selected as 7 in our case where delay-times are chosen at equal interval logarithmically between 1us and 10ms. The SSIM and MSE between the reconstruction of $D_B^\delta $ and background $D_B^\delta $ without extracerebral perturbation do not vary much after the 7th iteration as shown in Fig. 10. The plots are shown for iterative multi delay-time reconstruction with and without filter.

 figure: Fig. 10.

Fig. 10. The SSIM and MSE plots show that seven iterations over delay-time is optimal beyond which there is negligible change in the reconstruction.

Download Full Size | PDF

In addition, we have compared the reconstruction quantitatively in terms of MSE to show that using a high-density of detectors improves the reconstruction. For this, we have used two spherical inhomogeneities of 0.3cm diameter with 1cm separation between their centers at a depth of 1cm from surface and reconstructed the inhomogeneity using multiple sets of detectors. As the distance between the detectors increases, the number of detectors decreases. The Fig. 11 shows that as the distance between detectors increases, the MSE also increases, indicating that the quality of reconstruction can be improved by using a high-density of detectors.

 figure: Fig. 11.

Fig. 11. The MSE of the reconstruction increases as the distance between detector increases explaining the need for a high-density measurement for tomography.

Download Full Size | PDF

Multi-layer tissue phantom experiment:

The reconstruction results for 3-layer tissue mimicking flow phantom are shown in Fig. 12, where the entire top layer is perturbed and a spatial perturbation inside a pipe of diameter 0.4cm is given in the third layer at a depth of 1cm from the surface. The estimated fitted ${D_B}$ for first and third layer are 0.92 × 10−8 cm2/s and 1.87 × 10−8 cm2/s. The estimated fitted ${D_B}$ for the top layer perturbation is 1.21 × 10−8 cm2/s and for the perturbation in the pipe is 2.94 × 10−8 cm2/s. The ${D_B}$ values of the tissue mimicking phantoms are measured for individual layers by using conventional DCS system and the analytical solution for the semi-infinite solution to CDE is used for fitting the ${g_1}$ measurement obtained from the DCS system. As evident from the results, the effect of iterative multi delay-times based reconstruction along with filter and spatial weight W(z) gives reconstruction of deep layer perturbation with better spatial localization and less interferences from surface layers.

 figure: Fig. 12.

Fig. 12. The reconstructed $D_B^\delta $ with and without filter for a 3-layer flow phantom having two different flow rates at different depths.

Download Full Size | PDF

3.2 Human experiment

We perform different experiments for perturbing the flow in the extracerebral layer and deep tissue either simultaneously or one at a time as discussed in section 2-5-3. The autocorrelation measurements obtained for different cases from the M-DCT system at a SD separation of 2 cm are shown in Fig. 13(a). The mean and standard deviation of three trials of a single subject for difference in autocorrelation measurements for cognitive number task, heating pad and both of them together from the baseline autocorrelation measurement are shown in Fig. 13(b). The rightmost curve represents the baseline flow and leftmost curve represents the number task along with heating pad which perturbs the blood flow in both brain and at the surface. The curves for the perturbation of the scalp and the brain individually lie in between. It can be seen that the effect of perturbing the blood flow at scalp by the heating pad can be seen in the later part of the correlation delay-time more than the earlier part. Also, perturbing the blood flow in the prefrontal cortex by the number task affects the earlier part of the correlation delay-time more than the later part. The Fig. 13(c) shows the effect of filter on autocorrelation measurements obtained by M-DCT and conventional DCT for number task at a SD of 1.5cm.

 figure: Fig. 13.

Fig. 13. Measurement obtained from M-DCT for different perturbations: (a) autocorrelation measurements, (b) mean and standard deviation of three trials of a single subject for difference in autocorrelation measurements from baseline and (c) effect of filter on autocorrelation measurement by M-DCT and conventional DCT for number task

Download Full Size | PDF

To estimate the depth of the brain, the mean and standard deviation of the quantity ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$ and the sensitivity S for SDs ranging from 0.5cm to 2cm for five subjects (3 trials each) are plotted in Fig. 14(a) and (b) respectively. The spatial variance of S shows a peak at mean SD of 1.65cm as depicted in Fig. 14(c). This gives us an approximate estimation of brain depth ${\cong} $ 0.8cm.

 figure: Fig. 14.

Fig. 14. Estimation of depth of brain from autocorrelation measurement for five subjects: (a) the quantity ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$ (b) the sensitivity S and (c) the spatial variance of S gives a maxima at 1.65cm estimating the brain depth to be ${\cong} $0.8cm.

Download Full Size | PDF

The reconstruction of $D_B^\delta $ for the cognitive number task is shown in Fig. 15(a) under different conditions. The artefacts and extracerebral interferences are reduced and the spatial localization of the perturbation in the prefrontal cortex of the brain improves due to the filter with spatial weight along with iterative multi delay-time scheme. The Fig. 15(b) and (c) show the reconstruction of $D_B^\delta $ in a human-atlas [36], rather than a rectangular geometry, using conventional single delay-time and multiple delay-time with filter and W(z) respectively.

 figure: Fig. 15.

Fig. 15. Reconstruction of $D_B^\delta $ during the cognitive number task: (a) at different depths for different cases, (b) human-atlas representation using single delay-time and (c) multiple delay-time with filter and spatial weight W(z)

Download Full Size | PDF

The reconstruction of $D_B^\delta $ between number task along with heating pad and baseline is shown in Fig. 16(a) for single delay-time and multiple delay-time with filter having spatial weight W(z). The results show that even if there is an external perturbation in the scalp due to the heating pad, the filter is able to remove it and the final result becomes comparable to the result in Fig. 15(a) where there is no external perturbation.

 figure: Fig. 16.

Fig. 16. Reconstruction of $D_B^\delta $ at different depths during the cognitive number task along with heating pad: (a) depth wise reconstruction (b) mean and standard deviation of speckle contrast from LSCI of blood flow with and without heating pad and difference between them measured at FP1 region of brain for a single subject.

Download Full Size | PDF

In order to show the perturbation of surface blood flow by using the heating pad, we have shown the mean and standard deviation of speckle contrast κ2 obtained from LSCI on a single subject in Fig. 16(b). The speckle contrast κ2 decreases when the heating pad is used due to increase in surface blood flow.

Finally, we show the mean and standard deviation of rSBF and rCBF for 3 trials of intra-subject (n = 3) and inter-subject (n = 21) for number task and number task along with heating pad in Fig. 17(a) and (b) respectively. The baseline DB is estimated to be 1.1 × 10−8 cm2/s for scalp, 0 for skull and 2.8 × 10−8 cm2/s for brain. For the estimation of DB in scalp and brain, we have fitted ${g_1}$ from 10−6s to 10−2s at a small SD = 0.6cm and ${g_1}$ from 10−6s to 10−4s at a long SD = 2.5cm to analytical solution of CDE for semi-infinite geometry respectively. The DB for skull is fixed at 0 from Ref. [26,27]. The rSBF is around 1% when the filter is used for reconstruction, indicating the efficiency of the filter in removing the extracerebral blood flow. The mean rCBF difference is around 5% when the filter is not used for the cases of number task with and without heating pad indicating an increase in error in estimating the rCBF. The use of the filter reduces the mean rCBF difference to only 2% between the cases of number task with and without heating pad.

 figure: Fig. 17.

Fig. 17. Mean and standard deviation of rSBF and rCBF indicating the effect of filter at a depth of 0.2cm and 1cm respectively for 3 trials of single human subject and 2l trials of all human subjects for: (a) number task and (b) number task along with heating pad

Download Full Size | PDF

In order to find the statistical significance of the blood flow due to the number task and the effect of filter on the reconstruction of blood flow at extracerebral layer, Student’s t-test was performed for 3 trials of 7 subjects each. A two-tailed paired t-test (n = 21) for rCBF before and after the task at a depth of 1cm from the surface was performed. It reveals that the cognitive number task produces a significant increase in rCBF with t = 7.84, p < 0.000001. For illustrating the effect of filter, another t-test (n = 21) was performed for rSBF with and without filter at a depth of 0.2cm from surface for before and after the number task. The t values are found out as t = 2.075, p > 0.050 and t = 8.05, p < 0.000001 respectively with and without filter. This shows that there is no significant difference in the rSBF before and after the number task when the filter is used.

4. Discussion

The measurements for DCS or DCT in the human brain are often corrupted by instrumental artefacts and interferences from extracerebral layers [15,4143], resulting in incorrect rCBF in the reconstruction results. These are usually corrected by using appropriate multi-layer modelling [26,43] or partial volume fractions [44] to remove or reduce the influence of extracerebral layers. In this paper, to minimize the interferences due to extracerebral layers, we introduce a filter designed using parameters determined from the measurements at shorter SD separation and subsequently applied to both measurement and the Jacobian. In reflection geometry, the reconstructions are biased towards the shallower region due to the inherent higher sensitivity of Jacobians at shorter SD separation, which is corrected by using a spatial weight along with the above-mentioned filter in the Jacobian matrix. While the filter removes the artefacts in reconstruction from the extracerebral layers, the spatial weight incorporated in it takes care of the inherent biasness of the inversion algorithm towards shallower regions.

It is well known that inappropriate selection of delay-time may lead to a large error in reconstruction. In case of multi-layer tissue, the selection of a single delay-time for tomographic reconstruction may not contain the entire information of the decorrelation measurement. In addition, the autocorrelation at a single delay-time may be corrupted by noise [8,25]. In this paper, we use a iterative multi delay-time algorithm that uses the entire autocorrelation measurement, rather than a single delay-time, resulting in better reconstruction. The results in Figs. 8, 12 and 15 indicate improvement in reconstructions when iterative multi delay-time scheme is employed in contrast to the single delay-time as commonly used in DCT systems.

In Ref. [24], an integral of the autocorrelation over multiple delay-times to give a single value is used for the reconstruction of blood flow rather than using multiple delay-times. The authors have also reconstructed a non-linear form of mean square displacement for each delay-time in order to extract the storage and loss moduli. In Ref. [45], a Nth order linear DCT has been proposed which also uses an integral form of ${g_1}(\tau )$ with a Nth order Taylor polynomial which does not seek the solution for CDE, rather utilizes the autocorrelation function at fifty delay-times to find the blood flow based on linear regression. The researchers show that the blood flow obtained by the iterative multi delay-time scheme have better accuracy and robustness when compared to a single delay-time based blood flow. In order to use the entire autocorrelation measurement, the above-mentioned research works have employed the integral of the full autocorrelation measurement for the tomography, wherein the fine details of the autocorrelation measurement are integrated out, which are crucial in a multi-layer media. In our paper, we use the entire autocorrelation at multiple delay-times, rather than using an integrated quantity. The autocorrelation delay-times used in this paper is in the range of 1us to 10ms distributed equally on a logarithmic scale. The use of a greater number of delay-times (currently this paper employs seven delay-times) does not make any significant change in the reconstruction results as shown in Fig. 10. However, the range of delay-times used for reconstruction depends on the type of tissue and may vary for tissues having different optical properties. We would like to note that in time-domain DCS systems, by selecting appropriate multi-exponential model, time of arrival of photons and multi correlation delay-time, it is possible to separate the extracerebral blood flow from deep tissue blood flow [33]. However, such approaches are not possible in continuous wave DCS which necessities for the proposed filter.

The framework proposed for the filter with spatial weight for three-dimensional reconstruction of blood flow has been extensively validated using simulation studies as shown in Fig. 8. For the experimental validation, we use a multi-speckle DCT system based on our recently introduced M-DCS system [9]. M-DCS uses a low frame rate camera for measuring intensity autocorrelation detection rather than using expensive high sampling rate detectors with multi-channel hardware auto correlators, which reduces the cost by several folds leading to a high-density DCS system. In this paper, we have extended M-DCS to tomography (M-DCT) resulting in a high-density DCT system with relatively inexpensive instrumentation. For simultaneous speckle detection using array detectors for blood-flow measurements, any other previously reported methods [4648] can also be used. The filter design and the associated algorithm presented here are equally applicable to all such high-density speckle-based blood flow imaging methods.

Decoupling the effect of perturbations in the scalp from reconstructions using data from shorter SD separations has been reported in diffuse optical spectroscopy/tomography where spatial filtering methods like principal component analysis (PCA) [4951] or general linear model (GLM) [15] or its combination (PCA-GLM) [12] were employed. However, the corrections are applied only to the intensity measurements or to the derived oxygen saturation. Our proposed method combines both filter and the spatial weight to a single function which is applied to both measurement and the Jacobian matrix to reconstruct the blood flow. The simulation studies shown in section 2.5.1 shows the effect of reconstruction for several cases of filter and spatial weight along with single or multiple delay-times. As evident from the results in Figs. 8, 12 and 15, the reconstruction using single delay-time with spatial weight introduces surface artefacts and interference from the scalp. The multiple delay-time reconstruction scheme with spatial weight and without filter gives better spatial localization but with surface artefacts contributed by flow in the scalp. Similarly, the reconstruction using multiple delay-time scheme and filter but without any spatial weight results in minimal surface artefacts but poorly localized flow. Therefore, we need a combination of filter along with spatial weight using iterative multi delay-time scheme for reconstruction in DCT to remove all the extracerebral interferences and artefacts and get a good spatial localization for the perturbation present in the brain.

One of the demerits of our method is the need for prior information on the depth of the brain and skull from the surface. Although we have come up with a new method to approximately determine the depth of the brain, it may be more appropriate to extract the information about the thickness of the layers from CT/MRI [52]. However, an error analysis is also performed on the data based on the forward model with various depths ranging between 0.6 cm to 1 cm (-25% to 25% of the actual depth), while the inverse model assumes the depth of brain at 0.8 cm, the rCBF would give a maximum of 4% error with filter and a maximum of 9% error without filter. Additionally, we would like to note that a generic human head atlas is used to mimic the effects of curvature of human head. However, the mesh needs to be further customized according to the individual subject head size for better reconstruction.

The current in-vivo experimental studies involve the acquisition of M-DCT data based on 100 images of 10 exposure times, at 900 × 900 pixels, which takes around 25s per source. This could be further reduced by optimizing the source-detector pairs [53], acquisition of data with lesser resolution and by utilizing synthetic exposure datas [54]. By adapting our recently proposed algorithm [40], which avoids the need of inverting the Jacobian matrix, the total computational time of M-DCT system can be considerably reduced thus moving towards real-time tomography.

Funding

Wadhwani Research Centre for Bioengineering (WRCB), Indian Institute of Technology Bombay; Department of Biotechnology, Ministry of Science and Technology, India (Ramalingaswamy Fellowship-2016); Department of Science and Technology, Ministry of Science and Technology, India (SERB–Early career research award); Indian Institute of Technology Bombay (SEED GRANT).

Acknowledgments

Portions of this work were presented at the Biophotonics Congress: Biomedical Optics 2022 conference titled, ‘A Spatially Weighted Filter Design for High Density Diffuse Correlation Tomography’(JM3A.22) [55].

Disclosures

The authors have no conflict of interest.

Data availability

The data that support the findings of this study are available upon reasonable request from the authors.

References

1. T. Durduran and A. G. Yodh, “Diffuse correlation spectroscopy for non-invasive, micro-vascular cerebral blood flow measurement,” NeuroImage 85, 51–63 (2014). [CrossRef]  

2. D. A. Boas, L. E. Campbell, and A. G. Yodh, “Scattering and Imaging with Diffusing Temporal Field Correlations,” Phys. Rev. Lett. 75(9), 1855–1858 (1995). [CrossRef]  

3. C. Huang, Y. Lin, L. He, D. Irwin, M. M. Szabunio, and G. Yu, “Alignment of sources and detectors on breast surface for noncontact diffuse correlation tomography of breast tumors,” Appl. Opt. 54(29), 8808–8816 (2015). [CrossRef]  

4. S. Han, A.R. Proctor, J. B. Vella, D. S. Benoit, and R. Choe, “Non-invasive diffuse correlation tomography reveals spatial and temporal blood flow differences in murine bone grafting approaches,” Biomed. Opt. Express 7(9), 3262–3279 (2016). [CrossRef]  

5. C. Zhou, G. Yu, D. Furuya, J. H. Greenberg, A. G. Yodh, and T. Durduran, “Diffuse optical correlation tomography of cerebral blood flow during cortical spreading depression in rat brain,” Opt. Express 14(3), 1125–1144 (2006). [CrossRef]  

6. J. P. Culver, T. Durduran, D. Furuya, C. Cheung, J. H. Greenberg, and A. G. Yodh, “Diffuse optical tomography of cerebral blood flow, oxygenation, and metabolism in rat during focal ischemia,” J. Cereb. Blood Flow Metab. 23(8), 911–924 (2003). [CrossRef]  

7. S. Han, J. Johansson, M. Mireles, A. R. Proctor, M. D. Hoffman, J. B. Vella, D. S. W. Benoit, T. Durduran, and R. Choe, “Non-contact scanning diffuse correlation tomography system for three-dimensional blood flow imaging in a murine bone graft model,” Biomed. Opt. Express 6(7), 2695 (2015). [CrossRef]  

8. Y. Lin, C. Huang, D. Irwin, L. He, Y. Shang, and G. Yu, “Three-dimensional flow contrast imaging of deep tissue using noncontact diffuse correlation tomography,” Appl. Phys. Lett. 104(12), 121103 (2014). [CrossRef]  

9. H. M. Varma and K. Murali, “Multi-speckle diffuse correlation spectroscopy to measure cerebral blood flow,” Biomed. Opt. Express 11(11), 6699–6709 (2020). [CrossRef]  

10. K. Murali, A. K. Nandakumaran, T. Durduran, and H. M. Varma, “Recovery of the diffuse correlation spectroscopy data-type from speckle contrast measurements: towards low-cost, deep-tissue blood flow measurements,” Biomed. Opt. Express 10(10), 5395–5413 (2019). [CrossRef]  

11. J. Ren, S. Han, A. R. Proctor, D. E. Desa, G. A. Ramirez, V. R. D. Ching-Roa, J. B. Majeski, I. A. Dar, N. E. Barber, A. M. Forti, D. S. W. Benoit, and R. Choe, “Longitudinal 3D Blood Flow Distribution Provided by Diffuse Correlation Tomography during Bone Healing in a Murine Fracture Model,” Photochem. Photobiol. 96(2), 380–387 (2020). [CrossRef]  

12. F. Zhang, D. Cheong, A. Khan, Y. Chen, L. Ding, and H. Yuan, “Correcting physiological noise in whole-head functional near-infrared spectroscopy,” J. Neurosci. Methods 360, 109262 (2021). [CrossRef]  

13. W. Zhou, M. Zhao, O. Kholiqov, and V. J. Srinivasan, “Multi-exposure interferometric diffusing wave spectroscopy,” Opt. Lett. 46(18), 4498–4501 (2021). [CrossRef]  

14. T. Yamada, S. Umeyama, and M. Ohashi, “Removal of motion artifacts originating from optode fluctuations during functional near-infrared spectroscopy measurements,” Biomed. Opt. Express 6(12), 4632–4649 (2015). [CrossRef]  

15. L. Gagnon, K. Perdue, D. N. Greve, D. Goldenholz, G. Kaskhedikar, and D. A. Boas, “Improved recovery of the hemodynamic response in diffuse optical imaging using short optode separations and state-space modeling,” NeuroImage 56(3), 1362–1371 (2011). [CrossRef]  

16. N. M. Gregg, B. R. White, B. W. Zeff, A. J. Berger, and J. P. Culver, “Brain specificity of diffuse optical imaging: improvements from superficial signal regression and tomography,” Front. Neuroenerg. 2, 14 (2010). [CrossRef]  

17. J. Virtanen, T. E. Noponen, and P. Meriläinen, “Comparison of principal and independent component analysis in removing extracerebral interference from near-infrared spectroscopy signals,” J. Biomed. Opt. 14(5), 054032 (2009). [CrossRef]  

18. L. Gagnon, M. Yücel, D. A. Boas, and R. J. Cooper, “Further improvement in reducing superficial contamination in NIRS using double short separation measurements,” NeuroImage 85, 127–135 (2014). [CrossRef]  

19. H. Zhao and E. M. Buckley, “Influence of source–detector separation on diffuse correlation spectroscopy measurements of cerebral blood flow with a multilayered analytical model,” Neurophotonics 9(03), 035002 (2022). [CrossRef]  

20. S. R. Arridge and J. C. Schotland, “Optical tomography: forward and inverse problems,” Inverse Problems 25(12), 123010 (2009). [CrossRef]  

21. B. Pogue, T. McBride, J. Prewitt, U. L. Österberg, and K. D. Paulsen, “Spatially variant regularization improves diffuse optical tomography,” Appl. Opt. 38(13), 2950–2961 (1999). [CrossRef]  

22. H. Niu, F. Tian, Z. J. Lin, and H. Liu, “Development of a compensation algorithm for accurate depth localization in diffuse optical tomography,” Opt. Lett. 35(3), 429–431 (2010). [CrossRef]  

23. B. White, “Developing High-Density Diffuse Optical Tomography for Neuroimaging,” (2012). https://www.proquest.com/docview/926209544?pq-origsite=gscholar&fromopenview=true

24. H. M. Varma, A. Nandakumaran, and R. M. Vasu, “Study of turbid media with light: Recovery of mechanical and optical properties from boundary measurement of intensity autocorrelation of light,” J. Opt. Soc. Am. A 26(6), 1472–1483 (2009). [CrossRef]  

25. J. Zuo, X. Zhang, J. Lu, Z. Gui, and Y. Shang, “Impact of reconstruction algorithms on diffuse correlation tomography blood flow imaging,” IEEE Access 8, 31882–31891 (2020). [CrossRef]  

26. H. Zhao, E. Sathialingam, and E. M. Buckley, “Accuracy of diffuse correlation spectroscopy measurements of cerebral blood flow when using a three-layer analytical model,” Biomed. Opt. Express 12(11), 7149–7161 (2021). [CrossRef]  

27. J. Li, L. Qiu, C. S. Poon, and U. Sunar, “Analytical models for time-domain diffuse correlation spectroscopy for multi-layer and heterogeneous turbid media,” Biomed. Opt. Express 8(12), 5518–5532 (2017). [CrossRef]  

28. O. Kholiqov, W. Zhou, T. Zhang, M. Zhao, S. Ghandiparsi, and V. J. Srinivasan, “Scanning interferometric near-infrared spectroscopy,” Opt. Lett. 47(1), 110–113 (2022). [CrossRef]  

29. K. Murali and H. M. Varma, “Laser speckle simulation using stochastic differential equation for blood flow imaging applications,” in Microscopy Histopathology and Analytics (Optica Publishing Group, 2022), p. JM3A. 27.

30. K. Murali, A. K. Nandakumaran, and H. M. Varma, “On the equivalence of speckle contrast-based and diffuse correlation spectroscopy methods in measuring in vivo blood flow,” Opt. Lett. 45(14), 3993–3996 (2020). [CrossRef]  

31. Y. Hoshi, M. Shimada, C. Sato, and Y. Iguchi, “Reevaluation of near-infrared light propagation in the adult human head: implications for functional near-infrared spectroscopy,” J. Biomed. Opt. 10(6), 064032 (2005). [CrossRef]  

32. K. Verdecchia, M. Diop, A. Lee, L. B. Morrison, T. Y. Lee, and K. S. Lawrence, “Assessment of a multi-layered diffuse correlation spectroscopy method for monitoring cerebral blood flow in adults,” Biomed. Opt. Express 7(9), 3659–3674 (2016). [CrossRef]  

33. S. Samaei, P. Sawosz, M. Kacprzak, Ż. Pastuszak, D. Borycki, and A. Liebert, “Time-domain diffuse correlation spectroscopy (TD-DCS) for noninvasive, depth-dependent blood flow quantification in human tissue in vivo,” Sci. Rep. 11(1), 1817 (2021). [CrossRef]  

34. R. Choe, “Diffuse optical tomography and spectroscopy of breast cancer and fetal brain,”, (2005), https://www.proquest.com/docview/305412841?pq-origsite=gscholar&fromopenview=true.

35. A. Kienle, F. Foschum, and R. Michels, “Optical properties of fat emulsions,” Opt. Express 16(8), 5907–5925 (2008). [CrossRef]  

36. “NEVA EM | Electromagnetic Modeling and Simulation,” https://www.nevaelectromagnetics.com/high-resolution-head-models.

37. C. Hock, F. Müller-Spahn, S. Schuh-Hofer, M. Hofmann, U. Dirnagl, and A. Villringer, “Age dependency of changes in cerebral hemoglobin oxygenation during brain activation: A near-infrared spectroscopy study,” J. Cereb. Blood Flow Metab. 15(6), 1103–1108 (1995). [CrossRef]  

38. D. Irwin, L. Dong, Y. Shang, R. Cheng, M. Kudrimoti, S. D. Stevens, and G. Yu, “Influences of tissue absorption and scattering on diffuse correlation spectroscopy blood flow measurements,” Biomed. Opt. Express 2(7), 1969–1985 (2011). [CrossRef]  

39. J. Feng, Q. Sun, Z. Li, Z. Sun, and K. Jia, “Back-propagation neural network-based reconstruction algorithm for diffuse optical tomography,” J. Biomed. Opt. 24(05), 1 (2018). [CrossRef]  

40. R. Paul, K. Murali, S. Chetia, and H. M. Varma, “A simple algorithm for diffuse optical tomography without Jacobian inversion,” Biomed. Phys. Eng. Express 8(4), 045001 (2022). [CrossRef]  

41. D. Milej, A. Abdalmalak, A. Rajaram, and K. St. Lawrence, “Direct assessment of extracerebral signal contamination on optical measurements of cerebral blood flow, oxygenation, and metabolism,” Neurophotonics 7(04), 045002 (2020). [CrossRef]  

42. M. M. Wu, S.-T. Chan, D. Mazumder, D. Tamborini, K. A. Stephens, B. Deng, P. Farzam, J. Y. Chu, M. A. Franceschini, J. Z. Qu, and S. A. Carp, “Improved accuracy of cerebral blood flow quantification in the presence of systemic physiology cross-talk using multi-layer Monte Carlo modeling,” Neurophotonics 8(01), 015001 (2021). [CrossRef]  

43. J. Li, G. Dietsche, D. Iftime, S. E. Skipetrov, G. Maret, T. Elbert, B. Rockstroh, and T. Gisler, “Noninvasive detection of functional brain activity with near-infrared diffusing-wave spectroscopy,” J. Biomed. Opt. 10(4), 044002 (2005). [CrossRef]  

44. A. G. Yodh, C. Zhou, G. Yu, J. Wang, J. H. Greenberg, J. A. Detre, M. G. Burnett, and T. Durduran, “Diffuse optical measurement of blood flow, blood oxygenation, and metabolism in a human brain during sensorimotor cortex activation,” Opt. Lett. 29(15), 1766–1768 (2004). [CrossRef]  

45. X. Zhang, Z. Gui, Z. Qiao, Y. Liu, and Y. Shang, “Nth-order linear algorithm for diffuse correlation tomography,” Biomed. Opt. Express 9(5), 2365–2382 (2018). [CrossRef]  

46. J. Dong, K. Lee, and R. Bi, “Deep tissue flowmetry based on diffuse speckle contrast analysis,” Opt. Lett. 38(9), 1401–1403 (2013). [CrossRef]  

47. H. M. Varma, C. P. Valdes, A. K. Kristoffersen, J. P. Culver, T. Durduran, A. Devor, S. Sakadži, V. Srinivasan, M. Yaseen, K. Nizar, P. Saisan, P. Tian, A. Dale, S. Vinogradov, M.A. Franceschini, M. J. Leahy, J. G. Enfield, N. T. Clancy, J. ODoherty, P. McNamara, and G. E. Nilsson, “Speckle contrast optical tomography: A new method for deep tissue three-dimensional tomography of blood flow,” Biomed. Opt. Express 5(4), 1275–1289 (2014). [CrossRef]  

48. C. Huang, D. Irwin, Y. Lin, Y. Shang, L. He, W. Kong, J. Luo, and G. Yu, “Speckle contrast diffuse correlation tomography of complex turbid medium flow,” Med. Phys. 42(7), 4000–4006 (2015). [CrossRef]  

49. T. J. Huppert, S. G. Diamond, M. A. Franceschini, and D. A. Boas, “HomER: a review of time-series analysis methods for near-infrared spectroscopy of the brain,” Appl. Opt. 48(10), D280–D298 (2009). [CrossRef]  

50. T. Sato, I. Nambu, K. Takeda, T. Aihara, O. Yamashita, Y. Isogaya, Y. Inoue, Y. Otaka, Y. Wada, M. Kawato, M. aki Sato, and R. Osu, “Reduction of global interference of scalp-hemodynamics in functional near-infrared spectroscopy using short distance probes,” NeuroImage 141, 120–132 (2016). [CrossRef]  

51. Y. Zhang, D. H. Brooks, M. A. Franceschini, and D. A. Boas, “Eigenvector-based spatial filtering for reduction of physiological interference in diffuse optical imaging,” J. Biomed. Opt. 10(1), 011014 (2005). [CrossRef]  

52. M.M. Wu, K. Perdue, S. T. Chan, K. A. Stephens, B. Deng, M. A. Frasnceschini, and S.A. Carp, “Complete head cerebral sensitivity mapping for diffuse correlation spectroscopy using subject-specific magnetic resonance imaging models,” Biomed. Opt. Express 13(3), 1131–1151 (2022). [CrossRef]  

53. M. E. Eames, B. W. Pogue, P. K. Yalavarthy, H. Dehghani, S. Srinivasan, and B. W. Pogue, “An efficient Jacobian reduction method for diffuse optical image reconstruction,” Opt. Express 15(24), 15908–15919 (2007). [CrossRef]  

54. A. M. Safi, C. Hernandez-Isidro, S. Cini, S. Moka, M. Harrah, C. L. Passaglia, and A. B. Parthasarathy, “Quantitative cerebral blood flow imaging with synthetic single-shot multi-exposure laser speckle imaging,” Biophotonics Congress2021, p. BW3B.4.

55. R. Paul, K. Murali, and H. M. Varma, “A spatially weighted filter design for high density diffuse correlation tomography,” in Optics and the Brain (Optica Publishing Group, 2022), p. JM3A.22.

Data availability

The data that support the findings of this study are available upon reasonable request from the authors.

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 (17)

Fig. 1.
Fig. 1. Mesh to simulate 3-layer head model: (a) scalp, skull and brain with two inhomogeneities and (b) positions of sources and detectors indicated by red stars and green dots respectively
Fig. 2.
Fig. 2. Experimental setup for flow reconstruction using two different flow rates at different layers in a 3-layer phantom.
Fig. 3.
Fig. 3. Experimental setup for tomographic imaging of blood flow in prefrontal cortex: (a) protocol timeline, (b) M-DCT set-up for the experiment and (c) sources and detectors indicated by red and green dots respectively.
Fig. 4.
Fig. 4. Choice of $S{D_{opt}}:$ (a) ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$ is plotted for each SD and (b) the sensitivity S plotted as function of mean SD.
Fig. 5.
Fig. 5. Filter design and choice of $\varepsilon $: (a) $\tilde F(\omega )$ with $S{D_{opt}} = 0.7cm$ and $\varepsilon = 1.5$ and (b) MSE for different values of $\varepsilon $ ranging from 1.3 to 2.1.
Fig. 6.
Fig. 6. Spatial weight function W(z) designed based on depth of brain and SNR of the measurement data.
Fig. 7.
Fig. 7. Effect of filter: (a) plot of filtered and unfiltered $g_1^\delta $ and (b) plot of filtered and unfiltered $G$. The filtered curves shift towards earlier part of τ, which corresponds mostly to the deep tissue.
Fig. 8.
Fig. 8. Reconstruction of $D_B^\delta $ in XY plane for different depths (Z axis) of 3-layer model a$\to $ original data, b$\to $ single delay-time (${\tau _c}$ at single SD), c$\to $ single delay-time (${\tau _c}$ at single SD) and W(z), d$\to $ single delay-time (different ${\tau _c}$ at different SDs) and W(z), e$\to $ single delay-time with filter on measurement only and W(z), f$\to $ multi delay-time and W(z), g$\to $ multi delay-time and filter without W(z), h$\to $ multi delay-time and filter with W(z).
Fig. 9.
Fig. 9. Simulation results: (a) rSBF comparison for reconstruction of the extracerebral perturbation at a depth of 0.2cm (b) rCBF comparison for reconstruction of the extracerebral perturbation at a depth of 1cm and (c) MSE for reconstruction of $D_B^\delta $
Fig. 10.
Fig. 10. The SSIM and MSE plots show that seven iterations over delay-time is optimal beyond which there is negligible change in the reconstruction.
Fig. 11.
Fig. 11. The MSE of the reconstruction increases as the distance between detector increases explaining the need for a high-density measurement for tomography.
Fig. 12.
Fig. 12. The reconstructed $D_B^\delta $ with and without filter for a 3-layer flow phantom having two different flow rates at different depths.
Fig. 13.
Fig. 13. Measurement obtained from M-DCT for different perturbations: (a) autocorrelation measurements, (b) mean and standard deviation of three trials of a single subject for difference in autocorrelation measurements from baseline and (c) effect of filter on autocorrelation measurement by M-DCT and conventional DCT for number task
Fig. 14.
Fig. 14. Estimation of depth of brain from autocorrelation measurement for five subjects: (a) the quantity ${}_\omega ^{max}|{\varDelta \widetilde {{g_1}}({{r_s},{r_d},\omega } )} |$ (b) the sensitivity S and (c) the spatial variance of S gives a maxima at 1.65cm estimating the brain depth to be ${\cong} $0.8cm.
Fig. 15.
Fig. 15. Reconstruction of $D_B^\delta $ during the cognitive number task: (a) at different depths for different cases, (b) human-atlas representation using single delay-time and (c) multiple delay-time with filter and spatial weight W(z)
Fig. 16.
Fig. 16. Reconstruction of $D_B^\delta $ at different depths during the cognitive number task along with heating pad: (a) depth wise reconstruction (b) mean and standard deviation of speckle contrast from LSCI of blood flow with and without heating pad and difference between them measured at FP1 region of brain for a single subject.
Fig. 17.
Fig. 17. Mean and standard deviation of rSBF and rCBF indicating the effect of filter at a depth of 0.2cm and 1cm respectively for 3 trials of single human subject and 2l trials of all human subjects for: (a) number task and (b) number task along with heating pad

Tables (1)

Tables Icon

Table 1. Optical properties and thickness of scalp, skull and brain used for simulation, phantom and in-vivo human studies

Equations (12)

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

( ( D ( r ) ) μ a ( r ) 1 3 μ s k 0 2 Δ r 2 ( τ ) ) G 1 ( r , τ ) = Q ( r ) .
g 1 δ ( r s , r d , τ ) G 1 δ ( r s , r d , τ ) G 1 0 ( r s , r d , 0 ) = 2 u s k 0 2 G 1 0 ( r s , r d , 0 ) G ( r d , r , τ ) G ( r , r s , τ ) D B δ ( r ) τ d r ,
Δ g 1 ~ ( r s , r d , ω ) = g 1 ~ δ , s b g 1 ~ δ , s ,
S = r [ ω m a x | Δ g 1 ~ ( r s , r d , ω ) | ] .
F ~ ( ω ) = 1 Q ( ω ) max ( Q ( ω ) ) + ε ,
Q ( ω ) F T ( e k τ ) F T ( i = 1 N a i e b i τ ) ,
g 1 ~ δ , F ( r s , r d , ω ) = g 1 ~ δ ( r s , r d , ω ) × Q ( ω ) .
G ~ F ( r 1 , r 2 , ω ) = G ~ ( r 1 , r 2 , ω ) × Q ( ω )
W ( z ) = 2 C D F ( 1 2 π σ 2 e α ( z μ ) 2 2 σ 2 ) 1 2 π σ 2 e ( z μ ) 2 2 σ 2 ,
g 1 δ , F ( r s , r d , τ ) = 2 u s k 0 2 G 1 0 ( r s , r d , 0 ) τ G F ( r d , r , τ ) G F ( r , r s , τ ) W ( r ) D B δ d r ,
M S E = [ ( i = 1 N | x t r u e ( i ) x r e c o n ( i ) | ) / N ] 2 + ( i = 1 N | x t r u e ( i ) x ¯ r e c o n | ) / N ,
S S I M = ( 2 x ¯ t r u e x ¯ r e c o n + c 1 ) ( 2 σ t r u e , r e c o n + c 2 ) ( x ¯ t r u e 2 + x ¯ r e c o n 2 + c 1 ) ( σ t r u e 2 + σ r e c o n 2 + c 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.