<figure><img src="image_1.png" /></figure>

HAL
open science

# Choice of numerical implementation of spatial contrast calculation impacts microcirculation quantitation in laser speckle contrast imaging

Marc Chammas, Frédéric Pain

<figure><img src="image_2.png" /></figure>

## ► To cite this version:

Marc Chammas, Frédéric Pain. Choice of numerical implementation of spatial contrast calculation impacts microcirculation quantitation in laser speckle contrast imaging. Journal of Biomedical Optics, 2025, 30 (4), pp.046006. 10.1117/1.JBO.30.4.046006. hal-05037804

HAL Id: hal-05037804

https://hal.science/hal-05037804v1

Submitted on 17 Apr 2025

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_3.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Choice of numerical implementation of spatial contrast calculation impacts microcirculation quantitation in laser speckle contrast imaging

Marc Chammas and Frédéric Pain*

Université Paris-Saclay, Institut d'Optique Graduate School, CNRS, Laboratoire Charles Fabry, Palaiseau, France

**ABSTRACT.** **Significance:** Laser speckle contrast imaging (LSCI) allows noninvasive imaging of microcirculation. Its scope of clinical applications is growing, yet the literature lacks a comparison of the accuracy of methods used to compute the spatial contrast $K_s$ from which the blood flow index is derived.

**Aim:** We aim to evaluate the impact on flow quantitation of different computational approaches used to derive $K_s$.

**Approach:** We compare numerical calculation of $K_s$ in Python and ImageJ applied to noise-free simulated data and to experimental data acquired *in vivo* in anesthetized mice. The estimation of the decorrelation time $\tau_c$, inversely proportional to the blood flow index, is carried out following two approaches: LSCI asymptotic estimation and fitting the multiple exposure speckle imaging (MESI) model to $K_s(T)$.

**Results:** For simulation data, we found variations of up to 58% for the blood flow index in the LSCI approach. Nonlinear fitting of the MESI model was less affected with discrepancies of only a few percent. Considering experimental data, the LSCI approximation led to $K_s$ with relative differences (up to 35%) depending on the calculation methods. The noise and limited exposure time strongly limited the accuracy of the LSCI asymptotic estimation. Adjustment of the MESI model to the data led to consistent values of $\tau_c$ in the 0.05 to 1 ms range with significant variations depending on the method used to calculate $K_s$.

**Conclusions:** Numerical methods used to calculate $K_s$ should be precisely acknowledged and validated against direct calculation to ensure accuracy. *Uniform* filter approach leads to accurate $K_s$ values and is 100 times more computationally efficient than the *Direct* calculation. Other investigated methods lead to various levels of errors in flow index estimation using LSCI. Errors are minimized using larger kernels. MESI derivation of $\tau_c$ is not immune but less affected by such methodological biases.

© The Authors. Published by SPIE under a Creative Commons Attribution 4.0 International License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI. [DOI: 10.1117/1.JBO.30.4.046006]

**Keywords:** laser speckle contrast imaging; spatial contrast; 2D filters; convolution; blood flow

Paper 240335GR received Nov. 18, 2024; revised Feb. 12, 2025; accepted Mar. 25, 2025; published Apr. 16, 2025.

## 1 Introduction

Laser speckle contrast imaging (LSCI) relies on the illumination of a sample with coherent light from a laser. Speckles are interferometric patterns that are observed when coherent light is used to

<footnote>*Address all correspondence to Frédéric Pain, frederic.pain@universite-paris-saclay.fr</footnote>
illuminate objects. When a sample contains moving scatterers such as red blood cells in a blood vessel, the pattern evolves over time. Recording these time-varying speckle patterns with a camera for a given exposure time leads to a local blur in the image. The local statistics of the resulting grainy speckle pattern contain information on the movements of the scatterers within the sample. The most common metric used is the speckle contrast, defined as the local ratio of the standard deviation to the mean of the gray levels of the pixels $K_s = \sigma/\langle I \rangle$. It can be computed in spatial, temporal, or spatio-temporal ensembles of neighboring pixels.¹ For computing efficiency, the images are processed numerically using sliding convolution filters. In most implementations, the spatial speckle contrast is calculated over windows of 5 × 5 or 7 × 7 pixels. This allows a reasonable spatial resolution and adequate statistical accuracy.² The speckle contrast values quantify the blur of the local speckle pattern, which is in turn related to the local decorrelation time $\tau_c$ and ultimately to the speed of the scatterers movements.³ As the exact relationship between $\tau_c$ and the speed is not known, a relative metric is used, the blood flow index (BFI), estimated as $1/\tau_c$. In LSCI, the BFI can be derived as $\text{BFI} = 1/\text{TK}^2$ if $T/\tau_c > 100$.⁴ Multiple exposure speckle imaging (MESI) is an improved variant of LSCI, where several values of $K_s$ are obtained at different exposure times. A mathematical model can then be fitted to the resulting $K_s(T)$ providing a more robust estimation of the local decorrelation time $\tau_c$.⁵ The relative blood flow changes are then estimated through the relative changes of the inverse correlation time (relative ICT) defined as $1/\tau_c$.⁶ Regarding quantitation and comparison of measurements between systems, several instrumental and methodological parameters have been pointed out and investigated throughout the years, including the effect of the speckle to pixel size ratio, the contribution of specular reflected light, and the amount of system noises.¹,⁷ Although the accuracy of the metrics detailed above depends on the accuracy of $K_s$ calculation, the effect of the numerical methods implemented to calculate the speckle contrast was less investigated. A computationally efficient approach was proposed based on the calculation of the variance as the mean of squared local intensities values minus the square of the mean of these values.⁸ However, the actual numerical implementation of $K_s$ calculation is seldom provided in the literature. Studying the accuracy of $K_s$ numerical calculations and its potential impact on BFI or ICT is relevant in light of the current efforts toward improved standardization and relative and absolute quantitation.⁷,⁹⁻¹² The aim of this work is to compare the accuracy of different numerical approaches for spatial contrast calculations and the corresponding $\tau_c$ evaluations. The computational speed and accuracy in both LSCI and MESI implementations for imaging of microcirculation have been investigated on simulated and experimental raw speckle data.

