## 1.

## Introduction

Laser speckle contrast imaging (LSCI) is a full-field optical technique used for imaging blood flow changes *in vivo* without scanning.^{1, 2, 3, 4, 5, 6, 7} The speckle pattern is generated by random interference of the scattered laser light and is modulated by motion of the scatterers. LSCI accesses the speed information of the scatterers by calculating the speckle contrast
$\left(K\right)$
,^{8} which is defined as the ratio of the standard deviation of intensity to the mean value of intensity in the speckle patterns.

As illustrated in Fig. 1
, the laser speckle contrast analysis can be performed based on spatial statistics,^{2} temporal statistics,^{9, 10} or a combination of both.^{3, 11, 12} The laser speckle spatial contrast analysis (LSSCA) method^{2} performs speckle contrast calculation in the spatial domain using a spatial window. LSSCA achieves high temporal resolution with the loss of spatial resolution, impeding its application on monitoring blood flow changes in small vessels. The laser speckle temporal contrast analysis (LSTCA) method,^{9, 10} which is based on temporal statistics, computes speckle contrast images using a sequence of speckle images acquired along a few time points instead of using a spatial window. LSTCA preserves the original spatial resolution by sacrificing the temporal resolution, making it inappropriate in applications where video frame rate visualization of blood flow is required. To compromise spatial and temporal resolution, several methods that attempted to combine spatial and temporal statistics have been presented. For each of these methods, a total number of
${N}_{s}\times {N}_{s}\times N$
pixels is used for calculating a
$K$
value, where
${N}_{s}$
is the length of a square matrix in the spatial domain and
$N$
is the number of frames. The spatial-based laser speckle contrast analysis (SLASCA) method^{3} is an improvement of LSSCA, in which
$N$
frames of speckle contrast images computed by LSSCA with a
${N}_{s}\times {N}_{s}$
spatial window are temporally averaged. By performing a temporal average, SLASCA displays much lower image noise than LSSCA does. The temporal-based laser speckle contrast analysis (TLASCA) method^{11} is a modification of LSTCA, in which a speckle contrast image computed by LSTCA from
$N$
frames of raw speckle images is spatially averaged using a
${N}_{s}\times {N}_{s}$
square matrix (such as
$3\times 3\phantom{\rule{0.3em}{0ex}}\text{pixels}$
). By performing a spatial average, TLASCA obtains a speckle contrast image with acceptable signal-to-noise ratio (SNR) with fewer number of frames than LSTCA does. Recently, a spatiotemporal laser speckle contrast analysis (STLASCA) method^{12} has been presented, in which speckle contrast is calculated directly from an
${N}_{s}\times {N}_{s}\times N$
pixel cube. Because all the
${N}_{s}\times {N}_{s}\times N$
pixel are used to compute a
$K$
value, STLASCA achieves high *SNR* without a further averaging procedure as TLASCA or SLASCA does.

However, the focus of previous studies has been on the processing time,^{11, 12} but the statistical accuracy of these methods, which is important in LSCI, has not yet been thoroughly compared. Since the introduction of some parallel data processing devices such as the graphics processing unit into LSCI,^{13} the computational efficiency has been significantly promoted, thus the data processing time should no longer be the primary consideration. In contrast, the close correlation of speckle contrast value with speed of motion determines that the statistical accuracy of a method for computing speckle contrast is very important. Through numerical simulation and *in vivo* rat cortical blood flow imaging, this study investigates the mean speckle contrast values and the relative noises of the speckle contrast images computed by SLASCA, TLASCA, and STLASCA.

## 2.

## Methods and Materials

## 2.1.

### Numerical Simulation of Time-Integrated Dynamic Speckles

Numerical simulation of time-integrated dynamic speckle images is performed to compare the prior three methods. The simulating program is written in MATLAB language. A sequence of statistically independent
$Z$
values following Gaussian distribution have being generated by^{14}

## Eq. 1

