Abstract
Laser speckle-based blood flow imaging is a well-accepted and widely used method for pre-clinical and clinical applications. Although it was introduced as a method to measure only superficial blood flow (< 1mm depth), several recently introduced variants resulted in measuring deep tissue blood flow (a few cm) as well. A means of simulating laser speckles is often necessary for the analysis and development of these imaging modalities, as evident from many such attempts towards developing simulation tools in the past. Such methods often employ Fourier transforms or statistical tools to simulate speckles with desired statistical properties. We present the first method to use a stochastic differential equation to generate laser speckles with a pre-determined probability density function and a temporal auto-correlation. The method allows the choice of apriori gamma distribution along with simple exponential or more complex temporal auto-correlation statistics for simulated speckles, making it suitable for different blood flow profiles. In contrast to the existing methods that often generate speckles associated with superficial flow, we simulate both superficial and diffuse speckles leading to applications in deep tissue blood flow imaging. In addition, we have also incorporated appropriate models for noise associated with the detectors to simulate realistic speckles. We have validated our model by comparing the simulated speckles with those obtained from in-vivo studies in mice and healthy human subject.
© 2022 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement
1. Introduction
Non-invasive imaging of microvascular blood flow using coherent laser source illumination and associated scattered speckle intensity [1,2] measurements have been studied over several decades with promising results. Conventionally, as in laser speckle contrast imaging (LSCI) [3] systems, the sample is illuminated uniformly and the speckle contrast of the measured speckles is used to quantify the surface blood flow up to a depth of 1mm. A spatio-temporal speckle contrast measured either at single exposure time [3] or multiple exposure time [4] of the camera are often employed for this purpose. With the advent of high-speed cameras and single-photon avalanche diode (SPAD) arrays, utilizing intensity autocorrelation instead of the time-averaged speckle contrast to quantify blood flow is gaining much attention [5–8].
On the other hand, diffuse correlation spectroscopy (DCS) uses a focused laser illumination to measure the deep tissue blood flow (order of a few cm) using spatially resolved diffuse speckle measurements. Speckle intensity measured at multiple source-detector (SD) separations are used to quantify blood flow, with the propagation model being the Correlation diffusion equation (CDE) [2,9]. Most often, a high-density measurement of speckles at several SD pairs is required either for better signal-to-noise ratio (SNR) [9] or to implement so-called called diffuse correlation tomography (DCT) [10]. However, in diffuse speckle contrast analysis (DSCA) [11], an inexpensive alternative to high-density autocorrelation measurements, speckle contrast from multi-speckle measurements from low frame rate cameras was used to measure deep tissue blood flow. Concurrently, speckle contrast optical spectroscopy (SCOS) [12] and its tomographic counterpart called speckle contrast optical tomography (SCOT) [13] were introduced with a better quantification of flow using spatially resolved speckle contrast measurements.
There have been other variants of the speckle based blood flow imaging such as interferometric DWS (i-DWS) [14], speckle contrast diffuse correlation tomography (sc-DCT) [15], multi-speckle DCS (M-DCS) [16,17], multi-exposure interferometric diffusing wave spectroscopy (mi-DWS) [18], Interferometric speckle visibility spectroscopy (ISVS) [19], Integrated Diffuse Speckle Contrast Spectroscopy (iDSCS) [20], holographic Fourier domain DCS [21], and inverse of the numerical integration of squared g1 (INISg1) [22]. In addition to the blood flow measurements, laser speckles have been employed to measure viscoelastic properties as in laser speckle rheology [23] and ultrasound-assisted DWS [24,25].
Although many such modalities use laser speckles to measure blood flow, a rigorous methodology for the simulation of laser speckles is lacking. It is not only essential for improving the quantification of blood flow by optimizing the source and detector characteristics but also for having a quantitative comparison among these modalities. Some of the existing methods used to generate laser speckles in the past are reported as follows. One common approach is using a fast Fourier transform of the phase matrix, which generates independent speckle patterns, but without any pre-determined correlation statistics [26]. In Ref [27], Duncan et al., proposed a method based on a copula algorithm for generating speckle sequences for a given auto-correlation model that utilizes the concepts of the quantile function and Fourier transform. However, this method was not utilized for generating spatially varying diffuse speckle patterns. The concept of coherent imaging principles was used to generate spatially varying speckle images using pre-defined correlation coefficient matrices in Ref [28]. A comprehensive review of the above-mentioned methods used for generating speckles can be found in chapter 3 of Ref [29]. Very recently, in Ref [30], James et al., have used a method based on the Karhunen–Loève expansion of the electric field to simulate time-integrated speckle patterns. All the above-mentioned work stress the need for laser speckle simulation tools for blood flow imaging applications.
A new approach for the simulation of speckle intensity based on stochastic differential equation (SDE) is introduced in this paper. A SDE, in simpler terms, is an ordinary differential equation with a stochastic source term, whose solution is a stochastic process. SDE has several applications in various fields that include, but not limited to, chemistry, finance, epidemiology, mechanics and microelectronics [31,32]. Here, for the first time, we apply the SDE to simulate laser speckles for blood flow imaging. The laser speckles originating from a turbid medium are characterized by probability density function (PDF) which is either negative exponential or gamma and with a simple exponential or complex temporal autocorrelation [33]. The coefficients of the SDE are determined in such a way that the above-mentioned statistics of the speckles can be incorporated into the simulation of speckles. In addition, spatially varying diffuse speckles are also simulated by appropriately modifying the coefficients of SDE using the solution to photon diffusion equation with a potential for application in deep-tissue blood flow imaging. In order to make this into a comprehensive laser speckle simulation tool, we have also applied realistic noise models to the simulated speckle pattern based on noise associated with the common speckle measuring detectors like array detectors (cameras) and avalanche photodiodes (APD).
We have compared the simulated speckles using SDE with the $in-vivo$ experimental speckles generated from the cerebral blood flow in mice and healthy adult humans. The presented speckle simulation tool incorporating both conventional and diffuse speckles with the flexibility to incorporate pre-determined PDF and complex temporal autocorrelations could guide researchers in choosing appropriate models and techniques to achieve better blood flow quantification.
2. Simulation model for speckle intensity using SDE
2.1 Theory
Consider a SDE which consists of an ordinary differential equation (ODE) along with a stochastic term [34–36], given as
Here a and b are called the drift and diffusion terms, respectively, and W is the Wiener process (or Brownian motion) such that $dW(t) \approx \sqrt {\delta t} N(0,1)$, where N denotes a Gaussian random variable with zero mean and unit variance. It is well known that functions $a$ and $b$ satisfy the Fokker-Planck equation [31], which is given as
Here $p$, PDF of random variable $x$, is assumed to be independent of time, as it is a stationary random process. Note that $a$ and $b$ implicitly depend on time. For $x(t)$ to represent speckle intensity, it should have a negative exponential PDF, i.e., $p(x)=\frac {1}{\mu } exp(-x/\mu )$, where the mean and the standard deviation are equal to $\mu$ [33]. Additionally, the auto-correlation of the speckle intensity obeys exponential decay [1,37]. The usual procedure to enforce this in the solution of Eq. (1) is to choose the parameters $a$ and $b$ using Eq. (2) appropriately. It can be seen that the terms $a(x(t)) = -\alpha (x-\mu )$ will ensure a correlation decay of the form $e^{-\alpha \tau }$ , which results in the equation $dx(t) = -\alpha (x - \mu )dt + b(x(t))dW(t)$ [34,35]. This predetermined $a$ is now substituted into Eq. (2) along with the exponential PDF to deduce $b$ as follows:
2.2 Computational implementation
Although there are several numerical methods to solve Eq. (3), we chose to use the well known Milstein scheme [32,34,35,38,39]. The Eq. (3) solved by Milstien scheme is given by
The SDE algorithm used to generate speckles, $I(t)$, is shown as a pseudocode below in Algorithm 1.
$*-$Input parameters; $\mu$ is mean intensity based on photon diffusion equation and $\tilde {\alpha }$ is from the auto-correlation function, based on the flow.
$\#-$randn - MATLAB command - generate random numbers
$*\#-$gives speckles with mean $\mu$ and correlation decay $1/\tilde {\alpha }$.
The autocorrelation $g_2$ was calculated using the standard MATLAB function $xcorr$ (with an unbiased estimate of the cross-correlation), wherein $g_2$ is defined as $g_2(\tau ) =\frac {(<I(t)I(t+\tau )>)}{((N-\mid m \mid )<I(t) >^2)}$, $N$ is the number of samples, and $m$ is the index of the lag $\tau$. This unbiased estimate of $xcorr$ is essential for obtaining a proper estimate of intensity autocorrelation, while biased and other forms of autocorrelation results in distorted or undesired autocorrelation for longer time integration. The effect of using biased/ unscaled autocorrelation is shown in Fig S3 of the supplemental document, wherein the biased auto-correlation does not follow the desired auto-correlation at a longer delay time, while unbiased auto-correlation follows the desired auto-correlation. Here $g_2^N=\frac {(<I(t)I(t+\tau )>)}{max(<I(t)I(t+\tau )>)}$ and $g_2^N (Biased/unscaled) =\frac {(g_2-1)}{\beta }$. This is due to the insufficient number of samples leading to an inherent bias at larger values of $\tau$, resulting in distortion of autocorrelation.
2.3 Simulation parameters
Depending on the application, we can take auto-correlation to be independent of spatial coordinates $r$ (in a sample with uniform illumination as in LSCI/ MESI) or dependent on $r$ as in diffuse speckles (DCS, DSCA, SCOS, sc-DCT, M-DCS and i-DWS). We have used $t_{min}=10^{-6} s,t_{max}=1s$ and $N=10^7$ samples to generate the speckle intensity.
2.3.1 Surface blood flow
For surface blood flow, we have normalized intensity auto-correlation $g_2(\tau )=1+\beta |g_1^2|$, where $g_1=e^{-t/\tau _c}$ and $\tau _c$ is the characteristic decay time (based on Siegert’s relation [2,40], $\beta$ is a constant depending on the collection optics of the experiment). In order to demonstrate the capability of the proposed method in simulating complicated autocorrelation models (henceforth referred to as complex autocorrelations) other than the simple exponential based one, we consider the autocorrelation model proposed in Ref [6] given by ,
Here $\beta$ is the instrumentation constant, $X$ and $d$ dictate the dynamic scattering regime, namely, multiple unordered (for $X=0.5$), single unordered or multiple ordered regime (for $X=1$) and single ordered regime (for $X=2$) and $\rho$ is the fraction of dynamic scatterers. Three representative tissue types were used for the simulations: large vessels, small vessels and parenchymas. The values of the above parameters for modelling were chosen based on Ref [6]. For the large vessels, $\tau _c=50 \mu s$ was used along with $\rho =0.98, d=0$ and $X=1$. For small vessels, we use $\tau _c=250 \mu s$, $\rho =0.98, d=1$ and $X=2$ and for the parenchyma $\tau _c=250 \mu s, \rho =0.8$, $d=1$ and $X=0.5$ was used.
2.3.2 Simulation of complex auto-correlation models
Any complex exponential function (such as $e^{k\sqrt {\tau }}$ or $e^{k \tau ^2}$ or its combinations) can be approximated as a finite weighted sum of exponentials [41], as given by $g_2=\sum _{i=0}^N (1/\theta _i) e^{a_i \tau }$. We solve Eq. (3) to obtain $I(t)$, for every $0 \leq i \leq n$, such that with mean $\mu =\tilde {\mu }/N$ and auto-correlation of the form $e^{a_i\tau }$. The resulting $I(t)'s$ were scaled by $1/\theta _i$ and added to get the final intensity of speckle pattern. Here we perform a non-linear least-squares fit to determine $a_i^{'}s$ for a fixed set of $\theta _i^{'}s$, where $a_i$ is the co-efficient of the fit, $\theta _i$ is the scaled factor and $n$ is the number of terms ($n=4$ is selected and $\theta _i$ were predetermined for these simulations). The terms $\tilde {\mu }$ and $\mu$ denote the mean intensity used in individual simulations (for every term $'n'$) and desired mean intensity. It is to be noted here that the individual simulations gives a random variables of the negative exponential distribution, while the resultant summation of random variable gives rise to gamma distribution [33,42] along with the desired complex autocorrelation.
While generating speckle intensity with the sum of such exponentials, the SDE in Eq. (3) can be modelled as, $dI(t)= -\alpha (I(t)- \alpha.p)dt +\sqrt {(2\alpha \mu I(t))} dW(t)$, which generates the speckles $I(t)$, with a PDF of the form $\frac {(1/\mu )^p I^{p-1} e^{-I/\mu }}{(J(p))}$; $J(p)=(p-1)\!$ and auto-correlation of the form $e^{-\alpha t}$. Note that the above PDF is a gamma distribution and a special case with $p=1$ leads to an exponential distribution. Here $p=1,2{\ldots }n,$ with $p$ being a positive integer. The Feller’s condition (i.e. $\frac {2\alpha \mu }{\sigma ^2} \geq 1$) for the above equation can be written as $\frac {2\alpha \mu p}{(\sqrt {2\alpha \mu })^2} = p$. As $p$ is a positive integer and $p\geq 1$, Feller’s condition [31] is always satisfied, thus generating positive speckle intensity for complex models.
2.3.3 Deep tissue blood flow
For deep tissue blood flow, we use a semi-infinite medium solution of the correlation diffusion equation (CDE) [2,37]. The above mentioned method for the generation of complex autocorrelation is used for the generation of the speckles. For the simulation, $\mu _a=0.1 cm^{-1}$, $\mu _s^{'}=12 cm^{-1}$ and $D_B=10^{-8} cm^{2}/s$ was used. In many cases, it may be sufficient to use a binomial approximation to the complex autocorrelation instead of adopting the approximation using weighted sum of simple exponential functions. The complex autocorrelation based on the infinite solution of CDE may be approximated (based on Binomial Approximation) to $g_1=e^{-\alpha _D \tau }$, where $\alpha _{D} =\frac {6 \mu _{s}^{'^2} k_{0}^2 D_{B}r}{2\sqrt {3\mu _{a}\mu _{s}^{'} }}$. Here, $\mu _a,\mu _s^{'},k_0$ and $D_B$ denote the absorption coefficient, reduced scattering coefficient, wave number and particle diffusion coefficient respectively. Although the results based on binomial approximation is not shown in this paper, it can be employed when the error in approximation is within the required accuracy.
2.3.4 In vivo animal experiment
To validate the SDE based algorithm, we have performed an in-vivo animal experiment in mice (Balp/C, 36 g, M, 4 months). This work was undertaken with the approval of the Institute Animal Ethical Committee (IAEC) at the APT research centre, Pune, India. The animal was anaesthetized using Ketamine and Xylazine and the skull of the mice was removed. For the uniform illumination, the laser beam (785 nm) was illuminated via a diffuser and for obtaining a point source (diameter $\leq 0.3 mm$), a focusing lens was used, as shown in Fig S1(a-b) in supplemental document. The detector used was a CCD camera (Basler acA 640-120 um) and the scattered intensity was measured at multiple exposure times. Using our recent multi-speckle DCS algorithm [16,43,44], the auto-correlation was measured and fitted to obtain the blood flow at an SD separation of 0.4 cm. The $\mu _a$ and $\mu _s^{'}$ were fixed as $0.1 cm^{-1}$ and $12 cm^{-1}$, respectively [2,45].
2.3.5 In vivo human experiment
A DCS system comprising a 785 nm laser (CW Thorlabs, L785P90) coupled to fibre that delivers light to the human forehead. At an SD separation of 2 cm, a single-mode fibre coupled to a single photon avalanche diode (SPAD, SPCM-AQRH-14 FC, Excelitas technologies) and a hardware autocorrelator (correlator.com) measures the intensity autocorrelation. The autocorrelation measured was used as the input to the SDE algorithm to generate the speckle intensity. The $\mu _a$ and $\mu _s^{'}$ were fixed as $0.12 cm^{-1}$ and $12 cm^{-1}$ respectively and the $D_B$ value obtained by fitting to the semi-infinite solution of CDE was $7.2 \times 10^{-9} cm^2/s$.
2.4 Noise model for speckle intensity pattern
Usually, for the deep tissue blood flow imaging (as in DCS), the auto-correlation is directly simulated either using Monte-Carlo simulation [46] or by solving CDE [37]. Instead, the proposed method simulates the speckle intensity (I(t)) using SDE. Hence, the appropriate noise model used in this work is the noise associated with the detection of intensity. For this, we use a camera-based noise model as follows: The noise model takes into account of quantum efficiency - QE, full well capacity - FWC, number of bits - n, pixel size - ps, dark current - $I_d$, exposure time - T, gain - g for a given sensor and irradiance - $I_r$. The model accounts for photon shot noise, dark noise, dark current shot noise and analogue to digital conversion (ADC) approximation based on Ref. [47]. The pseudocode for the same is given in Algorithm S1 of Supplemental document .
For the camera as a detector, the following specifications were used based on our recent work Ref: [16,45], Photometrics Prime BSI, with $QE-72\%$, pixel size $ps=5.6 \mu m$, $I_d$=68.6 $e^-$, FWC=45000$e^-$ and n=16). When an APD was used as a detector, the following parameters were used, QE=20$\%$, $Id=10^{-8} e^-$, gain=50. For APD-based detection, we did not integrate the speckle intensity for a given exposure time, and there was no quantization or thresholding based on FWC. For both camera and APD, note that we have predominantly used effects of quantum efficiency, dark noise, shot noise and gain of the detector, while we have neglected other types of noise such as readout noise, trans-impedance noise, fixed pattern noise, and source follower noise.
3. Results
3.1 Simulation of laser speckles for surface blood flow imaging
The simulated speckle intensity with the SDE algorithm for superficial blood flow are shown in Fig. 1. The sample is assumed to be uniformly illuminated by a laser source with a flow which gives rise to dynamic speckles with spatially independent simple exponential autocorrelation [1] given as $e^{-t/\tau _c}$, where $t$ is the correlation delay time and $\tau _c$ is the characteristics or decorrelation decay time. A sample discretised to $100 \times 100$ grid (representing pixels) with a high flow of $\tau _c=10^{-5} s$ is simulated in the centre channel and a $\tau _c=10^{-3} s$ is simulated on the sides, as shown in Fig. 1(a). We apply the SDE algorithm to every pixel to generate the speckle intensity both as a function of spatial variable (pixels) and time, as shown in Figs. 1(b) and 1(c), respectively. To confirm that the generated speckles have predetermined PDF and autocorrelation, we plot the histogram and intensity autocorrelation in Figs. 1(d) and (e), respectively. It is evident from the plots that the histogram fits a negative exponential PDF with an appropriate mean, and the autocorrelation fits a simple exponential decay with appropriate $\tau _c$ as desired. Here, the normalized intensity auto-correlation $g_2$ was computed based on the Siegert’s relation $g_2=1+\beta |g_1^2|$, where $g_1$ is the normalized field autocorrelation function and $\beta =1$ is the instrumentation constant. To confirm that every pixel obeys the desired autocorrelation, we fit the intensity autocorrelation $g_2$ for $\tau _c$ as shown in Fig. 1(g) and its inset figure, which agrees with the actual $\tau _c$ as shown in Fig. 1(f). The speckle intensity integrated to an exposure time of 1 ms are shown in Visualization 1 of the supplementary section as a video. To check the repeatability, we simulated speckles for 50 times, for a predetermined flow and mean intensity, which resulted in a mean percentage error of $2.8\%$ and $3\%$ for $\tau _c$ and mean intensity respectively. In addition to auto-correlation, we also validate the intensity of speckles in terms of speckle contrast. Figure 1(h) and Fig. 1(i) shows the speckle contrast produced spatially [1] (using a $3 \times 3$ window) and temporally [12] (over 100 frames), respectively, for the data generated at $100 \mu s$ exposure time, both indicating the speckle contrast has lower values when the flow is high at the centre.
It is well known that the sum of random variables with negative exponential distributions gives rise to a resultant random variable with gamma distribution [33]. The same concept applies to speckles where we sum up several speckles to get a resultant speckle with gamma distribution which is also associated with a decrement of instrumentation factor $\beta$ [26] to be less than one. In such a scenario, the autocorrelation of the resultant speckles will be a weighted sum of individual autocorrelations (of component speckles) [38]. Using the SDE algorithm, we simulate several speckle intensity and perform a summation of speckles, as shown in Fig. 2. The decrement in the instrumentation factor $\beta$ also summarizes the importance of the appropriate speckle to pixel mapping in the laser speckle based blood flow measurement systems, and our results of decrement in $\beta$ when the speckles from multiple pixels are summed up is consistent with earlier literature [26].
3.1.1 Simulation of different flow models for various tissue types in LSCI
Different flow models having more complex field autocorrelation (in contrast to simple exponential relations) that may arise in different parts of the tissues, as explained in Ref [6], are simulated by the SDE method. Using the above-mentioned principles based on the addition of speckles, as depicted in Fig. 2, the complex auto-correlation models (as given in Eq. (5) in the Methods section) were generated for different cases, such as with and without static layers [29], with different $n$ and $d$ values based on the tissue type. The autocorrelation of the speckle intensity generated using SDE is compared against the desired auto-correlation models as shown in Fig. 3 with a reasonably good fit, thus indicating the capability of the proposed method to model complex blood flow models.
3.1.2 In-vivo experimental validation of laser speckle contrast imaging
For in-vivo validation of the proposed algorithm, we measured the surface blood flow in the mice brain using our recently proposed Multi-speckle imaging system [16,43] as shown in Fig. 4(a). The system uses a low-frame rate camera to measure the speckle intensity at multiple exposure times and compute the intensity autocorrelation, which is related to blood flow. The schematic of the experimental setup is shown in Fig S1 (a) of the supplemental section. For the simulation of speckles, the intensity image of the region of interest in the brain at 1.1 ms exposure time was used as the mean intensity at every pixel. The characteristic decay time $\tau _c$ determined experimentally is shown in Fig. 4(b). For two representative pixels, the multi exposure speckle contrast data obtained by the experiment is compared with the simulated speckle contrast data. The corresponding $\tau _c$ values are reported in the legend of inset of Fig. 4(b) and are comparable to one another. The SDE based intensity data at $1.1 ms$ exposure time is shown in Fig. 4(c). Visualization 2 (Supplemental section) shows the simulated speckle intensity fluctuations in mice brain. The temporal speckle contrast obtained from the simulated SDE-based speckle intensity is shown for a given exposure time of $T=1.1 ms$ is compared with the experimental data.
3.2 Simulation of laser speckles for deep tissue blood flow imaging
The schematic of the DCS system with a pointed source used for simulation is shown in Fig. 5(a). A laser source is pointed to the centre of the object, and the light diffuses symmetrically in all directions. The sample is assumed to be homogenous, and the optical properties used for the simulation are $\mu _a=0.1 cm^{-1}$, $\mu _s^{'}=0.12 cm^{-1}$ and the dynamic property (particle diffusion coefficient) used is $D_B=10^{-8} cm^2/s$. An infinite solution of Diffusion equation [2,48] and semi-infinite CDE [2,37] are used to generate the intensity and autocorrelation, respectively, for the simulated speckles. The SDE-based algorithm was used to simulate the speckle intensity and the integrated intensity $(\phi )$ at an exposure time of T=1ms as shown in Fig. 5(c), which is comparable to the desired intensity as shown in Fig. 5(b). The log of the intensity is represented in Fig. 5(b, c, e) for better clarity. The intensity profile as a function of SD separation $(log(\phi.r))$ is shown in Fig. 5(d) along with a linear fit. It is well known that the slope of the linear fit should approximately follow $k=\sqrt {3\mu _a \mu _s^{'}} =2.08$, which is the desired value, whereas the slope obtained from the SDE-based intensity speckle was 2.1. The speckle intensity at an exposure time of $10 \mu s$ is shown in Fig. 5(e) (Refer to the supplementary video Visualization 3 for the intensity of speckle fluctuations as a function of time and space). To validate the proposed method, we have compared the histogram (PDF) and the autocorrelation of the speckles obtained by the simulation against their desired values. Figure 5(f) compares the SDE-based intensity autocorrelation at two representative SD separations (r=1cm and r=3 cm) against the desired auto-correlation, which are in reasonable comparison with one another. Figure 5(g) shows the histogram of the speckle intensity obtained at a given representative SD separation which fits a gamma distribution with the desired scale parameter (Desired $\mu =2.26 \times 10^{-15}$ and the fitted $\mu =1.85 \times 10^{-15}$). Figures 5(h) and 5(i) shows the intensity auto-correlations obtained at two delay times $\tau =10 \mu s$ and $\tau =100 \mu s$ respectively. In addition, we also plot the computed speckle contrast as a function of multiple exposure times for two different SD separations ($r=1 cm$ and $3 cm$) in Fig. 5(j). These results show that the SDE algorithm can be used for the generation of both uniform as well as diffusing speckle intensity with applications in conventional laser speckle imaging and deep-tissue imaging, respectively.
3.2.1 In-vivo experimental validation of diffusing speckles generation
The capability of the SDE algorithm to generate speckle intensity comparable to $in-vivo$ experimental studies is shown in Fig. 6. The results show the in-vivo validation performed in both mice and adult human brain for the speckles generation using the SDE algorithm. The experimental method is detailed in the methods section. Figure 6(a) shows the comparison between the experimental intensity autocorrelation at an SD separation of 0.4 cm in the brain of a mouse in alive and sacrificed condition against the SDE-based intensity autocorrelation. The $g_2$ curves obtained experimentally were fitted to the semi-infinite solution of CDE, which resulted in $D_B$ of $1.57 \times 10^{-9} cm^2/s$ and $1.49 \times 10^{-11} cm^2/s$ in alive and sacrificed condition. The auto-correlation of speckle intensity generated by SDE resulted in $D_B$ of $1.45 \times 10^{-9} cm^2/s$ and $1.53 \times 10^{-11} cm^2/s$. Figure 6(b) shows the experimental autocorrelation in an adult human head at an SD separation of 2 cm compared to the SDE-based intensity speckle simulation. These results indicate that SDE-based algorithms can be used to obtain speckle intensity for deep tissue blood flow comparable with the experimental results based on $in-vivo$ animal and human studies.
3.3 Modelling of speckle intensity with realistic noise models
The speckle intensity obtained using SDE were further added with detector noise for more realistic analysis. For this, we assume that our detectors were (a) a camera and (b) an APD [2]. If the detector is a camera, we use the camera noise model for the speckle intensity spatially. For the APD, we apply the noise to temporal intensity fluctuations at a given pixel. Figure 7(a) shows the speckle intensity with a camera noise model with dark and shot noise effects. The speckle intensity deviated from the linear path when dark noise was added. Figure 7(b) shows the effect of the noise on an APD, wherein the noise was added to the speckle intensity for a given SD pair. It can be seen that SDE-based autocorrelation without noise follows the desired autocorrelation, while speckle intensity with added noise has the autocorrelation deviating slightly from the desired value. It is also to be noted that the effect of noise is evident at larger values of delay time $\tau$ than that of earlier times as the number of samples available is large at earlier delay times which may results in averaging of noise. The effect of the noises added to speckle intensity is more realistic than that of the noise models used earlier in DCS, as the DCS usually follows the noise models for auto-correlation, while the noises should be applied to the fundamental quantity which is speckle intensity. Figure 7(c) shows the camera noise effects on normalized speckle contrast $(\kappa _N)$. The speckle contrast fits the initial exposure times; however, as the exposure time increases, due to inadequate dynamic range or full-well capacity of the camera, it saturates, which is evident from the plot.
Although we have used a $3 \times 3$ window for obtaining spatial speckle contrast as shown in Fig. 1(h); as reported in Ref [49], spatial window of different sizes may used for optimal computation of speckle contrast. In case of calculation of temporal speckle contrast as shown in Fig. 1(i), it is well known that, more the number of pixels/ frames used, better would be the SNR of the speckle contrast, as shown in Ref [50] in the noise model for speckle contrast. Figure 7(d) shows the mean and the variance of the multi-exposure speckle contrast calculated by varying the number of speckles. As the number of speckles used to calculate speckle contrast increases, the standard deviation decreases and mean speckle contrast approaches the analytical solution.
4. Discussion
This paper proposes a comprehensive tool for laser intensity speckle simulation that uses a novel stochastic differential equations-based algorithm. The proposed algorithm generates speckle intensity that follow a desired PDF and autocorrelation function for both surface and deep tissue blood flow. The speckle intensity produced by SDE based algorithm were validated by extensive simulation studies and in-vivo experiments measuring speckles from mice brain and an adult human head.
For the surface blood flow, we simulate the speckle intensity with a simple field autocorrelation given by $g_1=exp(-t/\tau _c)$ for different values of $\tau _c$ as given in Fig. 1. We validate the SDE algorithm by comparing the two main characteristics of the speckle intensity [33]: the PDF of the speckles and autocorrelation, against their respective desired values as shown in Fig. 1(d,e). In addition, we generate speckle contrast [1] (Fig. 1(h)), which is one of the widely used quantifiers in measuring surface blood flow. We also validate the LSCI system for different tissue types, such as large vessels, small vessels and parenchyma types for different values of static scatterers and different orders of the diffusing regime based on different flow models. This shows that SDE could be a potential tool for the simulation of speckle intensity for different flow models, as it is established in Ref [6] that different tissue types follow different flow models in measuring surface blood flow. An in-vivo validation of the proposed SDE algorithm in mice brain is shown in Fig. 4 by comparing simulated speckles against the in-vivo experimental speckle measurements. Thus, the SDE algorithms could be used as a robust tool for generating speckle intensity for LSCI with different flow models.
In order to simulate the deep tissue blood flow speckles, we combine the SDE with the solution of the photon diffusion equation along with appropriate PDF and autocorrelation. The infinite solution of the photon diffusion equation is used for the generation of the mean intensity (or initial intensity $\mu$ in Eq. (3)). Similarly, the semi-infinite solution of CDE was used to select the normalized field autocorrelation function [2], which is not a simple autocorrelation function of exponential type. In order to accommodate complicated (not simple exponential) autocorrelation, we adopt the method of approximating such functions using weighted series of simple exponentials. The speckles simulated for each such simple exponentials are appropriately added to get the speckles with desired PDF and autocorrelation. It can be deduced that such speckles obey gamma distribution, which is the case for the sum of several spatially averaged speckles due to the finite size of the detector [50,51]. The simulation results as shown in Fig. 5 and the in-vivo results on mice and healthy human adults as shown in Fig. 6 validate the proposed algorithm. Additionally, the multi-exposure speckle contrast, a derived quantity used by several modalities was also generated from the SDE based speckle intensity and its validation is shown in Fig. 5(k).
In order to simulate speckles in a more realistic manner, we have also incorporated the noise models associated with commonly employed detectors for speckle measurement. The noise model considers critical detector parameters such as quantum efficiency, dark noise, pixel size / active detection area, gain, full well capacity, dynamic range, and quantization noise. The proposed noise models are directly associated with the speckle intensity rather than the existing noise models that have been developed for either the autocorrelation or the speckle contrast, which are derived quantities from the speckle as in LSCI/DCS [50,52]. We have not included the effects laser coherence length, line width of the laser and numerical aperture of the camera lens in the model. It is shown in Ref [53] that, volume holographic grating stabilized laser diode and other laser sources whose spectral width is less than 1 nm improves the speckle contrast. In addition, a complete theoretical study on speckle contrast and laser coherence length and imaging characteristics of the optical system is reported in Ref [54]. These models could be incorporated to model the effect of coherence length and line width in order to access its effect on the speckle intensity produced by the proposed method.
As mentioned in Section 2.3, we have used $t_{min}=10^{-6} s,t_{max}=1s$ and $N=10^7$ samples to generate the speckle intensity. In this time discretization, within $5\%$ error in blood flow obtained using intensity autocorrelation, we can simulate the flow or $\tau _c$ in the range of $5\mu s$ to $40 ms$. In principle, the SDE algorithm could simulate speckle intensity without any constraints on blood flow. However, considering the practical computational constraints such as the numerical instability of the algorithm and the requirement of a sufficient number of samples for the calculation of intensity autocorrelation, the above range is found to be optimal. However, we can increase the velocity range by varying the time range of the simulations. Additionally, to check the SDE algorithm’s repeatability, we produced speckle intensity for a predetermined flow and mean intensity over 50 trails. We found that the mean error in the flow and mean intensity produced is around $2.8\%$ and $3\%$, respectively, which shows that the algorithm is repeatable. The noise characteristics of the camera, such as, sensor full well capacity, quantum efficiency, and signal-to-noise ratio, are added later after the simulation of speckle intensity. The above-mentioned parameters determines the usable range of flow thus mimicking the realistic case. For instance, for the given camera (Basler acA-640-120 um), based on its quantum efficiency, full-well capacity and noise characteristics, $\tau _c$ should be within the range of $0.1$ to $10 ms$, beyond which the intensity saturates due to the finite dynamic range of the camera. It is crucial to use the unbiased auto-correlation as given in the Methods section for auto-correlation calculation. In case of biased/scaled auto-correlation models are used, it could lead to bias in the larger values of $\tau$, as shown in Fig S3 in supplemental document. The bias is caused by the insufficient number of samples for the auto-correlation at larger values of $\tau$, whose effects are removed when an unbiased auto-correlation is used. Additionally, it is evident that we have generated temporal intensity speckles at every pixels and not spatially across multiple pixels. Under the assumption of ergodicity, the temporal speckle contrast and spatial speckle contrast are equivalent. However, further mathematical analysis, as given in Ref [42,55], would be required to ensure the same.
Although there are a few methods to produce speckle intensity, as given in Ref [26], most do not provide a means to apriori select complex (non-exponential) autocorrelation. Even those which could incorporate autocorrelation statistics in blood flow are usually a simple exponential function as in Ref [28,30,56]. We have introduced the complex autocorrelation models for both surface as well as deep tissue blood flow, by approximation using finite weighted exponential series. This allows us to model speckles from any tissues with multiple layers of optical and dynamical properties. This is essential, as the inappropriate forward modelling leads to significantly large errors in the reconstruction of blood flow [57,58] as encountered in modelling the flow in the human head and animal brain with multiple layers [59].
The proposed algorithm could simulate speckle intensity of any desired auto-correlation and speckle statistics by appropriately selecting the parameters $a$ and $b$ of Eq. (3). Although the proposed method in this paper, shows the applicability of producing speckle intensity for laser speckle contrast imaging, it could also be easily extended for other speckle statistics, which are not exponential or gamma distribution [60]. Moreover, although we have shown the commonly used reflection geometry for validation, simulations can also be implemented for transmission geometry which gives better signal-to-background ratio [61]. Additionally, in near future, the proposed laser simulation tool could be potentially used to compare multiple laser speckle imaging modalities and could optimize several parameters such as exposure time, number of frames/ speckles, optode positions, which are key factors for both LSCI and DCS based blood flow measurement systems.
The SDE-based speckle simulation for turbid medium like tissue is a new approach based on stochastic process whereas the most of the existing methods use statistical tools, ODE’s or Fourier-based approaches [26–28]. We have validated the proposed method using extensive simulation and in-vivo studies, thus proving it to be a robust, comprehensive tool to simulate speckle intensity for both surface and deep tissue flow imaging applications. We believe this approach would eventually help adopt stochastic differential equations to model light propagation in a dynamic turbid medium like tissue.
Funding
Science and Engineering Research Board (Early Career Research Award -2016); Department of Biotechnology, Ministry of Science and Technology, India (Ramalingaswamy Fellowship-2016); Indian Institute of Technology Bombay (Seed Grant, Wadhwani Research Center for Bioengineering (WRCB)).
Disclosures
The authors have no conflict of interest.
Data availability
Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.
Supplemental document
See Supplement 1 for supporting content.
References
1. D. A. Boas and A. K. Dunn, “Laser speckle contrast imaging in biomedical optics,” J. Biomed. Opt. 15(1), 011109 (2010). [CrossRef]
2. T. Durduran, R. Choe, W. Baker, and A. G. Yodh, “Diffuse optics for tissue monitoring and tomography,” Rep. Prog. Phys. 73(7), 076701 (2010). [CrossRef]
3. J. Briers and A. Fercher, “Retinal blood-flow visualization by means of laser speckle photography,” Investigative ophthalmology & visual science 22(2), 255–259 (1982).
4. A. B. Parthasarathy, W. J. Tom, A. Gopal, X. Zhang, and A. K. Dunn, “Robust flow measurement with multi-exposure speckle imaging,” Opt. Express 16(3), 1975–1989 (2008). [CrossRef]
5. W. Liu, R. Qian, S. Xu, P. Chandra Konda, J. Jönsson, M. Harfouche, D. Borycki, C. Cooke, E. Berrocal, Q. Dai, H. Wang, and R. Horstmeyer, “Fast and sensitive diffuse correlation spectroscopy with highly parallelized single photon detection,” APL Photonics 6(2), 026106 (2021). [CrossRef]
6. D. D. Postnov, J. Tang, S. E. Erdener, K. Kılıç, and D. A. Boas, “Dynamic light scattering imaging,” Sci. Adv. 6(45), eabc4628 (2020). [CrossRef]
7. J. D. Johansson, D. Portaluppi, M. Buttafava, and F. Villa, “A multipixel diffuse correlation spectroscopy system based on a single photon avalanche diode array,” J. Biophotonics 12(11), e201900091 (2019). [CrossRef]
8. T. Dragojević, D. Bronzi, H. M. Varma, C. P. Valdes, C. Castellvi, F. Villa, A. Tosi, C. Justicia, F. Zappa, and T. Durduran, “High-speed multi-exposure laser speckle contrast imaging with a single-photon counting camera,” Biomed. Opt. Express 6(8), 2865–2876 (2015). [CrossRef]
9. T. Durduran and A. G. Yodh, “Diffuse correlation spectroscopy for non-invasive, micro-vascular cerebral blood flow measurement,” NeuroImage 85, 51–63 (2014). [CrossRef]
10. J. P. Culver, T. Durduran, D. Furuya, C. Cheung, J. H. Greenberg, and A. 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]
11. R. Bi, J. Dong, and K. Lee, “Deep tissue flowmetry based on diffuse speckle contrast analysis,” Opt. Lett. 38(9), 1401–1403 (2013). [CrossRef]
12. C. P. Valdes, H. M. Varma, A. K. Kristoffersen, T. Dragojevic, J. P. Culver, and T. Durduran, “Speckle contrast optical spectroscopy, a non-invasive, diffuse optical method for measuring microvascular blood flow in tissue,” Biomed. Opt. Express 5(8), 2769–2784 (2014). [CrossRef]
13. H. M. Varma, C. P. Valdes, A. K. Kristoffersen, J. P. Culver, and T. Durduran, “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]
14. W. Zhou, O. Kholiqov, S. P. Chong, and V. J. Srinivasan, “Highly parallel, interferometric diffusing wave spectroscopy for monitoring cerebral blood flow dynamics,” Optica 5(5), 518–527 (2018). [CrossRef]
15. 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]
16. K. Murali and H. M. Varma, “Multi-speckle diffuse correlation spectroscopy to measure cerebral blood flow,” Biomed. Opt. Express 11(11), 6699–6709 (2020). [CrossRef]
17. R. Paul, K. Murali, S. Chetia, and H. M. Varma, “High-density diffuse correlation tomography with enhanced depth localization and minimal surface artefacts,” Biomed. Opt. Express 13(11), 6081–6099 (2022). [CrossRef]
18. W. Zhou, M. Zhao, O. Kholiqov, and V. J. Srinivasan, “Multi-exposure interferometric diffusing wave spectroscopy,” Opt. Lett. 46(18), 4498–4501 (2021). [CrossRef]
19. J. Xu, A. K. Jahromi, J. Brake, J. E. Robinson, and C. Yang, “Interferometric speckle visibility spectroscopy (isvs) for human cerebral blood flow monitoring,” APL Photonics 5(12), 126102 (2020). [CrossRef]
20. A. Biswas, A. Takshi, and A. B. Parthasarathy, “A method for improving the dynamic range of integrated diffuse speckle contrast spectroscopy,” in Optics and the Brain, (Optical Society of America, 2021), pp. BTh1B–4.
21. E. James and S. Powell, “Fourier domain diffuse correlation spectroscopy with heterodyne holographic detection,” Biomed. Opt. Express 11(11), 6755–6779 (2020). [CrossRef]
22. M. Seong, Y. Oh, K. Lee, and J. G. Kim, “Blood flow estimation via numerical integration of temporal autocorrelation function in diffuse correlation spectroscopy,” Comput. Methods Programs Biomed. 222, 106933 (2022). [CrossRef]
23. Z. Hajjarian, H. T. Nia, S. Ahn, A. J. Grodzinsky, R. K. Jain, and S. K. Nadkarni, “Laser speckle rheology for evaluating the viscoelastic properties of hydrogel scaffolds,” Sci. Rep. 6(1), 37949 (2016). [CrossRef]
24. H. M. Varma, K. P. Mohanan, N. Hyvönen, A. K. Nandakumaran, and R. M. Vasu, “Ultrasound-modulated optical tomography: recovery of amplitude of vibration in the insonified region from boundary measurement of light correlation,” J. Opt. Soc. Am. A 28(11), 2322–2331 (2011). [CrossRef]
25. R. S. Chandran, S. Sarkar, R. Kanhirodan, D. Roy, and R. M. Vasu, “Diffusing-wave spectroscopy in an inhomogeneous object: local viscoelastic spectra from ultrasound-assisted measurement of correlation decay arising from the ultrasound focal volume,” Phys. Rev. E 90(1), 012303 (2014). [CrossRef]
26. S. J. Kirkpatrick, D. D. Duncan, and E. M. Wells-Gray, “Detrimental effects of speckle-pixel size matching in laser speckle contrast imaging,” Opt. Lett. 33(24), 2886–2888 (2008). [CrossRef]
27. D. D. Duncan and S. J. Kirkpatrick, “The copula: a tool for simulating speckle dynamics,” J. Opt. Soc. Am. A 25(1), 231–237 (2008). [CrossRef]
28. L. Song, Z. Zhou, X. Wang, X. Zhao, and D. S. Elson, “Simulation of speckle patterns with pre-defined correlation distributions,” Biomed. Opt. Express 7(3), 798–809 (2016). [CrossRef]
29. H. J. Rabal and R. A. Braga Jr, Dynamic Laser Speckle and Applications (CRC press, 2018).
30. E. James, S. Powell, and P. Munro, “Simulation of statistically accurate time-integrated dynamic speckle patterns in biomedical optics,” Opt. Lett. 46(17), 4390–4393 (2021). [CrossRef]
31. C. A. Braumann, Introduction to Stochastic Differential Equations with Applications to Modelling in Biology and Finance (John Wiley & Sons, 2019).
32. C. Gardiner, “Stochastic methods: a handbook for the natural and social sciences 4th ed. (2009),”.
33. J. W. Goodman, Speckle Phenomena in Optics: Theory and Applications (Roberts and Company Publishers, 2007).
34. G. N. Mil’shtein, “Approximate integration of stochastic differential equations,” Teoriya veroyatnostei i ee primeneniya 19, 583–588 (1974).
35. G. Mil’shtejn, “Approximate integration of stochastic differential equations,” Theory Probab. & Its Appl. 19(3), 557–562 (1975). [CrossRef]
36. Z. Yin and S. Gan, “An improved milstein method for stiff stochastic differential equations,” Adv. Differ. Eqs. 2015(1), 369 (2015). [CrossRef]
37. D. A. Boas, L. Campbell, and A. G. Yodh, “Scattering and imaging with diffusing temporal field correlations,” Phys. Rev. Lett. 75(9), 1855–1858 (1995). [CrossRef]
38. R. Zárate-Mi nano and F. Milano, “Construction of sde-based wind speed models with exponentially decaying autocorrelation,” Renewable Energy 94, 186–196 (2016). [CrossRef]
39. D. J. Higham, “An algorithmic introduction to numerical simulation of stochastic differential equations,” SIAM Rev. 43(3), 525–546 (2001). [CrossRef]
40. E. M. Buckley, A. B. Parthasarathy, P. E. Grant, A. G. Yodh, and M. A. Franceschini, “Diffuse correlation spectroscopy for measurement of cerebral blood flow: future prospects,” Neurophotonics 1(1), 011009 (2014). [CrossRef]
41. B. M. Bibby, I. M. Skovgaard, and M. Sørensen, “Diffusion-type models with given marginal distribution and autocorrelation function,” Bernoulli 11(2), 191–220 (2005). [CrossRef]
42. O. Ibe, Markov Processes for Stochastic Modeling (Newnes, 2013), Chapter 1, 2 and chapter 9.
43. K. Murali, A. 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]
44. K. Murali, A. 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]
45. 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]
46. Q. Fang and D. A. Boas, “Monte carlo simulation of photon migration in 3d turbid media accelerated by graphics processing units,” Opt. Express 17(22), 20178–20190 (2009). [CrossRef]
47. M. Konnik and J. Welsh, “High-level numerical simulations of noise in ccd and cmos photosensors: review and tutorial,” arXiv, arXiv:1412.4031 (2014). [CrossRef]
48. S. R. Arridge, M. Cope, and D. Delpy, “The theoretical basis for the determination of optical pathlengths in tissue: temporal and frequency analysis,” Phys. Med. Biol. 37(7), 1531–1560 (1992). [CrossRef]
49. D. D. Duncan, S. J. Kirkpatrick, and R. K. Wang, “Statistics of local speckle contrast,” J. Opt. Soc. Am. A 25(1), 9–15 (2008). [CrossRef]
50. S. Yuan, “Sensitivity, noise and quantitative model of laser speckle contrast imaging,” Ph.D. thesis, Tufts University (2008).
51. J. C. Dainty, Laser Speckle and Related Phenomena, vol. 9 (Springer science & business Media, 2013).
52. 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]
53. D. D. Postnov, X. Cheng, S. E. Erdener, and D. A. Boas, “Choosing a laser for laser speckle contrast imaging,” Sci. Rep. 9(1), 2542 (2019). [CrossRef]
54. H. Fujii and T. Asakura, “A contrast variation of image speckle intensity under illumination of partially coherent light,” Opt. Commun. 12(1), 32–38 (1974). [CrossRef]
55. E. Parzen, “Conditions that a stochastic process be ergodic,” Ann. Math. Statist. 29(1), 299–301 (1958). [CrossRef]
56. X. Cheng, H. Chen, E. J. Sie, F. Marsili, and D. A. Boas, “Development of a monte carlo-wave model to simulate time domain diffuse correlation spectroscopy measurements from first principles,” J. Biomed. Opt. 27(08), 083009 (2022). [CrossRef]
57. L. Gagnon, M. Desjardins, J. Jehanne-Lacasse, L. Bherer, and F. Lesage, “Investigation of diffuse correlation spectroscopy in multi-layered media including the human head,” Opt. Express 16(20), 15514–15530 (2008). [CrossRef]
58. 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]
59. M. M. Wu, K. Perdue, S.-T. Chan, K. A. Stephens, B. Deng, M. A. Franceschini, 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]
60. N. Bender, H. Yılmaz, Y. Bromberg, and H. Cao, “Customizing speckle intensity statistics,” Optica 5(5), 595–600 (2018). [CrossRef]
61. D.-Y. Li, Q. Xia, T.-T. Yu, J.-T. Zhu, and D. Zhu, “Transmissive-detected laser speckle contrast imaging for blood flow monitoring in thick tissue: from monte carlo simulation to experimental demonstration,” Light: Sci. Appl. 10(1), 241 (2021). [CrossRef]