# 2 Methods

## 2.1 Speckle Images Datasets

### 2.1.1 *Simulated data*

Dataset #1 is composed of simulated images of time-integrated speckle patterns at 15 different exposure times. The simulated camera exposure times are 0.05, 0.1, 0.2, 0.5, 0.8, 1, 2, 3.5, 5, 10, 15, 20, 30, 40, and 60 ms. These exposure times span two orders of magnitude allowing to sample the decrease of the speckle contrast with increasing exposure time in MESI studies of mice brain microcirculation.¹³ Images have a bit depth of 12 bits and dimensions of 1000 × 1000 pixels. The simulations were carried out following the implementation described by James et al.¹⁴ The method allows the simulation of time-integrated speckle patterns for custom adjustable decorrelation functions and parameters. The geometry is a homogeneous, purely diffusive media with no static scatterers. A decorrelation function corresponding to multiple photon scattering in ordered flow was simulated. Other flow and scattering regimes can be simulated, but we chose specifically this one as *in vivo* measurements of $g_1(\tau_c)$ have shown that it is relevant for arterioles in functional areas of the mice brain.¹⁵ Two sets of data with decorrelation times of 0.1 and 1 ms were simulated in the range of $\tau_c$ values that were observed experimentally in mice brain arterioles.¹⁵ A ratio of 3 pixels per speckle was set in the simulation to ensure speckle grain sampling was better than the Shannon–Nyquist criteria.²
### 2.1.2 *Experimental in vivo data*

Dataset #2 is composed of experimental raw MESI speckle images acquired *in vivo* at the surface of exposed cortical tissues of a mouse. These data have been obtained previously.¹⁶ In this previous study, the number of exposure times was limited to 6 with a maximum of 20 ms to reduce the duration of an MESI sequence and improve the overall temporal resolution of the technique. *In vivo* imaging was carried out in accordance with the European Directive 2010/63/UE (Ethical project #CEB-03-2018) regarding the care and use of laboratory animals. *In vivo* recordings were carried out in 4-month-old, anesthetized mice (male, C57BL6, mean weight = 48 g) after a craniotomy was performed over the barrel cortex. Anesthesia was carried out with a cocktail of ketamine (100 mg/kg Imalgen) and medetomidine (0.5 mg/kg), injected intraperitoneally. The anesthesia level was adjusted as necessary throughout the experiment. Surgery begins when the mouse no longer responds to the hindpaw pinch. During the entire experiment, the animal was placed on a thermostatically controlled heating pad. Body temperature was maintained rectally at 37°C. After the animal was secured in a stereotactic frame, imaging was carried out over regions of ~3 × 2 mm² over the barrel cortex. A low-temperature agarose layer (1.2%) and a microscope coverslip were placed over the dura matter. A stabilized laser diode (SGL, Shanghai, China, 634 nm 200 mW laser diode) in free space was modulated in time and intensity using an acousto-optic modulator (AOM) (AA Opto-Electronic, Orsay, France) driven by an arbitrary waveform generator (Keysight 336001 Series, Penang, Malaysia). The AOM allows modulating the first diffraction order of the signal. A fixed aperture diaphragm was used to select the first order (modulated) of diffraction while filtering the non-diffracted light as well as eventual higher diffraction orders. After 15 min of laser warm-up, light was shone on the exposed cortical tissue at an oblique incidence to minimize specular reflection contributions. Cycling of the modulated laser pulses is synchronized with image acquisition. The camera was a CMOS Orca Flash 2.8 (Hamamatsu Photonics, Massy, France) with a frame rate of 46 fps at full frame. It was c-mounted on a MZ16 microscope (Leica Microsystems, Nanterre, France) providing a × 3.3 magnification resulting in a field of view of ~3 × 2 mm². A ratio of 2.3 pixels per speckle grain was measured experimentally. Dataset #2 is composed of 774 frames cycling through 6 exposure times (1, 2, 5, 10, 15, and 20 ms) corresponding to 129 MESI sequences.

## 2.2 Numerical Implementations of Speckle Contrast Calculations