$$Z\left(k\right)=\sqrt{-2\phantom{\rule{0.2em}{0ex}}\mathrm{ln}\phantom{\rule{0.2em}{0ex}}{X}_{1}}\phantom{\rule{0.2em}{0ex}}\mathrm{cos}(2\pi {X}_{2}+\frac{\pi}{2}\cdot \frac{k-1}{n-1}),$$^{14}for $Z$ values, 50 frames of uniformly distributed $T$ values on the unit interval that are statistically correlated are obtained. Subsequently, 50 frames of statistically correlated fully developed speckle images are obtained by performing fast Fourier transformation of the random phasor, which is

^{14}where $m$ is a multiplicative factor that affects the correlation time of the fully developed speckle image sequence, which will then affect the decorrelation time ${\tau}_{c}$ of the time-integrated dynamic speckle images. Due to the fact that most of the speckle contrast values of a speckle contrast image are less than 0.6 in biomedical applications of LSCI,

^{15}in this study $m$ is set to be 12 to shorten the correlation time. The minimum speckle size of the speckle images is chosen to be two pixels to maximize the speckle contrast.

^{16}Then one frame of a time-integrated dynamic speckle image is obtained by averaging 25 consecutive frames of the generated fully developed speckle images. To generate 30 frames of time-integrated dynamic speckle images, the prior procedures are repeated. It is worth to note that what the random number generators produce are not random numbers, but pseudorandom numbers. If the seed of a generator is the same as the previous one, then the produced pseudorandom numbers will be the same as the previous ones. To ensure statistical independence among the time-integrated dynamic speckle images, for each time a time-integrated dynamic speckle image frame is simulated, and different seeds for producing ${X}_{1}$ and ${X}_{2}$ are used.

## 2.2.

### Data Analysis of the Simulated Time-Integrated Dynamic Speckle Sequence

One frame of the simulated time-integrated dynamic speckle images is shown in Fig. 2
. Due to the effect of time integration, the speckle image is to some degree blurred. Figure 2 shows the intensity probability density function (PDF) of one of the simulated fully developed speckle images (green solid line), the theoretical negative exponential intensity PDF curve for a fully developed speckle image^{8} (black dashed line), and the intensity PDF of one of the simulated time-integrated dynamic speckle images (red solid line). As is shown, the intensity PDF of the simulated fully developed speckle image agrees with the theoretical negative exponential curve very well, demonstrating the validity of the simulation method. As is also shown, intensity PDF of the simulated time-integrated dynamic speckle image is no longer negative-exponential shaped due to time integration.

$K$ values for these simulated speckle images are then computed using SLASCA, TLASCA, and STLASCA. The ${N}_{s}$ and $N$ values are variable. Each method obtains one frame of $K$ image using a certain ${N}_{s}$ and an $N$ value. For each method, the ${\mu}_{K}$ and the ${\sigma}_{K}$ values, which are the mean value and standard deviation of $K$ , respectively, are calculated from the whole $K$ image. Subsequently, the ${\sigma}_{K}\u2215{\mu}_{K}$ value, which is used to quantify the noise level of a speckle contrast image, is obtained. Global speckle contrast values, which are calculated over all the pixels of the speckle image realizations, are obtained as references of speckle contrast values.

## 2.3.

###
*In Vivo* Rat Cortical Blood Flow Imaging

