A new simulation method for laser speckle imaging to investigate hemodynamics

Speckle simulation is a powerful protocol to investigate the properties of speckle and evaluate image processing method. However, only static speckle images can be simulated by available methods without considering time-integrated effect of CCD. A time–integrated dynamic speckle simulation method basedon coherent imaging was developed. Through the new simulation method, the effect of speckle size on LSCI was investigated. The smaller the speckle size is, the higher the spatial resolution become.But the one-dimensional speckle size should exceed two pixels to sample the speckle pattern. The characteristics of existing speckle contrast imaging methods were studied based on spatial statistics, and optimal parameters are given to obtain accurate and less noisy image. In general, the new simulation method for laser speckle imaging is a powerful tool to monitor blood flow in vivo and lay a solid foundation for the study of hemodynamics.


Introduction
Port-wine stains (PWSs) are discolorations of the human skin caused by vascular anomalies [1]. The histological features of PWSs are characterized by enlarged and expanded capillaries and posterior venules with diameters of 30-300 μm and are present at depths of 200-800 μm in the superficial dermis. PWSs require early treatment given that they will influence the physical and psychological health of patients and tend to darken and thicken with age and growth. Pulsed dye laser (PDL) irradiation at 585/595 nm has become the standard choice for the clinical treatment of PWSs [2][3][4][5]. PDL is based on the theory of selective photothermolysis that was proposed by Anderson and Parrish [6]. This theory states that laser energy is selectively absorbed by hemoglobin in red blood cells. In PDL, the target lesion is heated and destroyed without damaging the surrounding tissue. The online monitoring of blood vessel morphology, blood flow, and thrombosis formation is essential for the precise treatment of vascular proliferative skin diseases. Nevertheless, conventional invasive measurement methods, such as micro-particle image velocimetry, are difficult to employ in the treatment of PWSs.
Various methods have been proposed to simulate speckle phenomenon. The simulation of speckle phenomenon, which ignores the physical process of light transmission in the tissue, but makes some assumptions about the statistical properties of the scattering field, such as uniform phase distribution, and then uses scalar diffraction theory to perform simulation analysis. Duncan et al. proposed a dynamic speckle simulation based on copula transformation [7]. The core idea is to use a static speckle field as the initial speckle field and the opposite scattering field as the target field, using the copula transformation theory to generate correlated field sequences between the two to achieve the purpose of simulating a dynamic speckle field. Based on this, Qiu proposed time-integrated speckle simulation method [8]. However, the copula method can only simulated temporal correlated speckle patterns and mix the spatial correlation. Song et al. introduce a simple method based on the concept of the coherent imaging principle, to generate speckle images with a pre-defined correlation matrix to mimic situations with complex correlations across speckle images [9]. But in this method the timeintegrated effect of CCD is not taken into account. In this work, a time-integrated dynamic speckle simulation method based on coherent imaging was developed for the simulation of speckle image to investigate hemodynamics.

Simulation method
Laser speckle is a random interference pattern that is produced by the addition of scattered coherent laser light with slightly different path lengths. The dynamic lightscattering theory and approximate model state that speckle intensity is modulated by the velocity, amplitude, and frequency of particle scattering. Since laser speckle contrast imaging (LSCI) is based on the relationship between the contrast value and the correlation time of the light intensity autocorrelation function to obtain the velocity of the scattering particles such as red blood cells, it is necessary to simulate the time-integated dynamic laser speckle images and study the speckle size on contrast images and invastigate the performance of some laser speckle image processing methods. The simulating program is written in MATLAB language. A sequence of statistically independent W values with a specific correlation r to 1 i e   have being generated by 1 2 where 1  and 2  are two two-dimensional statistically independent random variables that follow uniform distribution on the interval -  ( ， ) . 1  and 2  can be produced by setting two different seeds for a random number generator. The optical field in the object plane is W, then the statistically correlated speckle images which is the intensity distribution in the image plane can be expessed as: where F stands for Fourier thansform. H(x, y) is the aperture function of the lens. It acts as a low pass filter and the central region of H(x, y) that contains non-zeros values is called clear aperture. The ratio of the size of the clear aperture to the size of matrix W sets the minimum size of the speckles [10]. However, when the size of the clear aperture is equal to that of W, the imaging system is equivalent to free-space light propagation. In this case the speckle imges should be expressed as: In both cases, the relationship between the input correlation coefficients r and the correlation of the simulated speckle patterns  has been deduced and proved to be power 2.
Therefore, a sequence of speckle images with correlation  can be generated by adjusting r which stands for the the correlation coefficient of W and initial optical field 1 i e   . The fully developed speckle contrast images with predefined correlation distribution are generated through the above simulation program. However, every image is still a static speckle image essentially without considering the time-integrated effect of CCD. To simulate the timeintegrated process of CCD, the first n frames of the fully developed speckle image stack were averaged to get one frame of time-integrated speckle image, while larger n indicates longer exposure time.