The Python and ImageJ implementations of $K_s$ calculations are summarized in Table 1. Python was chosen as it is open source and has become one of the common languages for teaching and research in the biomedical engineering community. Versions of software and Python modules are summarized in Table 2. Data were processed using a MacBook Air with an M1 microchip and 16 GB of random access memory (RAM). The reference calculation of $K_s$ is referred to as the *Direct* method. It was implemented using nested “for” loops based on the mathematical definitions of the mean and variance. Other methods rely on vectorized approaches using the “sums algorithm”⁸ and the calculation of local squared, summed, and averaged values using convolution filters. The corresponding Python and ImageJ scripts are available in a Zenodo repository.¹⁷ For both ImageJ macro and Scipy uniform_filter Python function, the data are processed as 3D stacks, which is convenient for LSCI and MESI data consisting of stacks (i.e., movies) of 2D image frames. On the contrary, Scipy correlated2d (adapted from Ref. 18) and SimpleITK convolution methods require 2D arrays as inputs, which implies the use of “for” loops

<table><caption>**Table 1** Numerical methods implemented to calculate $K_s$.</caption><thead><tr><td>Name</td><td>Language</td><td>Libraries/modules</td><td>Input</td><td>Parameters</td></tr></thead><tbody><tr><td>Uniform_filter</td><td>Python</td><td>Scipy.ndimage</td><td>Images stack</td><td>2D or 3D kernel</td></tr><tr><td>Correlated2d</td><td>Python</td><td>Scipy.signal</td><td>Single 2D image</td><td>2D kernel</td></tr><tr><td>SimpleITK</td><td>Python</td><td>SimpleITK</td><td>Single 2D image</td><td>2D circular kernel</td></tr><tr><td>ImageJ_Macro</td><td>FIJI</td><td>Built-in Filters</td><td>Images stack</td><td>2D circular kernel</td></tr></tbody></table>
<table><caption>**Table 2** Versions of software and Python modules.</caption><thead><tr><th>Python</th><th>Scipy</th><th>Numpy</th><th>SimpleITK</th><th>ImageJ</th></tr></thead><tbody><tr><td>3.10.14</td><td>1.14.0</td><td>2.0.0</td><td>2.3.1</td><td>2.14.0 1.54f</td></tr></tbody></table>

to process 3D stacks corresponding to frame sequences. By default, ImageJ uses circular convolution kernels, whereas *Uniform_filter* and *Scipy correlated2d* use square kernels. *SimpleITK* convolution is implemented with 2D circular kernels that are defined in a custom separate function. Square kernels have a side dimension of 5 or 7 pixels. For circular kernels, radii of 2 or 3 pixels are set, resulting in kernels of 5 or 7 pixel diameter. The default edge management for the different methods is self-explanatory and named “reflect,” “constant padding,” and “zero padding” for *Uniform_filter*, *Correlated2d*, and *SimpleITK* convolutions, respectively. They can be customized to some extent for *Uniform_filter* and *Correlated2d* to “mirror” or “symmetrical” values around the edges. The edge management affects only the values on the sides of the images. Consequently, to avoid eventual bias, we have cropped the borders of the final speckle contrast images to remove the affected edge pixels. For each dataset, the corresponding $K_s$ images are calculated using the different methods listed in Table 1. For dataset #1, as the imaging field is homogeneous, the resulting mean $K_s$ is calculated as the average of all pixels for each contrast frame. For dataset #2, the spatial contrast frames are calculated for all 774 frames, and then, for each exposure time, the corresponding 129 frames are averaged, providing an average MESI sequence with six exposure times. Then, three regions of interest (ROIs) are drawn, over two arterioles and a parenchyma region [see Fig. 2(a)]. Care was taken not to include the vessel walls in the arterial ROIs. The speckle contrast values within each ROI are averaged to provide the mean $K_s$ as a function of the exposure times.

## 2.3 Evaluation of Decorrelation Time

BFI is defined as $1/\tau_c$. For LSCI data, it is estimated as $BFI = 1/TK^2$, following the asymptotic approximation.⁴ For MESI, $\tau_c$ values are estimated using a nonlinear regression of the MESI model on the data⁵,¹⁹

$$K(T, \tau_c) = \left\{ \beta \rho^2 \left( \frac{e^{-2x} - 1 + 2x}{(2x)^2} \right) + 4\beta\rho(1-\rho) \left( \frac{e^{-x} - 1 + x}{x^2} \right) + v_{ne} + v_{noise} \right\}^{\frac{1}{2}}, \quad (1)$$

$T$ is the exposure time, $\beta$ is a constant depending on instrumental parameters, $\rho$ is the proportion of dynamic scatterers involved in the speckle patterns, $x$ is $T/\tau_c$, $v_{ne}$ accounts for nonergodic contribution, and $v_{noise}$ accounts for noise contributions. To fit the MESI model to the dataset #1, the parameter $\rho$ was set to 1, $v_{noise}$ and $v_{ne}$ were set to 0 corresponding to the simulation’s “ideal” conditions, i.e., an imaging field composed only of dynamic scatterers, and no experimental noises. For experimental dataset #2, the nonergodic and experimental noises’ contributions were pooled in a single parameter $v_{noise}$. The parameters were constrained within the following physical boundaries: $0.05 < \beta < 0.5$, $0.1 < \rho < 1$, $0 < \tau_c < 10$, $0 < v_{noise} < 0.30$.