An *in vivo* cortical blood flow imaging experiment is also performed to compare these three methods. The schematic setup for the experiment is shown in Fig. 3
. An adult male Wistar rat weighing around
$200\phantom{\rule{0.3em}{0ex}}\mathrm{g}$
was anesthetized and fixed in a stereotaxic instrument. An approximately
$5.0\times 5.0\text{-}\mathrm{mm}$
cranial window with intact dura was formed by removing the skull overlying one side of the parietal cortex with a high speed dental drill (Fine Science Tools, USA) under constant saline cooling. A beam of He-Ne laser (Melles Griot, America;
$632.8\phantom{\rule{0.3em}{0ex}}\mathrm{nm}$
and
$15\phantom{\rule{0.3em}{0ex}}\mathrm{mW}$
) was expanded and collimated to illuminate the cranial window at about
$30\text{-}\mathrm{deg}$
incidence. 30 frames of statistically independent laser speckle images were acquired by a
$12\text{-bit}$
charge-coupled device (CCD) camera (PixelFly QE, PCO Computer, Germany; pixel
$\text{size}=6.45\times 6.45\phantom{\rule{0.3em}{0ex}}\mu \mathrm{m}$
) attached to a microscope (Z16 APO, Leica, Germany; working distance
$97\phantom{\rule{0.3em}{0ex}}\mathrm{mm}$
) for data processing. The CCD exposure duration was
$20\phantom{\rule{0.3em}{0ex}}\mathrm{ms}$
and the frame interval time is approximately
$87\phantom{\rule{0.3em}{0ex}}\mathrm{ms}$
. The system magnification is adjusted to
$3.15\times $
, and the aperture diaphragm is well controlled to ensure the average speckle size of the images to be approximately two pixels.^{16} A variable attenuator was used in the light path to ensure the light intensity within the dynamic range of the CCD camera. The whole setup was placed on a vibration-isolator table (VH3036W, Newport).

## 2.4.

### Data Analysis of Speckle Images Obtained from *In Vivo* Rat Cortical Blood Flow Imaging

The data were processed by programs written in MATLAB language. An image of the cortex under white light illumination is shown in Fig. 4
. A typical speckle contrast image computed by STLASCA over a
$3\times 3\times 20\phantom{\rule{0.3em}{0ex}}\text{pixels}$
stack is shown in Fig. 4. The broken vessels that can be identified in Fig. 4 disappear in Fig. 4, demonstrating again that laser speckle contrast analysis obtains the speed information of the scatterers, rather than simply the structure pattern of the vessels.^{17} Two small rectangular ROIs,
$R1$
with
$157\times 151\phantom{\rule{0.3em}{0ex}}\text{pixels}$
and
$R2$
with
$107\times 211\phantom{\rule{0.3em}{0ex}}\text{pixels}$
of the speckle images, are selected for analysis. The magnified images of
$R1$
and
$R2$
are shown in Figs. 4 and 4, respectively. The red lines depict the center lines of the selected vessels. The locations of the vessel center lines are obtained by performing cubic polynomial fit of the skeleton of the speckle contrast images with removed spurs. The number of pixels from point
$P1$
to
$P2$
of vessel 1 and from
$P3$
to
$P4$
of vessel 2 along the vessel center lines are both 130, which is approximately
$266\phantom{\rule{0.3em}{0ex}}\mu \mathrm{m}$
in distance. In such a short distance, the blood flow along the vessel center line can be considered as uniform, thus the speckle contrast values have no significant fluctuation. By setting
${N}_{s}$
to 3 and
$N$
to values ranging from 3 to 30, each method obtains two groups of
${\mu}_{K}$
and
${\sigma}_{K}\u2215{\mu}_{K}$
values. The first group of
${\mu}_{K}$
and
${\sigma}_{K}\u2215{\mu}_{K}$
values are calculated from the 130
$K$
values along vessel center line 1, and the second group of
${\mu}_{K}$
and
${\sigma}_{K}\u2215{\mu}_{K}$
values are calculated from the 130
$K$
values along vessel center line 2.

## 3.

## Results

The ${\mu}_{K}$ values as a function of number of frames $N$ are shown in Figs. 5 and 5 . The length of the square matrix ${N}_{s}$ is selected to be three pixels for Fig. 5 and nine pixels for Fig. 5. The ${\mu}_{K}$ values computed by STLASCA increase with either increasing $N$ or increasing ${N}_{s}$ . The ${\mu}_{K}$ values computed by TLASCA increase with increasing $N$ , but remain unchanged with increasing ${N}_{s}$ . The ${\mu}_{K}$ values computed by SLASCA increase with increasing ${N}_{s}$ , but remain unchanged with increasing $N$ . As is shown in either Fig. 5 or Fig. 5, the ${\mu}_{K}$ value computed by STLASCA is closest to the global contrast value for a given $N$ value. The global speckle contrast curve keeps almost a constant of 0.33. Given certain ${N}_{s}$ and $N$ values, the relative error between the mean speckle contrast value calculated by SLASCA, TLASCA, or STLASCA and the global speckle contrast value can be obtained. In Fig. 5, it is shown that the maximum relative error of ${\mu}_{K}$ computed by STLASCA is approximately 5% when $N$ is 3, while it is higher than 13% when using the other methods. In Fig. 5, the ${\mu}_{K}$ values computed by STLASCA almost coincide with global speckle contrast values, demonstrating that STLASCA utilizes thoroughly the number of pixels that are used for calculating one $K$ value.

