1. Introduction

Diffusing-wave spectroscopy (DWS), the extension of quasi-elastic light scattering (QELS) to turbid, multiple scattering media [1, 2], has become an important tool for investigating soft materials such as colloidal suspensions, gels, foams and granular materials. DWS provides access to the microscopic dynamics within these materials by measuring the temporal autocorrelation function g ⁽²⁾(r,x) = 〈I(r,t)I(r,t + τ)〉/〈I(r,t)〉² of the fluctuating intensity I(r,t) of the speckle pattern at the position r generated by illuminating the sample with light with large coherence length. Similar to the situation in an interferometer, motions of the scattering particles within the volume swept by the diffuse photon cloud lead to fluctuations in the phases and hence to intensity fluctuations. The field autocorrelation function g ⁽¹⁾(r,t) = 〈E ^*(r,t)E(r,t + τ)〉/〈|E(r,t)|²〉 is related to the measured g ⁽²⁾(r,τ) by the Siegert relation g ⁽²⁾(r,τ) = 1 + β |g ⁽¹⁾(r, t)|² , where β is the coherence factor. The field autocorrelation function yields detailed information on the mean-squared displacement 〈Δr ²(τ)〉 of scatterers within the time τ by the relation [3]

Here k ₀ = 2π/λ is the wavenumber of light in the medium, and P(r,s) is the normalized distribution of photon path lengths s at the position r for a source located at the origin. The quantity l ^* is the transport mean free path length which characterizes the length scale over which the photon direction is completely randomized. In contrast to QELS, the accumulated phase shift after s/l ^* » 1 scattering events leads to a decay of g ⁽¹⁾(r, t) already for typical displacements Δr much smaller than the wavelength of light, λ. As the decay of g ⁽¹⁾(r, t) is determined by the temporal evolution of the mean-square displacement, different types of motions within the turbid medium, such as laminar shear flow or Brownian motion, can be distinguished by the shape of the autocorrelation function [4, 5, 6]. This enables DWS to detect dynamical heterogeneities deep within a turbid medium even when scattering or absorption contrasts are entirely absent [7, 8, 9]. More recently, the dynamical contrast of DWS has been used to detect dynamical heterogeneities within biological tissue, e.g. skin burns [10], perfusion in tumors [11, 12, 13], exercise-dependent dynamics in muscle tissue [14] and functional activation in cortical tissue in rodents [15, 16, 17, 18] and humans [19, 20]. For systems with rapidly fluctuating speckles and stationary average speckle intensity 〈I(r,t)〉, the ensemble-average intensity autocorrelation function g ⁽²⁾(τ) equals the time average $g_t^{\left(2\right)}$(τ) = 〈I(r,t)I(r,t + τ)〉_t/ 〈I(r,t)${〉}_t^2$ which is typically computed by real-time digital or hardware correlators. Nevertheless, the time-average requires integration times T which are, depending on photon count rate, typically a factor 10² … 10⁵ larger than the decay time of g ⁽²⁾(r, t). In many soft materials, however, the dynamics is not stationary over the time required for measuring a time-averaged intensity autocorrelation function. Solutions to this problem are provided by multispeckle detection based on CCD arrays [21, 22] or echo techniques [23], which allow to measure ensemble-average g ⁽²⁾(r, τ) with integration times T not significantly longer than the decay time of g ⁽²⁾(r, t) . Likewise, biomedical applications of DWS are likely to be affected by transient and non-stationary dynamics due to e.g. pulsation, slow (0.1 Hz) vasomotion or, in the cerebral cortex, due to functionally enhanced perfusion after a stimulus, requiring schemes for time-resolved measurement of autocorrelation functions. However, the temporal resolution of standard CCD cameras is too low for DWS experiments probing dynamics in deep tissue, where field autocorrelation functions g ⁽¹⁾(r, t) decay within several 10 to 100 μs [20]. The first application of time-resolved QELS to biological tissue known to us is the pulsation-synchronized measurement of the blood flow in a rabbit femoral vein by Tanaka and Benedek [24]. Measurements of the laser Doppler linewidth from blood flowing in model capillaries showed a linear dependence on flow rate consistent with a linear velocity profile across the scattering volume. Measurements of the Doppler linewidth in vivo synchronized with the pulsation showed, on the other hand, no variation in the vein of the rabbit over the pulsation cycle. Laser Doppler experiments with integration times of 8 ms in the dog coronary artery by Nishihara et al. [25] using a spectrum analyzer showed, on the other hand, Doppler frequency shifts clearly following the pulsation. The limited bandwidth of analog laser Doppler does, however, severely restrict its ability to be used for non-invasive measurements of dynamics in deep tissue. On the other hand, most current types of multi-tau correlators used for sampling non-exponential DWS autocorrelation functions from deep tissue dynamics do not allow for direct gating which is required for pulsation-synchronized detection of DWS signals from tissue.