Speckle images simulation of Brownian motion
There are two motion types of the scattering particles: Brownian motion (unordered motion) and laminar flow (ordered motion). The correlation coefficient between the initial speckle image and following ones is based on the type of particle motion. According to Dunn's research [11], the correlation coefficients of speckle images of particle with Brownian motion follows a negative exponential function expressed as: where c  is correlation time and  is the delay time. The input correlation r is set to be a three-dimensional matrix ( 600 600 51   ) with space as the first two dimensions and time as the third one. The spatial distribution is that the 150 th to 450 th rows are Brownian motion and the others are stationary part. It has been stated that the relationship between the input correlation coefficients r and the correlation of the simulated speckle patterns  has been deduced and proved to be power 2. Therefore, the input correlation r is expressed as: To simulate the temporal changes of speckle patterns of Brownian motion, c  was set to be 370 μs according to V.Rajan's research [12] and  was varied from 0 μs to 1.85 ms with an increment of 37 μs. The simulation result is shown in Figure 1. It is clear that the correlation between the initial speckle image and the following ones remains constant in the stationary part while the correlation decreases as a negative exponential function in Brownian motion region. That the simulation results fit well with the expected curve confirm the validity of our simulation method. Fig. 1. The correlation coefficient between the initial speckle image and the following ones for Brownian motion (red diamond) and stationary part (blue asterisk). The black dotted line is a negative exponential function.
In this section, to simulate the time-integrated process of CCD, the first 25 frames of the fully developed speckle image stack were averaged to get one frame of timeintegrated speckle image. Laser speckle spatial contrast algorithm ( 7 7  window size) was used to obtain the contrast image of static speckle and dynamic speckle. The regional flow velocity difference of the timeintegrated dynamic imaging speckle obtained by the improved method can be clearly reflected.

Complex flow structures of blood flow in a rat dorsal skin window
Laser speckle contrast imaging is a powerful tool to monitor the dynamics of blood flow in vivo. Laser speckle simulation is a powerful protocol to investigate the properties of laser speckle and evaluate image processing method without sacrifice of animal. And it can change only one parameter and remain others constant.
Based on the above simulation method, it is available to simulate complex flow structures of blood flow in a dorsal skin window. The vessel structure is shown in Figure 2. The capillary network was removed by Photoshop to highlight main vessel branch. The correlation matrix r was set according to P. Zakharov's research [13]: where b is the correlation time, x is the delay time, c characterizes the ratio of static part to the total detected light intensity according to reference [13] and it is the lowest correlation coefficient. The parameter setting for the blood flow simulation is shown in Table 1. The white light image of vessels in a dorsal window is segmented into a 5 gray-scale levels image shown as Figure 3. The correlation for different gray level is shown in Figure 4.    The simulation results are shown in Figure 5. By controlling the ratio of the size of the clear aperture to the size of matrix W, the speckle size in Figure 6 is 2 pixels. The speckle contrast image corresponds well with the correlation distribution.
In this section, the simulated speckle images are used to study the effect of speckle size on speckle contrast. The laser speckle spatial contrast analysis (LSSCA) method performs speckle contrast calculation in the spatial domain using a spatial window. LSSCA achieves high temporal resolution with the loss of spatial resolution. As the spatial calculation window size s s N N  increases, the spatial contrast increases. The statistical noise is lowest when the speckle size is 2 and the larger the spatial calculation window size is, the lower the statistical noise is when the speckle size is larger than 2. Therefore, 2 pixels is the optimal speckle size to obtain accurate and less noisy image.

CONCLUSIONS
In this paper, a time-integrated dynamic speckle simulation method based on coherent imaging was developed. The method was proved to be effective by simulating the speckle imaging of red blood cells flowing in a straight vessel. Complex flow structures of blood flow in a rat dorsal skin window were simulated, and the effects of speckle size and spatial window length on speckle contrast was investigated. The effect of speckle size on LSCI was investigated through this method. The smaller the speckle size is, the higher the spatial resolution become. But the one-dimensional speckle size should exceed two pixels to sample the speckle pattern, which is the so-called Nyquist criterion. Therefore, 2 pixels is the optimal speckle size to obtain accurate and less noisy image. In general, the new simulation method for speckle imaging is a powerful tool to monitor blood flow in vivo and lay a solid foundation for the study of hemodynamics.