The relative noise of the $K$ images computed by the previously mentioned three methods can be quantified by ${\sigma}_{K}\u2215{\mu}_{K}$ . The ${\sigma}_{K}\u2215{\mu}_{K}$ values as a function of number of frames $N$ are shown in Figs. 5 and 5. ${N}_{s}$ values of three and nine pixels are selected for Figs. 5 and 5, respectively. As is shown, ${\sigma}_{K}\u2215{\mu}_{K}$ correlates with both the size of ${N}_{s}$ and $N$ . For a given ${N}_{s}$ value, ${\sigma}_{K}\u2215{\mu}_{K}$ values computed by all of the three methods decrease with increasing $N$ . On the other hand, for a given $N$ value, the ${\sigma}_{K}\u2215{\mu}_{K}$ values decrease with increasing ${N}_{s}$ . It is shown that the ${\sigma}_{K}\u2215{\mu}_{K}$ value computed by TLASCA is slightly higher than the one computed by the others when $N$ is less than 5. This is due to the fact that TLASCA is mainly based on temporal statistics, thus an insufficient number of frames will result in a high ${\sigma}_{K}$ value. With increasing $N$ , the difference in the ${\sigma}_{K}\u2215{\mu}_{K}$ values computed by the three methods becomes insignificant.

The ${\mu}_{K}$ values calculated from each group of 130 $K$ values along vessel center line 1 and vessel center line 2 as a function of number of frames $N$ are shown in Figs. 6 and 6 , respectively. ${N}_{s}$ value of 3 is used for computing $K$ values. As is shown, the ${\mu}_{K}$ values obtained from vessel 1 are lower than those obtained from vessel 2, which suggest that the flow speed of vessel 1 is faster than vessel 2. It is shown in Figs. 6 and 6 that the changes in the ${\mu}_{K}$ values as a function of $N$ computed by these three methods is consistent with the simulation results in Fig. 5. Likewise, the changes in the ${\sigma}_{K}\u2215{\mu}_{K}$ values as a function of $N$ in Figs. 6 and 6 are also consistent with the simulation results in Fig. 5. To save paper, results obtained with other ${N}_{s}$ values are not shown, in which changes in the ${\mu}_{K}$ and the ${\sigma}_{K}\u2215{\mu}_{K}$ values as functions of $N$ are similar to the simulation results.

The
${\sigma}_{K}\u2215{\mu}_{K}$
values computed by STLASCA as a function of number of pixels
${N}_{p}$
used for computing one speckle contrast value are then investigated, where
${N}_{p}={N}_{s}\times {N}_{s}\times N$
. The result is shown in Fig. 7
. The markers represent the
${\sigma}_{K}\u2215{\mu}_{K}$
values obtained by simulation and from vessels 1 and 2. The black solid line draws the expected tendency of
${\sigma}_{K}\u2215{\mu}_{K}$
values changing with
${N}_{p}$
. As is shown, the
${\sigma}_{K}\u2215{\mu}_{K}$
values are obtained from either simulation or experiment, and they scale with the number of pixels as
${N}_{p}^{-0.5}$
, which is consistent with previous reports where the
${\sigma}_{K}\u2215{\mu}_{K}$
value computed by LASCA or SLASCA as a function of
${N}_{p}$
is investigated.^{18, 19}

## 4.

## Discussion