In this publication we propose a technique for measuring transient or non-stationary dynamics in deep tissue with sub-second temporal resolution, using two gated single-photon counting detectors and a standard ungated multi-tau correlator. The distortions of the intensity autocorrelation functions induced by periodic or quasi-periodic gating of the detector are accurately corrected by the autocorrelation function of a non-fluctuating reference signal which is detected by a second detector gated synchronously with the one detecting the signal from the sample. We apply the technique to pulsation-synchronized measurements of field autocorrelation functions from the human forearm and wrist above the radial artery. While the tissue dynamics measured on the forearm does not significantly vary over the pulsation cycle, data from the wrist reveal a modulation of up to 90% during a pulsation cycle, much larger than the pulsation-induced variations of the transmitted intensity.

2. Method and Experimental Setup

As a light source we use an external-cavity diode laser operating at λ = 802nm (Toptica TA100) (see Fig. 1). The laser beam is divided by a beam splitter into two arms: the reference arm passes through an attenuator onto an avalanche diode (APD2) operating in Geiger mode (Perkin-Elmer SPCM-AQR-15-FC), the other arm is used to illuminate the sample (colloidal latex suspension or tissue). Multiple scattered light is collected at a distance of about 20mm from the source with a single-mode fiber with a cut-off wavelength of 1200nm [26] and detected with a second avalanche diode (APD1, same type as APD2). At the operating wavelength λ = 802nm the receiver fiber transports 6 transverse modes [26], which does not significantly improve the signal-to-noise ratio of the field autocorrelation function, but allows for a higher accuracy in the measurement of the transmission. An external TTL trigger signal M(t) is used to gate both APDs simultaneously; when the gate is open, M(t) = 1, otherwise M(t) = 0. The detectors with gains K_l produce photon count signals n_l (t) = M(t)K_l I_l (t), where I_l (t) (l = 1,2) are the intensities impinging on the detectors. The photon counts n ₁(t) and n ₂(t) of APD1 and APD2, respectively, are fed into the 2 channels of a digital correlator (ALV5000E) operating in dual mode.

As the correlator can not be synchronized with the APDs, it is operating also when the APDs produce no output. The autocorrelation function $g_t^{\left(2\right)}$(τ) of the sample photon counts n ₁(t) measured by the time-averaging correlator is

Assuming that (i) the detectors are linear (i.e. constant gains K_l ), (ii) during each of the time windows when M(t) = 1 the sample signal I ₁(t) is characterized by the same field autocorrelation function, and (iii) that M(t) and I ₁(t) are statistically independent, we can factorize the integral in Eq. (2), yielding

where g ⁽²⁾(τ) is the (ensemble-averaged) intensity autocorrelation function. The regularisation factor α(τ) relating the measured photon count autocorrelation function $g_t^{\left(2\right)}$(τ) to g ⁽²⁾(τ) is given by

 

Fig. 1. Experimental setup for time-resolved DWS. BS: beam splitter, APD: avalanche photodiodes.

Download Full Size | PDF

We first consider the case of a periodic gating signal M(t) with N cycles within the entire experimental duration T, each cycle being characterized by the open time T ₁ and the closed time T ₂. The periodicity of α(τ) (with period T ₂) leads to a maximal distortion of the photon count autocorrelation function by a factor 1 + T ₂/T ₁ ≥ 1 for τ = 0, and to a periodic modulation of $g_t^{\left(2\right)}$ (τ) for times much larger than the decay time of g ⁽²⁾(τ). When g ⁽²⁾(τ) decays completely within the open time T ₁, the entire information on the sample dynamics is contained within the first period of α(τ) , and the regularisation factor reduces, for 0 ≤ τ ≤ T ₁, to

Following tissue dynamics in different phases of the pulsation cycle requires quasi-periodic gating since pulsation is not strictly periodic. The modulation function M(t) derived from, e.g., a pulse oxymeter, is then characterized by a fixed open time T ₁, while the time T ₂ between open periods varies in different cycles. Defining N _open as the number of open-gate intervals within the total integration time T, we obtain for the regularisation factor