# 3 Results and Discussion

The local spatial speckle contrast for all frames in datasets #1 was computed using the numerical methods from Table 1. It corresponds to the time requested to calculate the contrast images for one MESI sequence of 1000 × 1000 pixels images with 15 different exposure times. We have first evaluated the computation times for each method and then the influence on the contrast values, resulting $\tau_c$ in LSCI and MESI.

## 3.1 Computing Time

Execution times were evaluated in Python using the time module. ImageJ macro was timed in batch mode (i.e., no screen output occurred during the calculation). For each Python approach listed in Table 1, an average timing was evaluated based on 100 iterations of the same calculation to smooth out system overheads. The absolute execution times are only indicative as they highly depend on the hardware used for calculation. The relative execution times are more relevant.
<table><caption>**Table 3** Execution times for different implementations of speckle calculation on dataset #1 ($\tau_c$ = 0.1 ms).</caption><thead><tr><th rowspan="2">Kernel dimension (pixels)</th><th colspan="5">Execution times (s)</th></tr><tr><th>Direct</th><th>Uniform_filter</th><th>Correlated2d</th><th>SimpleITK</th><th>ImageJ</th></tr></thead><tbody><tr><td>5</td><td>410</td><td>0.2</td><td>1.6</td><td>0.2</td><td>0.4</td></tr><tr><td>7</td><td>407</td><td>0.2</td><td>3.0</td><td>0.4</td><td>0.4</td></tr></tbody></table>

The results in Table 3 indicate that the direct estimation is ~100 to 1000 times slower than the other methods. Overall, *simpleITK* and *Uniform_filter* are the fastest. *Correlated2d* implementation is 8 to 15 times slower depending on the kernel dimension. *ImageJ* timing is twice that of *uniform*. Methods were applied using kernels of 5 or 7 pixel side or radius. Except for *Scipy_correlated2D* and *SimpleITK*, the execution times of which were almost doubled, the other timings were not significantly affected by the kernel dimensions.

## 3.2 Speckle Contrast Values for Simulated Data

All the different methods were applied to calculate the mean contrast for each frame of dataset #1. The resulting $K_s$ for kernels with a 5 or 7 pixel side or diameter are shown in Figs. 1(a) and 1(b) for the data with a simulated decorrelation time $\tau_c$ = 0.1 ms. Considering the *Direct* calculation as the reference, the relative errors for each method are presented in Figs. 1(c) and 1(d). The uniform method provides values that are almost identical to those of the *Direct* calculation. *ImageJ* calculation also presents values with less than 5% difference with the reference. *SimpleITK* and finally *Correlated2d* give the values with the strongest discrepancies, up to 20% for *Correlated2d* with a 5 pixel window. As expected, the calculations using larger kernels lead to lower discrepancies as the statistics of the speckle patterns are better sampled.² Similar trends are observed for data with a simulated decorrelation time $\tau_c$ = 1 ms (Fig. S1 in the Supplementary Material).

<figure><img src="image_4.png" /><figcaption>**Fig. 1** $K_s(T)$ derived for the different methods and relative errors for simulated data at $\tau_c$ = 0.1 ms. *Direct* calculation is considered the reference. (a) $K_s(T)$ for kernels of 5 × 5; (b) $K_s(T)$ for kernels of 7 × 7; (c) relative errors (kernels of 5 × 5); and (d) relative errors (kernels of 7 × 7).</figcaption></figure>
### 3.3 Estimation of the Decorrelation Times for Simulated Data

Table 4 shows the estimations of $\tau_c$ derived from the LSCI asymptotic approximation and from the nonlinear regression of the MESI model on $K_s(T)$ derived using all methods applied to dataset #1. For each exposure time, values of $K_s$ are calculated as the average of all pixels. In the LSCI asymptotic approximation, regardless of the method used to compute $K_s$, the estimated $\tau_c$ is closer to the expected simulated value when $T/\tau_c$ increases. As expected, the LSCI estimations for the shorter exposure times are poor in all cases. LSCI estimations at the longest exposure time (60 ms) are fair for the data simulated with $\tau_c = 0.1$ ms, but poor for the data simulated with $\tau_c = 1$ ms. Comparing the methods, for LSCI, the lowest accuracy is observed for Correlated2d in line with the discrepancies observed for $K_s$ calculation [Figs. 1(c) and 1(d)]. On the contrary, $\tau_c$ values derived from nonlinear regressions of the MESI model are less affected by the numerical method used. Considering the analysis with 7 × 7 kernels, all the MESI estimates of $\tau_c$ are