As shown in both Figs. 5 and 6, for the three methods SLASCA, TLASCA, and STLASCA, the
${\sigma}_{K}$
values correlate with both the
${N}_{s}$
and
$N$
values, while the
${\mu}_{K}$
values show different changes with
${N}_{s}$
and
$N$
. For SLASCA, the
${\mu}_{K}$
value is dependent on the
${N}_{s}$
value, but is independent of the
$N$
value. Therefore, SLASCA results in the same
${\mu}_{K}$
value as the one obtained by the basic LSSCA method. For TLASCA, the
${\mu}_{K}$
value is dependent on the
$N$
value, but is independent to the
${N}_{s}$
value. Therefore, TLASCA results in the same
${\mu}_{K}$
value as the one obtained by the basic LSTCA method. For STLASCA, on the contrary, the
${\mu}_{K}$
value correlates with both the
${N}_{s}$
and
$N$
values used, so increase in either the
${N}_{s}$
or
$N$
value results in higher
${\mu}_{K}$
value. The results shown in Figs. 5 and 6 demonstrate that STLASCA utilizes thoroughly the number of pixels
${N}_{p}$
, and thereby provide higher statistical accuracy in the laser speckle contrast value than SLASCA and TLASCA. For comparison in detail, suppose an
${N}_{s}$
value of 7 is used. According to the correlation between the
${\sigma}_{K}\u2215{\mu}_{K}$
value and
${N}_{p}$
, the relative noise is approximately 14% for a single speckle contrast image computed by LSSCA. If five frames of the computed speckle contrast images are averaged, then the relative noise is reduced to approximately 5%, which is more than 2.5 times lower than the relative noise of a single speckle contrast image, demonstrating the validation of SLASCA in reducing image noise. However, as is shown in our numerical simulation result, there remains approximately 5% relative error in
${\mu}_{K}$
for SLASCA. If TLASCA is used, with
${N}_{s}$
and
$N$
values the same as those used for SLASCA, then the relative error in
${\mu}_{K}$
for TLASCA is approximately 7%, which is a little higher than SLASCA. If STLASCA is used, the relative error in
${\mu}_{K}$
can be reduced to be within 0.5%, which is much smaller than those obtained by SLASCA and TLASCA. As has been pointed out, the statistical accuracy of a speckle contrast analysis method is important in LSCI.^{16} The method that achieves maximized
${\mu}_{K}$
means that it maximizes the variation in
$K$
, and thereby maximizes the variation in laser speckle contrast image. Our results demonstrate that STLASCA achieves more statistical accuracy than SLASCA and TLASCA do, given the same
${N}_{s}$
and
$N$
values for these methods. In the practical applications of LSCI, statistical accuracy of speckle contrast, spatial resolution, and temporal resolution are three important measures that should be taken into consideration comprehensively. For practical application of STLASCA, a
$3\times 3\times 15$
or
$5\times 5\times 5\phantom{\rule{0.3em}{0ex}}\text{pixel}$
stack is recommended, resulting in approximately 0.8% relative error in
${\mu}_{K}$
and approximately 8% relative noise in the speckle contrast image.

## 5.

## Conclusion

The statistical accuracy of a laser speckle contrast analysis method for blood flow imaging is important. This study quantitatively compares the statistical accuracy of the currently presented spatiotemporal laser speckle contrast analysis methods through both numerical simulation and *in vivo* laser speckle imaging of rat cortical blood flow. We demonstrate that STLASCA most effectively utilizes the number of pixels used for calculating one
$K$
value, and thereby achieves higher statistical accuracy in the mean speckle contrast than TLASCA and SLASCA do, which is especially meaningful for LSCI being a potentially quantitative tool of estimating blood flow.

## Acknowledgments

This work is supported by the National High Technology Research and Development Program of China (grant number 2007AA02Z303), the PhD Programs Foundation of Ministry of Education of China (grant number 20070487058), the National Natural Science Foundation of China (grant numbers 30711120171, 30970964, and 30801482) and the 111 Program.

## References

**,” Opt. Commun., 37 (5), 326 –330 (1981). https://doi.org/10.1016/0030-4018(81)90428-4 0030-4018 Google Scholar**