Thus once we know T ₁ and T ₂ for the periodic case or T ₁, T and N _open for the quasi-periodic case, we can use Eqs. (5) or (6) to calculate the regularisation factor α(τ), yielding the ensemble-average intensity autocorrelation function g ⁽²⁾(τ) for selected periods. The estimation of α(τ) from the intensity trace which is delivered by many correlators is, however, unreliable, since the intensity trace which is formed by averaging the predetection signal over some 10 – 100 ms is not able to capture the open-gate intervals with sufficient accuracy, in particular for small number of cycles N _open < 100. More directly and reliably, the regularisation factor α(τ) can be measured with a second APD (APD2) which detects a constant intensity I ₂(t) and which is gated synchronously with APD1 by the same modulation signal M(t). The normalized autocorrelation function of the photon count signal M(t)K ₂ I ₂(t) for a constant intensity I ₂(t) is then, by Eqs. (3) and (4), equal to the regularisation factor α(τ).

In order to minimize distortions of the photon count autocorrelation function by detector nonlinearities [27] and afterpulsing [28], experiments were performed with average photon count rates (measured with ungated detection) between 200 kHz and 400 kHz. Adjusting the count rates of both sample and reference detectors to similar values (to within 10%) resulted in an unbiased determination of the field autocorrelation function g ⁽¹⁾(τ) for both periodic and quasi-periodic gated detection.

3. Results

To illustrate the basic features of the proposed technique, we have carried out three sets of experiments: first on a colloidal suspension, using periodic and quasi-periodic gating, second on the forearm, and third on the wrist of a human subject (over the radial artery). In the last two human experiments we have synchronized the detection with the pulsation.

3.1. Colloidal suspension

The sample used in these measurements was an aqueous colloidal suspension of polystyrene spheres with 04μm diameter (volume fraction 0.55%) in a cuvette with dimensions 5.5cm × 6.7cm × 15 cm. Source and receiver fibers were placed on the same side of the cuvette at a distance of 2.3cm.

For the periodic gating, we used a periodic TTL signal to trigger the APDs with intervals T ₁ = 200 ms and T ₂ = 600 ms. The corresponding time-averaged autocorrelation function $g_t^{\left(2\right)}$(τ) (Fig. 2a) shows, after a first decay at about 100μs, an inclined plateau followed by a second decay starting at about 10ms. The amplitude of $g_t^{\left(2\right)}$(τ) with gated detection is much larger here than the value g ⁽²⁾ (τ = 0) = 1.167 obtained with ungated detection using a 6-mode fiber [26]. This is explained by the fact that for gated detection, the measured intensity autocorrelation function $g_t^{\left(2\right)}$(τ) differs from the ensemble-averaged autocorrelation function g ⁽²⁾(τ) by the regularisation factor α(τ) ≥ 1 (see Eq. (3)). The second decay in $g_t^{\left(2\right)}$(τ) is due to the periodic modulation of the average sample photon count rate, as is revealed by a direct measurement of the regularisation factor α(τ) directly using the reference signal which reproduces the second decay in $g_t^{\left(2\right)}$(τ) (Fig. 2b). The amplitude of the regularisation factor α(τ) obtained from fitting Eq. (6) to the data, is T/(N _open T ₁) = 4.080, close to the theoretical value α(τ = 0) = 1 + T ₂/T ₁ = 4. From the decay rate T/(N _open $T_1^2$) = (20.6 ±0.2) s^-1 and the amplitude of α(τ) one obtains a value for the open time T ₁ = 198μs, in very good agreement with the theoretical value T ₁ = 200μs. The difference between theoretical and experimental values is due to the slightly different numbers of closed and open intervals during the measurement. The autocorrelation function g ⁽²⁾(τ) calculated using Eq. (3) has a single decay reflecting the translational motion of the particles and is in excellent agreement with a measurement of g ⁽²⁾(τ) by ungated detection (see Fig. 2c).

For the quasi-periodic gating, we generated the gating signal M(t) by transforming the analog output of an optical pulse oxymeter (Nellcor N595), measuring the pulsation of a human subject, into a TTL signal with T ₁ = 200 ms and T ₂ varying between 630 and 700 ms. Measured time-averaged autocorrelation functions $g_t^{\left(2\right)}$(τ) again show two widely separated decays associated with translational motion and non-constant average photon-count rate (see Fig. 2d). From a fit of Eq. (6) to the measured α(τ) (Fig. 2e) we obtain the amplitude T/(N _open T ₁) = 4.050 and the decay rate T/(N _open T ₂) = (- 21.8 ± 0.3) s^-1, corresponding to an effective open time T ₁ = 186 ms which is slightly lower than the theoretical value T ₁ = 200 ms. Correcting the time-averaged photon-count correlation function $g_t^{\left(2\right)}$(τ) by the measured regularisation factor α(τ) using Eq. (3) then yields an ensemble-averaged autocorrelation function g ⁽²⁾(τ) which is practically indistinguishable from g ⁽²⁾(τ) measured without gating (Fig. 2f).

 