<table><caption>**Table 4** Estimated values of the decorrelation time $\tau_c$ (in ms) for dataset #1 obtained using the LSCI asymptotic approximation or by fitting the MESI model to the $K_s(T)$ derived using each numerical method and kernels of side or diameter of 5 and 7 pixels.</caption><thead><tr><th>Asymptotic approximation</th><th>Direct</th><th>Uniform_filter</th><th>Correlated2d</th><th>SimpleTK</th><th>ImageJ</th></tr></thead><tbody><tr><td colspan="6">Simulated τ<sub>c</sub> = 0.1 ms, kernel side or diameter of 5 pixels</td></tr><tr><td>Texp = 1 ms</td><td>0.060</td><td>0.060</td><td>0.040</td><td>0.047</td><td>0.056</td></tr><tr><td>Texp = 20 ms</td><td>0.069</td><td>0.069</td><td>0.045</td><td>0.054</td><td>0.065</td></tr><tr><td>Texp = 40 ms</td><td>0.079</td><td>0.079</td><td>0.052</td><td>0.062</td><td>0.074</td></tr><tr><td>Texp = 60 ms</td><td>0.095</td><td>0.095</td><td>0.061</td><td>0.074</td><td>0.091</td></tr><tr><td>MESI regression</td><td>0.106</td><td>0.106</td><td>0.091</td><td>0.101</td><td>0.105</td></tr><tr><td colspan="6">Simulated τ<sub>c</sub> = 0.1 ms, kernel side or diameter of 7 pixels</td></tr><tr><td>Texp = 1 ms</td><td>0.072</td><td>0.072</td><td>0.059</td><td>0.063</td><td>0.067</td></tr><tr><td>Texp = 20 ms</td><td>0.082</td><td>0.082</td><td>0.069</td><td>0.073</td><td>0.077</td></tr><tr><td>Texp = 40 ms</td><td>0.093</td><td>0.093</td><td>0.081</td><td>0.083</td><td>0.088</td></tr><tr><td>Texp = 60 ms</td><td>0.112</td><td>0.112</td><td>0.092</td><td>0.100</td><td>0.106</td></tr><tr><td>MESI regression</td><td>0.105</td><td>0.105</td><td>0.101</td><td>0.105</td><td>0.106</td></tr><tr><td colspan="6">Simulated τ<sub>c</sub> = 1 ms, kernel side or diameter of 5 pixels</td></tr><tr><td>Texp = 1 ms</td><td>0.321</td><td>0.321</td><td>0.228</td><td>0.252</td><td>0.299</td></tr><tr><td>Texp = 20 ms</td><td>0.626</td><td>0.626</td><td>0.407</td><td>0.489</td><td>0.585</td></tr><tr><td>Texp = 40 ms</td><td>0.646</td><td>0.646</td><td>0.417</td><td>0.504</td><td>0.605</td></tr><tr><td>Texp = 60 ms</td><td>0.655</td><td>0.655</td><td>0.425</td><td>0.511</td><td>0.612</td></tr><tr><td>MESI regression</td><td>0.935</td><td>0.935</td><td>0.838</td><td>0.858</td><td>0.918</td></tr><tr><td colspan="6">Simulated τ<sub>c</sub> = 1 ms, kernel side or diameter of 7 pixels</td></tr><tr><td>Texp = 1 ms</td><td>0.394</td><td>0.394</td><td>0.333</td><td>0.342</td><td>0.364</td></tr><tr><td>Texp = 20 ms</td><td>0.743</td><td>0.743</td><td>0.618</td><td>0.660</td><td>0.699</td></tr><tr><td>Texp = 40 ms</td><td>0.765</td><td>0.765</td><td>0.636</td><td>0.680</td><td>0.718</td></tr><tr><td>Texp = 60 ms</td><td>0.774</td><td>0.774</td><td>0.675</td><td>0.689</td><td>0.726</td></tr><tr><td>MESI regression</td><td>0.976</td><td>0.976</td><td>0.995</td><td>0.944</td><td>0.995</td></tr></tbody></table>
within a few percentages from the expected value. The values of $τ_c$ derived from the $K_s$ calculated with the 5 × 5 kernels are less accurate, with up to 17% discrepancy for the simulated data at $τ_c$ = 1 ms analyzed with *Scipy correlated2d*. As the *ICT* is computed as the inverse of $τ_c$, it is biased by the same relative error as $τ_c$. Regarding the kernel shape, most speckle imaging studies use square kernels as they are easy to implement in array processing. However, square kernels can introduce directional bias, whereas round kernels treat data isotropically. Nevertheless, the results show that the kernel shape is less critical than the numerical method and the dimension of the kernel used. Despite the use of kernels with different shapes, *Uniform* (square kernel) and *ImageJ* (round kernel) provide the best results considering the *Direct* method as the reference. All methods have better accuracy when larger kernels are used, regardless of the shape of the kernel.

## 3.4 Limitations of Simulation Data

Simulations from dataset #1 do not account for experimental noises that may affect significantly the speckle contrast values. These include camera noise (dark, read-out, fixed noise pattern, and shot noise) as well as potential instability of the illumination, contamination with noncoherent ambient light, or unwanted movements (vibrations or respiratory movements). Most of these contributions can be significantly reduced through careful design of the experimental set-up. However, they cannot be completely avoided. In addition, simulations in dataset #1 represent a purely dynamic, homogeneous media, whereas *in vivo* data are a mix of dynamic and static contributions. In Sec. 3.5, we have further evaluated the impact of the choice of the numerical implementation of $K_s$ calculation for real data acquired *in vivo* in mice.

## 3.5 Speckle Contrast Calculation and Decorrelation Times Estimation on *In Vivo* Data