*Flow visualization by means of single-exposure speckle photography***,” J. Biomed. Opt., 1 (2), 174 –179 (1996). https://doi.org/10.1117/12.231359 1083-3668 Google Scholar**

*Laser speckle contrast analysis (LASCA): a nonscanning, full-field technique for monitoring capillary blood flow***,” J. Cereb. Blood Flow Metab., 21 (3), 195 –201 (2001). https://doi.org/10.1097/00004647-200103000-00002 0271-678X Google Scholar**

*Dynamic imaging of cerebral blood flow using laser speckle***,” J. Cereb. Blood Flow Metab., 24 518 –525 (2004). https://doi.org/10.1097/00004647-200405000-00005 0271-678X Google Scholar**

*Spatiotemporal quantification of cerebral blood flow during functional activation in rat somatosensory cortex using laser-speckle flowmetry***,” Microvasc. Res., 68 (2), 143 –146 (2004). https://doi.org/10.1016/j.mvr.2004.04.003 0026-2862 Google Scholar**

*Laser speckle imaging for monitoring blood flow dynamics in the**in vivo*rodent dorsal skin fold model**,” J. Innov. Opt. Health Sci., 1 (2), 217 –226 (2008). https://doi.org/10.1142/S1793545808000121 Google Scholar**

*Tracing collateral circulation after ischemia in rat cortex by laser speckle imaging***,” J. Innov. Opt. Health Sci., 1 (2), 239 –256 (2008). https://doi.org/10.1142/S1793545808000236 Google Scholar**

*In vivo*mapping brain microcirculation by laser speckle contrast imaging: a magnetic resonance perspective of theoretical framework**,” J. Biomed. Opt., 8 (3), 559 –564 (2003). https://doi.org/10.1117/1.1578089 1083-3668 Google Scholar**

*Modified laser speckle imaging method with improved spatial resolution***,” Opt. Lett., 31 (12), 1824 –1826 (2006). https://doi.org/10.1364/OL.31.001824 0146-9592 Google Scholar**

*Imaging cerebral blood flow through the intact rat skull with temporal laser speckle imaging***,” IEEE Trans. Med. Imaging, 26 (6), 833 –842 (2007). https://doi.org/10.1109/TMI.2007.892643 0278-0062 Google Scholar**

*New insights into image processing of cortical blood flow monitors using laser speckle imaging***,” Proc. SPIE, 6858 685802 (2008). https://doi.org/10.1117/12.760514 0277-786X Google Scholar**

*Spatio-temporal algorithms for processing laser speckle imaging data***,” Opt. Express, 16 14321 –14329 (2008). https://doi.org/10.1364/OE.16.014321 1094-4087 Google Scholar**

*Fast blood flow visualization of high-resolution laser speckle imaging data using graphics processing unit***,” J. Opt. Soc. Am. A, 25 231 –237 (2008). https://doi.org/10.1364/JOSAA.25.000231 0740-3232 Google Scholar**

*The copula: a tool for simulating speckle dynamics***,” Opt. Express, 16 (5), 3197 –3203 (2008). https://doi.org/10.1364/OE.16.003197 1094-4087 Google Scholar**

*Impact of velocity distribution assumption on simplified laser speckle imaging equation***,” Opt. Lett., 33 (24), 2886 –2888 (2008). https://doi.org/10.1364/OL.33.002886 0146-9592 Google Scholar**

*Detrimental effects of speckle-pixel size matching in laser speckle contrast imaging***,” EPL, 82 (1), 18005 (2008). https://doi.org/10.1209/0295-5075/82/18005 0295-5075 Google Scholar**

*LASCA with a small number of scatterers: application for monitoring of microflow***,” Opt. Express, 13 (24), 9782 –9787 (2005). https://doi.org/10.1364/OPEX.13.009782 1094-4087 Google Scholar**

*Laser speckle imaging with an active noise reduction scheme***,” J. Opt. Soc. Am. A, 25 (1), 9 –15 (2008). https://doi.org/10.1364/JOSAA.25.000009 0740-3232 Google Scholar**

*Statistics of local speckle contrast*