Fig. 2. Circles: measured time-averaged photon-count autocorrelation functions $g_t^{\left(2\right)}$(τ) (a), (d), regularisation factors α(τ) (b), (e) and retrieved ensemble-averaged intensity autocorrelation functions g ⁽²⁾(τ) (c), (f) from a colloidal suspension. The data in the left column (a-c) were obtained with periodic gating with T ₁ = 200ms and T ₂ = 600ms. The right column (d-f) shows data obtained with quasi-periodic gating derived from the human pulsation signal with T ₁ = 200 ms and T ₂ varying between 630 and 700 ms. For comparison, g ⁽²⁾(τ) measured without gating is also shown in (c) and (f) (asterisks). The full lines in (b) and (e) are fits of Eq. (6) to the measured α(τ). Best-fit parameters are T/(N _open T ₁) = 4.080, T/(N _open T ₂) = (20.6 ± 0.2) s^-1 for (b) and T/(N _open T ¹) = 4.050, T/(N _open T ²) = (21.8±0.3)s^-1 for (e). Insets: differences between autocorrelation functions obtained with the gating technique and by an ungated experiment.

Download Full Size | PDF

3.2. Pulsation-synchronized DWS on the human forearm

In this experiment, we measured tissue dynamics at different phases of the pulsation. A sensor consisting of source and receiver fibers 2.3 cm apart was located on the left forearm of a human volunteer (see Fig. 3). The pulsation was monitored using the pulse oxymeter whose probe was placed at the tip of the ring finger. Its analog output was used for generating the gate control signal M(t) by an electronic circuit. In order to record the intensity autocorrelation function g (2)(τ) at different phases of the pulsation, we placed the detection window with duration T 1 = 200 ms with increasing delays with respect to the onset of the pulsation. The time span from the first to the last measurement point covered about two pulsation periods (at a heartbeat rate of 75/min), and the time between successive delays was 200 ms. Fig. 4 shows the pulsation waveform and the relative decay rate of the autocorrelation function which is given by the inverse of the decay time τd = ∫τ2 τ1 (g (2)(τ) - 1)dτ, normalized by τ d measured with an ungated experiment averaging over many pulsation cycles. We chose the integration times τ 1 = 0.8 × 10-6s and τ 2 = 1.024 × 10-3s which cover the full decay of g (2)(τ) without over-weighting the baseline.

 

Fig. 3. Experimental arrangement in the forearm experiment. The sensor consisting of source and receiver fibers is placed along the radius. Pulsation is monitored by a pulse oxymeter probe placed at the ring finger tip.

Download Full Size | PDF

We also analyzed measured autocorrelation functions using the solution of the correlation diffusion equation for a homogeneous half space [15], modelling the tissue dynamics by Brownian motion characterized by an effective diffusion coefficient. Using values for the tissue absorption coefficient μ _a = 0.23 cm^-1, reduced scattering coefficient μ′_s = 1/l ^* = 6.8cm^-1 and the refractive index n _m = 1.4 [29], we obtained the value D̄ = (2.95 ±0.16) × 10^-9cm²/s for the effective diffusion coefficient from an ungated measurement averaged over many pulsation cycles. Normalizing the effective diffusion coefficients obtained from time-resolved measurements with the average value D̄ we obtain a pattern which is very similar to the one of the relative decay rates (data not shown). The observation that at the forearm the DWS signal does not show any signature of pulsation is in line with pulsation-synchronized laser Doppler measurements inside the femoral vein of the rabbit [24] and reflects the absence of velocity and pressure waves in venous vasculature.

3.3. Pulsation-synchronized DWS on the human wrist

In this experiment, we placed the sensor on the left wrist of the same volunteer over the radial artery, with a distance between source and receiver of 2.3 cm (see Fig. 5). Except for the sensor location, the experiment protocol is the same as in Sec. 3.2. The measured autocorrelation function g ⁽²⁾(τ) shows a significantly faster decay during the systolic phase (see Fig. 6), indicating that the scattering particles in the tissue probed by the diffusing light move faster during the systolic phase than during the diastolic phase. Indeed, we observe that the autocorrelation functions differ most for relatively long lag times 20 μs ≤ τ ≤ 200 μs corresponding to photon paths which preferentially sample the superficial regions of the tissue where the artery is located. Unlike the situation at the forearm, we observe that the intensity autocorrelation function measured over the wrist cannot be described by a semi-infinite geometry with diffusing particles. The relative decay rate (see Fig. 7) increases steeply during the diastolic phase to a value 60 – 90% larger than the one measured during the diastolic phase of the pulsation. This pulsation-induced modulation of the decay rate is significantly larger than the one of the transmission measured which is less than about 5%.

 