Speckle contrast images derived from real data present inherently a nonhomogeneous field of view. Figure 2(a) shows, for all numerical methods, the speckle contrast images obtained for an exposure time of 1 ms and the corresponding relative differences with the contrast $K_s$ images obtained with the *Direct* method. Figure 2(b) shows the speckle contrast images and corresponding relative error maps for an exposure time of 20 ms. Although $K_s$ images look very similar for all methods, the relative error images show significant local differences in $K_s$ compared with the

<figure><img src="image_5.png" /><figcaption>**Fig. 2** (a) $K_s$ maps for all methods (kernel side or diameter of 7 pixels) at 1 ms exposure time, averaged over 129 MESI sequences (top row) and relative error maps compared with the direct method (bottom row). Three regions of interest are drawn over a large arteriole (ROI#1), a small arteriole (ROI#2) and parenchyma (ROI#3); (b) same as (a) for 20 ms exposure time.</figcaption></figure>
Direct calculation. Here, for the purpose of comparison, the same scale was applied to the relative error images. Figures S2 and S3 in the Supplementary Material show the $K_s$ and relative error images at all exposure times for all methods. The largest relative errors are observed for the *Correlated2D* (–25%) and *SimpleITK* (–20%) methods similar to what was observed for simulated data. *Uniform* methods provide almost error-free $K_s$ images, whereas *ImageJ* leads to intermediate errors in the ±5% range. The relative errors are relatively homogeneous inside the larger vessels, the smaller vessels, and parenchyma. For larger vessels, the gradients of integrated speckle values lead to larger errors in the wall regions, as observed in Fig. 2(b), particularly for *Correlated2d*. This bias is more pronounced at longer exposure times and can also be observed for *SimpleITK* and *ImageJ* for an exposure time of 20 ms. For all methods, the mean contrast $K_s(T)$ was calculated over the three ROIs, as shown in Fig. 2(a). The results are plotted in Figs. 3 and 4 as a function of exposure time for the different calculation methods and kernels of 5 × 5 and 7 × 7 pixels, respectively.

The absolute values of contrasts are much lower in the experimental data compared with the simulated dataset due to fast decorrelation and added experimental noises. Consistent with expectations, the lowest values are observed in the region with presumably the fastest flow (large arteriole, ROI1), and the highest $K_s$ are observed in the region with the slowest flow (ROI3). The small arteriole in ROI2 shows intermediate values. For kernels of 5 pixel side or diameter, the ranking and trends in terms of relative accuracy on $K_s(T)$ calculations are similar to those observed for calculations on simulated data [compare Figs. 1(a), 1(b), and 3]. The values of $K_s$ obtained with *Uniform* are almost identical to those obtained with the *Direct* calculation. Regarding other methods, *ImageJ* analysis remains close to the reference value followed by *SimpleITK* convolution and *Correlated2d*. For these two implementations, the relative errors are much smaller compared with those obtained on simulated data. This could be attributed to the

<figure><img src="image_6.png" /><figcaption>**Fig. 3** Evaluation of $K_s$ and relative errors for dataset #2 using all methods and kernels of side or diameter of 5 pixels. (a) $K_s$ values as a function of exposure times for all methods and kernels of 5 pixel side or diameter. (b) Relative errors for the $K_s$ values relative to the *Direct* method as a function of exposure times for all methods. Each row corresponds to one of the ROIs drawn in Fig. 2(a).</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Fig. 4** Evaluation of $K_s$ and relative errors for dataset #2 using all methods and kernels of side or diameter of 7 pixels. (a) $K_s$ values as a function of exposure times for all methods and kernels of 7 pixel side or diameter. (b) Relative errors for the $K_s$ values relative to the *direct* method as a function of exposure times for all methods. Each row corresponds to one of the ROIs drawn in Fig. 2(a).</figcaption></figure>

averaging of 129 MESI sequences for the *in vivo* dataset, whereas simulated data are single frames for each exposure time. In addition, the errors decrease slightly with exposure times for these methods. Indeed, the relative errors observed for *Scipycorrelated2d* are twice higher for ROI3 (parenchyma) compared with those for the large arteriole region (ROI1). This follows the values of the $K_s$ in these two regions, which are about twice higher in ROI3 than in ROI1. Regarding the influence of the kernel dimensions, as observed for the simulated data, using larger kernels leads to significantly lower errors (compare Figs. 3 and 4) notably for *Scipycorrelated2d* and *SimpleITK*. This can again be attributed to improved estimations of speckle pattern statistics. For 7 × 7 kernels and for all methods, the relative errors remain constant over the different exposures [Fig. 4(a)]. In the LSCI asymptotic approximation, $\tau_c$ is estimated directly from $1/K_s^2$, so in a first-order approximation, the relative uncertainty on $\tau_c$ is twice that on $K_s$. Considering the values of relative errors on $K_s$ for the maximum exposure time of 20 ms in Fig. 3, this theoretically leads to deviations on $\tau_c$ between 5% and 35% for kernels of 5 × 5 pixels. For kernels of 7 × 7 pixels, the deviations are limited to 5% to 12%. Finally, the graphs of $K_s(T)$ in Figs. 3 and 4 show that the maximal exposure time considered in this dataset does not allow reaching a clear asymptotic value of $K_s$, thus preventing an accurate estimation of $\tau_c$ using the LSCI approximation.

Table 5 presents the values of $\tau_c$ derived from the MESI model for all methods and mean $K_s(T)$ from the three ROIs, as shown in Fig. 2(a). The values obtained for kernels of 7 and 5 pixel side or diameter are shown. The values of $\tau_c$ estimated for both kernel dimensions show similar trends, but the $\tau_c$ values are systematically lower for kernels of 5 pixel side or diameter. Considering $K_s(T)$ obtained by the *Direct* and *Uniform_filter* methods using kernels of 7 × 7, the values of $\tau_c$ derived from the MESI model are very close and range between 0.499 ms for the parenchyma (ROI3) region down to 0.135 ms for the small arteriole (ROI2) and 0.094 ms in the largest arteriole (ROI1), in line with what was quantified *in vivo*.
<table><caption>**Table 5** Estimated values of the decorrelation time $\tau_c$ (in ms) obtained for each ROI in Fig. 2(a) by fitting the MESI model on $K_s(T)$ calculated with all the different numerical implementations with kernel side or diameter of 5 or 7 pixels.</caption><thead><tr><th>ROI</th><th>Arteriole diameter</th><th>Direct</th><th>Uniform_filter</th><th>SimpleITK</th><th>Correlated2D</th><th>ImageJ</th></tr></thead><tbody><tr><td colspan="7">Kernel side or diameter of 5 pixels</td></tr><tr><td>ROI1</td><td>178 ± 10 µm</td><td>0.076</td><td>0.076</td><td>0.061</td><td>0.064</td><td>0.074</td></tr><tr><td>ROI2</td><td>38 ± 10 µm</td><td>0.115</td><td>0.115</td><td>0.085</td><td>0.09</td><td>0.108</td></tr><tr><td>ROI3</td><td>Parenchyma</td><td>0.342</td><td>0.342</td><td>0.204</td><td>0.228</td><td>0.299</td></tr><tr><td colspan="7">Kernel side or diameter of 7 pixels</td></tr><tr><td>ROI1</td><td>178 ± 10 µm</td><td>0.094</td><td>0.094</td><td>0.083</td><td>0.082</td><td>0.087</td></tr><tr><td>ROI2</td><td>38 ± 10 µm</td><td>0.135</td><td>0.135</td><td>0.123</td><td>0.122</td><td>0.129</td></tr><tr><td>ROI3</td><td>Parenchyma</td><td>0.499</td><td>0.501</td><td>0.388</td><td>0.379</td><td>0.429</td></tr></tbody></table>

in arterioles in the mouse brain using direct decorrelation function measurement.¹⁵ Compared with these values, $\tau_c$ values derived from $K_s$ calculated using *Correlated2d*, *SimpleITK*, and *ImageJ* are systematically underestimated for all ROIs. Underestimations of $\tau_c$ are similar for *SimpleITK* and *Scipycorrelated2d* and significantly lower for *ImageJ*.

## 4 Conclusions

The study shows the strong impact on $\tau_c$ blood flow index quantitation of the numerical methods used to calculate the spatial contrast. The observed significant underestimations of $K_s$ can be partly minimized if larger kernels dimensions are used and if contrast frames' averaging is carried out. Among the methods tested, the *Uniform_filter* with a 2D square filter gave results almost identical to the *Direct* estimation but for a calculation time several orders of magnitude shorter. *ImageJ* circular filter implementation also provided fair results, within a few percentages of the *Direct* calculation. Other implementations (*Scipy correlated2d*, *SimpleTK*) gave much larger errors' event for large kernels and substantial contrast frames' averaging, which should be avoided. For the experimental data, the presence of noises and the lower values of contrast lead to a limited accuracy of the asymptotic LSCI approximation at a single exposure time. For all implementations, MESI analysis leads to physiologically consistent flow index estimation with $\tau_c$ decreasing with the increasing diameter of arterioles (associated to faster flow) with $\tau_c$ values in the same order of magnitude as absolute quantitation obtained in mice brain.¹⁵ The MESI model is fitted to the experimental curves of $K_s$ as a function of time, the shape of which is preserved by the different algorithms. This minimizes the impact of $K_s$ underestimations by the numerical methods. However, *ImageJ*, *SimpleITK*, and *Correlated2D* methods all led to significant underestimations of $\tau_c$. In conclusion, although the choice of the software implementation for $K_s$ calculation may seem trivial; the present data show that this often overlooked point definitely plays a role in the blood flow index quantitation in LSCI and MESI imaging. The Python method based on *Scipy Uniform* filter provides accurate calculations of $K_s$ calculations and subsequent estimation of $\tau_c$. It was also one of the fastest in terms of computation efficiency. In light of the continuous theoretical and experimental efforts and recent advances toward better relative and even absolute quantitation of blood flow with laser speckle imaging,⁷,⁹,¹¹,¹²,²⁰,²¹ we believe that the numerical implementation used in LSCI or MESI to calculate the speckle contrast should be systematically compared against the *Direct* calculation to avoid numerical bias.

## Disclosures

The authors have no relevant financial interests and no other potential conflicts of interest to disclose.
**Code and Data Availability**

The data and code supporting this study are available in the Zenodo open repository at https://doi.org/10.5281/zenodo.14178966.

**Acknowledgments**

Marc Chammas was funded as a PhD student by a scholarship from the French Ministère de l'Enseignement Supérieur et de la Recherche. He was affiliated with Ecole Doctorale Electrical, Optical, Bio-Physics and Engineering (EOBE, ED no. 575). We thank Dr. Hsin-Hon Lin from Chang Gung University in Taiwan for insightful discussions on simulations of time-integrated speckle patterns. We thank Dr. Hirac Gurden from Université Paris Cité, Laboratoire BFA, for help with small animal in vivo imaging experiments.

**References**

1. P. G. Vaz et al., "Laser speckle imaging to monitor microvascular blood flow: a review," *IEEE Rev. Biomed. Eng.* **9**, 106–120 (2016).

2. O. Thompson, M. Andrews, and E. Hirst, "Correction for spatial averaging in laser speckle contrast analysis," *Biomed. Opt. Express* **2**, 1021–1029 (2011).

3. J. D. Briers and S. Webster, "Laser speckle contrast analysis (LASCA): a nonscanning, full-field technique for monitoring capillary blood flow," *J. Biomed. Opt.* **1**, 174–179 (1996).

4. H. Cheng and T. Q. Duong, "Simplified laser-speckle-imaging analysis method and its application to retinal blood flow imaging," *Opt. Lett.* **32**, 2188–2190 (2007).

5. A. B. Parthasarathy et al., "Robust flow measurement with multi-exposure speckle imaging," *Opt. Express* **16**, 1975–1989 (2008).

6. S. M. S. Kazmi et al., "Chronic imaging of cortical blood flow using multi-exposure speckle imaging," *J. Cereb. Blood Flow Metab.* **33**, 798–808 (2013).

7. A. González Olmos et al., "Optimizing the precision of laser speckle contrast imaging," *Sci. Rep.* **13**, 17970 (2023).

8. W. J. Tom, A. Ponticorvo, and A. K. Dunn, "Efficient processing of laser speckle contrast images," *IEEE Trans. Med. Imaging* **27**, 1728–1738 (2008).

9. C. Wang et al., "Robust quantitative single-exposure laser speckle imaging with true flow speckle contrast in the temporal and spatial domains," *Biomed. Opt. Express* **10**, 4097–4114 (2019).

10. M. M. Qureshi et al., "Quantitative blood flow estimation in vivo by optical speckle image velocimetry," *Optica* **8**, 1092–1101 (2021).

11. D. A. Rivera and C. B. Schaffer, "Quasi-analytic solution for real-time multi-exposure speckle imaging of tissue perfusion," *Biomed. Opt. Express* **14**, 3950–3967 (2023).

12. X. Orlik, E. Colin, and A. Plyer, "Standardizing laser speckle orthogonal contrast imaging: achieving reproducible measurements across instruments," *Photonics* **11**, 585 (2024).

13. A. Santorelli, C. Smith, and A. K. Dunn, "Impact of logarithmic and linear spaced camera exposure durations for multi-exposure speckle imaging," *Proc. SPIE* **12378**, 1237803 (2023).

14. E. James, S. Powell, and P. Munro, "Simulation of statistically accurate time-integrated dynamic speckle patterns in biomedical optics," *Opt. Lett.* **46**, 4390–4393 (2021).

15. D. D. Postnov et al., "Dynamic light scattering imaging," *Sci. Adv.* **6**, eabc4628 (2020).

16. C.-Y. Yu et al., "Design and validation of a convolutional neural network for fast, model-free blood flow imaging with multiple exposure speckle imaging," *Biomed. Opt. Express* **14**, 4439–4454 (2023).

17. F. Pain and M. Chammas, "Python and ImageJ implementations of spatial speckle contrast calculations and associated simulated and experimental datasets," https://doi.org/10.5281/zenodo.14178966 (2024).

18. S. Prahl, "pyspeckle: python code for creation and analysis of laser speckle," https://doi.org/10.5281/zenodo.8311677 (2023).

19. C. Liu et al., "Choosing a model for laser speckle contrast imaging," *Biomed. Opt. Express* **12**, 3571–3583 (2021).

20. S. Zheng and J. Mertz, "Direct characterization of tissue dynamics with laser speckle contrast imaging," *Biomed. Opt. Express* **13**, 4118–4133 (2022).

21. B. Liu et al., "Dynamic light scattering and laser speckle contrast imaging of the brain: theory of the spatial and temporal statistics of speckle pattern evolution," *Biomed. Opt. Express* **15**, 579–593 (2024).

Marc Chammas received his PhD in biomedical optics from the University of Paris Saclay in 2023. His research interests include medical imaging, particularly methodologies for laser speckle imaging and the use of machine learning in biomedical image analysis. He is working on intellectual property as a patent engineer focusing on inventions in medical technologies and artificial intelligence for healthcare.
**Frédéric Pain** received his PhD in biomedical instrumentation from the University of Caen in 2000. He received his habilitation à Diriger des Recherches from Université Paris-Saclay in 2012. For the last 15 years, his research interests have been focused on a wide field *in vivo* optical functional imaging methods such as intrinsic optical imaging, autofluorescence imaging, and laser speckle contrast imaging.