Fig. 4. Oxymeter signal recorded at the ring finger tip (top) and relative decay rate of the intensity autocorrelation function g ⁽²⁾(τ) measured on the forearm for different phases of pulsation (bottom). The relative decay rate is the ratio between the decay rates ${\tau }_{\mathrm{d}}^{-1}$ measured with and without gating, respectively. During this experiment, the subject pulsation rate was about 75/min., thus the maximal delay is about two pulsation cycles. The error bars represent the standard deviation over 5 repeated measurements for each delay. In this experiment, the source-receiver distance is 2.3cm.

Download Full Size | PDF

 

Fig. 5. Experimental arrangement in the wrist experiment. Pulsation is monitored by a pulse oxymeter probe placed at the ring finger tip.

Download Full Size | PDF

 

Fig. 6. Intensity autocorrelation functions g ⁽²⁾(τ) recorded on the wrist for a systolic phase (time = [850,1050]ms, circles) and a diastolic phase (time = [1250,1450]ms, squares). The error bars represent the standard deviation over 5 repeated measurements. The source-receiver distance is 2.3cm.

Download Full Size | PDF

4. Discussion

As shown by the comparison of gated and ungated measurements on the colloidal suspension the separate measurement of the time-averaged autocorrelation function $g_t^{\left(2\right)}$(τ) and the regularisation factor α(τ) provides the ensemble-averaged intensity autocorrelation function g ⁽²⁾(τ) which reproduces the ungated measurement to within statistical errors.

The measurements on the forearm show that the dynamics of the sampled tissue do not significantly differ between different phases of pulsation (see Fig. 4). This observation might be due to that the tissue volume swept by the photon cloud, approximately up to 1 cm beneath the skin, contains no artery where red blood cells move with periodically modulated velocity. The good agreement of the measured g ⁽²⁾(τ) with the prediction for diffusive dynamics resembles the apparent Brownian dynamics in brain tissue [19, 20] and might have its origin in dynamical heterogeneities on length scales smaller than the transport mean free pathlength l ^* ≈ 1.5mm.

In contrast, the decay rates of the autocorrelation functions measured over the radial artery vary synchronously with the pulsation recorded at the ring finger tip (see Fig. 7). The transmitted intensity varies less than about 5%, consistent with the small variations of the partial blood volume during a pulsation cycle. Thus the large systolic increase of the decay rate of about 60 – 90% is likely to originate in the increased erythrocyte flow velocity during the systole, in line with the systolic increase of blood flow velocity observed in the dog coronary artery using Laser Doppler with an invasive probe [25]. Furthermore, the interpretation of the systolic increase of ${\tau }_{\mathrm{d}}^{-1}$ as an increase in arterial blood flow velocity is corroborated by the fact that in the radial artery, blood flow speed, pressure and local diameter vary in phase during the pulsation cycle [30, 31].

 

Fig. 7. Oxymeter signal recorded at the ring finger tip (top) and relative decay rate of the intensity autocorrelation function g ⁽²⁾(τ) measured on the wrist for different phases of pulsation (bottom). During this experiment, the subject heartbeat rate was about 75/min., thus the maximal delay is about two pulsation cycles. The error bars represent the standard deviation over 5 repeated measurements for each delay.

Download Full Size | PDF

5. Conclusions

In summary, we have presented a method for non-invasively probing transient dynamics in deep tissue with sub-second time resolution, using diffusing-wave spectroscopy with time-gated detection. This technique is compatible with standard multi-tau correlators without externally controllable gating function. While time-resolved measurements on the forearm show no significant difference for different phases of the pulsation cycle, measurements over the radial artery show that the autocorrelation function decays faster during the systolic phase than during the diastolic phase, consistent with the physiological pattern. The large variations of the DWS signal of up to 90% during a pulsation cycle are consistent with invasive measurements of arterial erythrocyte flow velocity and might thus allow to use pulsation-resolved DWS as a sensitive, non-invasive diagnostic tool for the assessment of deep tissue perfusion, e.g., after surgery or in intensive care units.