1.

Introduction

Optical methods have become increasingly important as diagnostic tools for noninvasive biomedical imaging, since they offer great advantages to established methods such as positron emission tomography and magnetic resonance imaging in terms of portability, price, and weakness of interaction (absence of radiation damage, etc.).^{1, 2} In particular, near-infrared light at wavelengths $\lambda$ between around 750 and $950\mathrm{nm}$ provides significant advantages over visible or infrared light, since absorption by the most abundant chromophores in tissue, such as melanin and hemoglobin, or water, is weak, while at the same time scattering mean free paths are long enough that light can penetrate up to a few centimeters into tissue. However, the strong scattering of light is a generic obstacle to direct imaging of subsurface features in the body. It has its origin in fluctuations of the dielectric constant on the length scale of the wavelength of light, associated with structural elements such as collagen fibers, cell membranes, and subcellular organelles. Rather than ballistic, light propagates diffusively deep inside tissue, with a photon diffusion coefficient $\mathcal{D}=c{l}^{*}∕3$ determined by the speed of light in the medium $c$ and the transport mean free path length $l^{*}$ .³ This quantity is roughly the distance over which scattering completely randomizes the direction of an incoming photon. Diffusive light propagation in tissue thus strongly reduces the spatial resolution of optical imaging methods, in particular for structures deep within the body.

Several methods are currently used to extract, from the distribution of photon fluxes measured at the surface, information on subsurface structures with a scattering or absorption contrast: 1. time-of-flight experiments with gated detection allow for the discrimination of the early arriving ballistic photons and the reconstruction of optical density distributions;⁴ 2. low-coherence tomography operates in reflection and allows to probe near-surface structures with a longitudinal resolution as small as several $\mu \mathrm{m}$ ;⁵ 3. photon migration imaging (PMI) uses a light source modulated at megahertz frequencies and measures the attenuation and retardation of the diffusely transmitted light intensity that is detected some distance away from the source.⁶ Photon migration imaging requires some means of calculating distributions of optical densities and absorption cross sections within the volume from an inversion of the measured data.⁷

PMI has recently been applied to noninvasive functional imaging of cerebral activity after stimulation of specific brain areas such as the visual and motor cortices. Using a back-projection method, low-resolution maps of cortical blood oxygenation following motor stimulation have been obtained with a temporal resolution of $160\mathrm{ms}$ .^{8, 9, 10} These experiments show that the total hemoglobin concentration increases within about $4\mathrm{s}$ after the stimulation by up to $0.3\mu \mathrm{M}$ above its average level of about $70\text{to}80\mu \mathrm{M}$ . At much shorter times after stimulation (tens of milliseconds), PMI experiments have revealed small stimulus-related signal changes that are interpreted as an increase in the cortical turbidity.^{11, 12, 13} This fast optical signal had previously been observed in reflectance measurements on exposed cortical areas.^{14, 15} Since the time scale of the fast optical signal is too short for significant optical changes due to increased blood flow, it is believed that it reflects the increased scattering from nerve cells following electrical activation, as has been directly observed in isolated neurons.¹⁶ Increased activation-related scattering, in turn, is attributed to mechanisms such as cell swelling, increased concentration of vesicles, and ion density gradients.¹⁷

In addition to the attenuation and retardation of the diffusely transmitted light intensity measured in PMI, a wealth of information is contained in the speckle pattern of multiple scattered light, which forms when a light source with a large coherence length is used. This speckle pattern is a complicated interferogram arising from the coherent superposition of multiple scattered electric fields at the detector. It contains information on the positions of all scatterers, such as cell membranes, erythrocytes, and mitochondria, in the volume swept by the diffuse photon cloud. Similar to the situation in a simple interferometer, motions of the scatterers (due to, e.g., shear flow or Brownian motion) impose changes in the phases of the scattered electric field $E(\mathbf{r},t)$ at the position $\mathbf{r}$ and time $t$ , which leads to temporal fluctuations of the detected intensity. The intensity $I(\mathbf{r},t)={\mid E(\mathbf{r},t)\mid }^2$ of the surface speckle pattern thus rapidly fluctuates at locations whose underlying dynamics are fast, while regions with slow dynamics will result in slow fluctuations of the surface speckle pattern. From spatially resolved measurements of the surface speckle contrast using, e.g., a CCD camera, time-resolved maps of the superficial blood flow in cortical tissue can be obtained.¹⁸ This technique has been used for the study of hemodynamic changes in exposed rat brains during cortical spreading depression^{19, 20} and somato-sensory activation.²¹

On the other hand, measurements of the temporal autocorrelation function $g^{\left(1\right)}(\mathbf{r},\tau )=⟨E^{*}(\mathbf{r},t)E(\mathbf{r},t+\tau )⟩∕⟨{\mid E(\mathbf{r},t)\mid }^2⟩$ of the fluctuating electric field amplitude $E(\mathbf{r},t)$ of the speckle pattern yields detailed information on the mean-squared displacement $⟨\Delta r^2\left(\tau \right)⟩$ of scatterers within the time $\tau$ by the relation^{3, 22}

Experiments using DWS for brain imaging are still relatively scarce. Aiming at mapping blood flow velocities in cortical vessels, Lohwasser and Soelkner have tried to account for the effects of multiple light scattering in laser Doppler spectra from brain phantoms.²⁹ DWS experiments by Cheung on a semi-exposed rat brain revealed that the field autocorrelation function shows a faster decay on hypercapnia.³⁰ The interpretation in terms of increased flow velocity of erythrocytes in the cortical vessels was also supported by parallel NIRS and DWS experiments on a rat brain during focal ischemia.³¹ Using motor cortex stimulation by finger opposition, Durduran ³² were able to show that a similar acceleration of $g^{\left(1\right)}(\mathbf{r},\tau )$ can be measured fully noninvasively over the human motor cortex. An analysis of the DWS data with a semi-infinite one-layer model, using an empirical correction factor accounting for the underestimation of the cortical diffusion coefficient by the one-layer model, resulted in functional increases of the cortical dynamics of about 40%, which were interpreted due to the functional enhancement of the cortical blood flow rate.

While the shape of the measured field autocorrelation functions reported by Cheung ³⁰ deviates from predictions for simple diffusion or random flow, the correction factor used by Durduran ³² suggests that using DWS to quantitatively determine cortical perfusion rates in an entirely noninvasive way is likely to be complicated by the heterogeneous optical and dynamical properties of the head consisting of scalp, skull, cerebrospinal fluid, and gray and white matter. Model calculations of autocorrelation functions from optically and dynamically heterogeneous multilayer media based on the solution of the correlation-diffusion equation do indeed show a complicated decay of $g^{\left(1\right)}(\mathbf{r},\tau )$ , which does not factorize into contributions from the individual layers.³³

We show that the DWS autocorrelation functions measured from the human motor cortex through intact scalp and skull can be analyzed quantitatively with a three-layer model allowing for different optical and dynamical parameters in the scalp, skull, and cortex. This allows for separating the cortical dynamics from tissue optical parameters and scalp dynamics reflecting peripheral perfusion. For right-handed subjects, we find a significant increase of the cortical diffusion coefficient when a finger opposition task is performed. In accordance with the enhanced activation of the contralateral hand area of the primary somato-sensory and primary motor cortex by the finger opposition task,^{34, 35, 36, 37} we find that the cortical diffusion coefficient for contralateral (right-hand) stimulation is significantly higher than the one for ipsilateral (left-hand) stimulation. An analysis of the possible origins of the observed increase in the cortical diffusion coefficient indicates that for the present motor stimulation protocol, the DWS signal is mainly dominated by functional increases of blood flow rate and less by increases in the blood volume, in accordance with the conclusions arrived at by Durduran ³²

2.

Methods

2.3.

Diffusing-Wave Spectroscopy

A diode laser operating at a wavelength ${\lambda }_0=802\mathrm{nm}$ (Toptica TA100) was coupled into a multimode optical fiber (core diameter $50\mu \mathrm{m}$ ; numerical aperture 0.22) serving as an illuminating source. Multiple scattered light was collected by a few-mode optical fiber guiding about six transverse modes³⁹ at a distance of $20\text{to}25\mathrm{mm}$ from the source and detected with an avalanche diode operating in Geiger mode (Perkin-Elmer SPCM-AQR-15-FC). The intensity autocorrelation function $g^{\left(2\right)}(\mathbf{r},\tau )=⟨I(\mathbf{r},t)I(\mathbf{r},t+\tau )⟩∕{⟨I(\mathbf{r},t)⟩}^2$ was then computed from the amplified and discriminated detector output by a digital multitau correlator (ALV5000E) at lag times between $200\mathrm{ns}$ and $3.1\mathrm{s}$ . Source and receiver fibers were mounted on a rigid fixation device that the subject could wear like a helmet. The intensity of the illuminating light was adjusted to be less than the maximum allowable intensity of $4\mathrm{mW}∕{\mathrm{mm}}^2$ .⁴⁰ Photon count rates of $20\text{to}80\mathrm{kHz}$ high enough for computation of intensity autocorrelation functions within the measurement period of $100\mathrm{s}$ were obtained by increasing both the incident laser power and the illuminated area on the scalp, thereby still keeping the peak intensity below the safety limit. This procedure was especially important for measurements on strongly pigmented subjects. The sensor, comprising light source and detector, was then fixed over the motor area on the left hemisphere [area C3 in the international 10–20 system for electroencephalography (EEG)⁴¹], using vertex and ears as external landmarks (see Fig. 1 ). The direction from source to detector was arranged perpendicular to the direction of the central sulcus to facilitate the photon cloud to cross the area activated most by the motor task. Before fixing the source and detector, the hair was combed away so that the illuminator and detector had good optical contact to the scalp. During both the resting baseline and the task periods, each lasting $100\mathrm{s}$ , we recorded ten normalized intensity autocorrelation functions, each with a duration of $10\mathrm{s}$ . To make sure that the amplitude of the measured intensity autocorrelation function was due to the few-mode receiver optics only (and not due to contributions of static scattering to the signal⁴²), we measured $g^{\left(2\right)}(\mathbf{r},\tau )$ , using the same setup, from a turbid latex sample immediately before and after the measurements on the subjects.

Fig. 1

Experimental geometry for the DWS experiments. N is nasion, I is inion, Cz is vertex, and SC is central sulcus. Source and receiver fibers (denoted by S and R, respectively) are placed over the C3 area of the somato-sensory cortex, approximately perpendicular to the central sulcus. For the right-handed subjects investigated here, finger opposition with the left hand corresponds to ipsilateral stimulation, and right-hand finger opposition to contralateral stimulation.

3.

Results

Figure 2 shows the reduced intensity autocorrelation functions $[g^{\left(2\right)}(\mathbf{r},\tau )-1]∕{\beta }_{\mathrm{coh}}$ from a selected subject obtained from the left hemisphere (C3) during baseline and activation, respectively. Data are averages over two blocks of ten runs each. The autocorrelation function recorded during contralateral stimulation shows a significantly faster decay than during the resting baseline period. Modeling the autocorrelation function by the solution of the correlation-diffusion equation for a semi-infinite one-layer geometry, assuming free diffusion as the mechanism for speckle decorrelation, gave rather poor agreement with the experimental data. Modeling the dynamics by a random flow $⟨\Delta r^2\left(\tau \right)⟩=⟨v^2⟩{\tau }^2$ with a Gaussian flow velocity $v$ yielded, due to the increased curvature of $g^{\left(1\right)}(\mathbf{r},\tau )$ , even poorer agreement than the free diffusion model. These discrepancies are resolved when we use the solution of the correlation-diffusion Eq. 2 for a three-layer medium, representing scalp and skull by finite layers and the cortex as the innermost, semi-infinite layer. With this dynamic three-layer model, the fitting procedure yields a cortical diffusion coefficient $D_3=4.9\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ during the resting period and $D_3=1.3\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ during the stimulation period. The best-fit value ${\Delta }_2=4.7\mathrm{mm}$ for the skull thickness obtained from the DWS data was, to within experimental error, identical for both baseline and stimulation data. ² Although the skull was modeled as a completely static layer, its incorporation was found to be necessary to account for the decay of the autocorrelation function at intermediate times $20\mu \mathrm{s}⩽\tau ⩽50\mu \mathrm{s}$ . The decay of the autocorrelation function at long times is dominated by the scalp diffusion coefficient $D_1$ , which is found to increase with stimulation (see below).

Fig. 2

Reduced intensity autocorrelation function $[g^{\left(2\right)}(r,\tau )-1]∕{\beta }_{\mathrm{coh}}$ measured over the C3 area of the motor cortex of a right-handed subject for a source-receiver distance of $20\mathrm{mm}$ . The blue circles are the resting baseline data; the red circles are the data for contralateral stimulation. The error bars are standard deviations from two blocks of ten runs each. Lines are fits of the analytical solution of the three-layer correlation-diffusion Eq. 2 to the data. Best-fit values for diffusion coefficients $D_1$ and $D_3$ of scalp and cortex, respectively, are $D_1=1.0\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $D_3=4.9\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ for the resting baseline period, and $D_1=1.5\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $D_3=1.3\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ for the stimulation period.

When measuring the autocorrelation function after the end of the motor task, the correlation functions were found to be identical with the baseline correlation functions before stimulation. Analyzing the individual autocorrelation functions from two blocks yields information on the distribution of the fitting parameters $D_1$ and $D_3$ . Due to the limited photon statistics and decay times in the $0.1\text{-}\mathrm{ms}$ range, individual autocorrelation functions recorded over $10\mathrm{s}$ are noisier than the averaged autocorrelation function. Nevertheless, the fitting routine converges, and, using a fixed skull thickness ${\Delta }_2=4.7\mathrm{mm}$ , the average best-fit values $D_1=(1.1\pm 0.3)\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $D_3=(4.7\pm 1.2)\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ for the baseline, and $D_1=(1.9\pm 0.5)\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $D_3=(1.4\pm 0.5)\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ for stimulation, agree well with the diffusion coefficients obtained from fitting the averaged data. Furthermore, this analysis shows that the differences in $D_3$ and $D_1$ on motor stimulation observed in the averaged autocorrelation functions are indeed significant.

When the stimulation is exerted by the left hand, ipsilateral to C3, the autocorrelation function recorded from C3 is found to decay significantly more slowly than for contralateral stimulation (see Fig. 3 ). This slower decay is mainly due to the fact that the cortical diffusion coefficient for contralateral stimulation $D_3=1.3\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ is larger than the one for ipsilateral stimulation $D_3=9.8\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ . In contrast, the scalp diffusion coefficients are not significantly different between ipsilateral and contralateral stimulation, although they are both increased with respect to the value for the resting period.

Fig. 3

Reduced intensity autocorrelation function $[g^{\left(2\right)}(r,\tau )-1]∕{\beta }_{\mathrm{coh}}$ measured over the C3 area of the motor cortex of a right-handed subject for a source-receiver distance of $20\mathrm{mm}$ . The magenta triangles are the data for ipsilateral (left-hand) stimulation; the red circles are the data for contralateral stimulation. The error bars are standard deviations from two blocks of ten runs each. Lines are fits of the analytical solution of the three-layer correlation-diffusion Eq. 2 to the data. Best-fit values for diffusion coefficients $D_1$ and $D_3$ of scalp and cortex, respectively, are $D_1=1.2\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $D_3=9.8\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ for ipsilateral stimulation and $D_1=1.5\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $D_3=1.3\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ for the contralateral stimulation.

Autocorrelation functions obtained from the 11 subjects measured during resting and stimulation periods showed equally good agreement with the predictions of the correlation-diffusion theory. The absolute values of the cortical diffusion coefficient were found to vary over more than a decade within the group of 11 subjects, both for the baseline and stimulation periods. However, normalizing the value of the cortical diffusion coefficient measured during stimulation periods with the one for the resting period for each subject, reveals that across the subject group the finger opposition task produces a significant increase of the cortical diffusion coefficient $D_3$ measured at C3 during contralateral stimulation by 38.6% [ $t$ -test (two-tailed): $t\left(10\right)=2.81$ , $p=0.018<0.050$ ], as shown in Fig. 4(a) . Ipsilateral (left hand) stimulation results in an increase of $D_3$ by about 20.4%, which is, however, not as significant as the one observed for contralateral stimulation[ $t$ -test (two-tailed): $t\left(10\right)=2.20$ , $p=0.052>0.050$ ]. Direct comparison of the relative cortical diffusion coefficients for contra- and ipsilateral stimulation shows a significant hemispheric asymmetry between the two responses [ $t$ -test (two-tailed): $t\left(10\right)=2.69$ , $p=0.023<0.050$ ].

Interestingly, the finger opposition task leads also to an increase of the scalp diffusion coefficient $D_1$ , albeit not as large as the increase of the cortical diffusion coefficient [see Fig. 4(b)]: contralateral stimulation leads to a 20.9% increase in $D_1$ above the baseline value [ $t$ -test (two-tailed): $t\left(10\right)=2.55$ , $p=0.029<0.050$ ], while an increase by about 16.7% can be seen for ipsilateral stimulation [ $t$ -test (two-tailed): $t\left(10\right)=2.86$ , $p=0.017<0.050$ ]. Direct comparison shows, on the other hand, that the difference in $D_1$ between contra- and ipsilateral responses is not significant [ $t$ -test (two-tailed): $t\left(10\right)=0.731$ , $p=0.480>0.050$ ].

4.

Discussion

The observed large spread of the value of $D_3$ between subjects might be due to the varying accuracy of localizing C3 using only external landmarks. A combination of DWS and independent measurements of the optical parameters by NIRS (such as in Ref. 32) or time-of-flight methods might help to resolve this issue.

In our experiments with the 11 subjects, in two subjects the response of the cortical diffusion coefficient to ipsilateral stimulation is larger than the one to contralateral stimulation, opposite to what would be expected from EEG data. This anomaly might result from a variant in the location of the contralateral somato-motor hand representation that may not be near the C3-scalp location in these subjects, or from an anomaly in the motor organization, which is known to occur in about 5% of the subjects.³⁷

As the theoretical expressions [Eqs. 16, 17] for the three-layer model used for extracting cortical and scalp diffusion coefficients contains ten parameters (layer thickness ${\Delta }_{1,2}$ , absorption path lengths $l_{1,2,3}^{\left(\mathrm{a}\right)}$ , transport mean free path lengths $l_{1,2,3}^{*}$ , and diffusion coefficients $D_{1,3}$ ) when the skull is treated as a static medium, it is important to know how robust the determination of $D_3$ , $D_1$ , and ${\Delta }_2$ is with respect to uncertainties in the optical parameters $l_n^{\left(\mathrm{a}\right)}$ and $l_n^{*}$ of scalp, skull, and cortex, and in the local thickness of the skin. In particular, an uncontrolled increase in the cortical transport mean free path could be interpreted as a decreased cortical diffusion coefficient. Conversely, a decrease in the cortical absorption coefficient could be interpreted as an increase of the diffusion coefficient. Weaker influences on the determination of the cortical diffusion coefficient are expected from changes in the optical properties of the scalp and skull, since these affect the decay of the autocorrelation function mainly at long times.

We have thus modeled autocorrelation functions with largely varying values of the optical and dynamical parameters using the analytical solution of the correlation-diffusion equation for a three-layer model of the human head, and compared these predictions with measured data. Figure 5(a) shows a plot of the square root of the sum of squared residuals $\chi$ between the measured and calculated intensity autocorrelation functions as a function of $D_3$ and $l_3^{\left(\mathrm{a}\right)}$ . Varying the cortical diffusion coefficient from $D_3=2.5\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ to $D_3=1.0\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ and the cortical absorption length from $l_3^{\left(\mathrm{a}\right)}=1.5\mathrm{cm}$ to $l_3^{\left(\mathrm{a}\right)}=6.0\mathrm{cm}$ , we find that $\chi$ shows a well-defined minimum at the best-fit values $D_3=4.9\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $l_3^{\left(\mathrm{a}\right)}=3.0\mathrm{cm}$ for the cortical diffusion coefficient and the cortical absorption length, respectively. However, fits with only slightly higher $\chi$ could be obtained by reducing $D_3$ to below $2.5\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and increasing $l_3^{\left(\mathrm{a}\right)}$ to about $5.0\mathrm{cm}$ . This correlation between $l_3^{\left(\mathrm{a}\right)}$ and $D_3$ reflects the fact that higher absorption cuts off long photon paths, which in turn might be interpreted as slowed-down dynamics (i.e., smaller $D_3$ ). An overestimation of $l_3^{\left(\mathrm{a}\right)}$ and a concurrent underestimation of $D_3$ can occur if the data quality is not sufficient to capture the rounding of $g^{\left(1\right)}(\mathbf{r},\tau )$ at short lag times, which is the hallmark of absorption.⁴⁸ Rather than corresponding to a physically meaningful value, the very large value of the cortical absorption length of $5.0\mathrm{cm}$ at the second minimum in $\chi (D_3,l_3^{\left(\mathrm{a}\right)})$ is more likely to reflect the noise in the data at short lag times.

Fig. 5

Plot of the sum of squared residuals $\chi$ between the prediction of the correlation-diffusion Eq. 2 for a three-layer model and experimental baseline data as a function of cortical diffusion coefficient $D_3$ and (a) cortical absorption length $l_3^{\left(\mathrm{a}\right)}$ , and (b) cortical transport mean free path length $l_3^{*}$ . Parameters for the calculation are $l_1^{\left(\mathrm{a}\right)}=3.97\mathrm{cm}$ , $l_1^{*}=0.0509\mathrm{cm}$ , ${\Delta }_1=0.10\mathrm{cm}$ for the scalp; $l_2^{\left(\mathrm{a}\right)}=3.97\mathrm{cm}$ , $l_2^{*}=0.0509\mathrm{cm}$ , ${\Delta }_2=0.47\mathrm{cm}$ for the skull; and $l_3^{\left(\mathrm{a}\right)}=2.96\mathrm{cm}$ , $l_3^{*}=0.060\mathrm{cm}$ , ${\Delta }_3=\infty$ for the cortex. Diffusion coefficients are $D_1=1.0\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ for the scalp, $D_2=0$ for the skull, and $D_3=4.9\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ for the brain. The asterisks mark the minimum of $\chi$ .

Bounds for the activation-related changes in $l_3^{\left(\mathrm{a}\right)}$ can be obtained from independent measurements using photon migration or time-of-flight methods. Recent experiments show that the concentration of deoxy-hemoglobin in the motor cortex decreases on activation by up to about $0.2\mu \mathrm{M}$ , while the increased oxygen consumption is reflected by an increase of the oxy-hemoglobin concentration by about $0.6\mu \mathrm{M}$ over the average total hemoglobin concentration of about $74\mu \mathrm{M}$ .¹⁰ At the wavelength ${\lambda }_0=802\mathrm{nm}$ used in our experiments, which is close to the isosbestic point of oxy-/deoxy-hemoglobin, relative changes in the cortical absorption coefficient $1∕l_3^{\left(\mathrm{a}\right)}$ are thus expected to be given by the changes in the total hemoglobin concentration, which are of the order of 1% around the baseline level. From our sensitivity plot, however, such changes would lead to changes in the cortical diffusion coefficient of no more than a few percent. On the other hand, we find that the cortical diffusion coefficient measured at C3 for contralateral stimulation is increased by 38.6% in the group average, which is clearly not explained by the perfusion-related increases in erythrocyte concentration or blood volume.

On the other hand, a plot of $\chi$ as a function of $D_3$ and $l_3^{*}$ [see Fig. 5(b)] shows that the minimum of $\chi (D_3,l_3^{*})$ around the optimal fit parameter $D_3=4.9\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ is quite pronounced. Varying the cortical transport mean free path length $l_3^{*}$ from $l_3^{*}=0.03\mathrm{cm}\text{to}{l}_3^{*}=0.12\mathrm{cm}$ does not significantly affect the value of $\chi$ for the cortical diffusion coefficient fixed at the best-fit value $D_3=4.9\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ . This shows that the determination of the cortical diffusion coefficient from fitting the full three-layer model to measured correlation functions is quite robust against uncertainties in the estimation of the cortical transport mean free path length.

Sensitivity plots of $\chi$ as a function of the cortical diffusion coefficient $D_3$ and the optical and dynamical parameters of the scalp show that the behavior of $g^{\left(1\right)}(\mathbf{r},\tau )$ at short times (governed by $D_3$ ) is sufficiently well decoupled from its long-time behavior that even sizeable changes in scalp parameters ${\Delta }_1$ , $l_1^{\left(\mathrm{a}\right)}$ and $l_1^{*}$ do not influence the determination of $D_3$ significantly. Conversely, the determination of the scalp diffusion coefficient $D_1$ was found to be largely unaffected by variations in cortical optical parameters.

Finally, the dependence of $\chi$ on $D_1$ and $D_3$ (Fig. 6 ) shows a single minimum for a wide range of diffusion coefficients $0.5\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}⩽D_1⩽2.0\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ and $2.5\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}⩽D_3⩽1.0\times {10}^{-8}{\mathrm{cm}}^2∕\mathrm{s}$ , which is rather flat in the direction of $D_1$ . This means that the estimate of the cortical diffusion coefficient $D_3$ is rather insensitive to uncertainties in value of the scalp diffusion coefficient $D_1$ .

It should be noted that the same qualitative behavior is shown when $\chi$ is calculated between theory and experimental data recorded during stimulation.

In addition, from the considerations on the robustness of the fitting procedure with respect to uncertainties in optical parameters, the accuracy in determining cortical diffusion coefficients also depends on the noise in the measured data arising from finite photon counting statistics. Noise in correlation functions is itself correlated and difficult to estimate.⁴⁹ In general, autocorrelation functions measured with low average photon count rate or fast decay times are particularly prone to errors. To assess the variance of $D_3$ arising from noise in the data, we analyzed the data for subject 3. These data show both a low average count rate of $21\mathrm{kHz}$ , and a short average decay time $\overline{\tau }=38\mu \mathrm{s}$ . ³ Analyzing the 20 autocorrelation functions from two blocks individually with fixed coherence factor ${\beta }_{\mathrm{coh}}=0.075$ ⁴ and fixed skull thickness ${\Delta }_2=5.2\mathrm{mm}$ yielded cortical diffusion coefficients $D_3=(6.7\pm 2.2)\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ during the resting baseline period and $D_3=(8.3\pm 1.7)\times {10}^{-9}{\mathrm{cm}}^2∕\mathrm{s}$ for contralateral stimulation. These values are in good agreement with the ones obtained from the analysis of the averaged autocorrelation functions. The analysis of these data shows that with even considerable noise in the measured data, the difference in cortical diffusion coefficients between resting baseline and stimulation periods is highly significant within the datasets from a single subject, at least for contralateral stimulation where the somato-sensory activation can be expected to be maximal.

In modeling our data with a three-layer geometry consisting of scalp, skull, and cortex, we have neglected the cerebrospinal fluid (CSF) layer between skull and cortex, whose absorption and transport mean free path lengths are both much larger than its thickness. Simulations have shown that the presence of the nearly transparent CSF can influence the spatial sensitivity profile or the path length distribution of diffusing photons;⁴⁶ for a thickness exceeding about $1\mathrm{mm}$ , photon propagation in the presence of the CSF is no longer well described by the diffusion approximation,⁵⁰ which is of particular importance for the accurate inversion of photon migration data. Light guiding by the CSF for large source-receiver distances reduces the photon flux through the cortex and, hence, the number of detected photons that can pick up its dynamics. The simulations show, on the other hand, that for the source-receiver distance used in our experiments, light-guiding effects should be small. The marked difference in the short-time behavior of $g^{\left(1\right)}(\mathbf{r},\tau )$ between resting and stimulation periods shows that with the present source-receiver distance, enough detected photons have indeed crossed the CSF and passed through the cortex for a measurable activation-related increase of the cortical diffusion coefficient. Surprisingly, the predictions for $g^{\left(1\right)}(\mathbf{r},\tau )$ from the three-layer model in the diffusion approximation show very good agreement with the data, even at the longest times where the effect of the CSF could be expected to be strongest. Indeed, Monte-Carlo simulations of DWS autocorrelation functions for four-layer systems, including a nonscattering CSF layer, show a shape that is very similar to the one for three-layer systems. An analysis of four-layer simulation data with a three-layer model shows that the obtained value of the cortical diffusion coefficient is underestimated only slightly.

Clearly our analysis of the sensitivity of the fitting procedure on uncertainties in the optical and dynamical parameters shows that within realistic bounds for the optical and dynamical tissue parameters, the faster decay of the autocorrelation function during a motor task can be assigned to an increase of the cortical diffusion coefficient due to cortical activation. It is, however, not quite clear whether this enhanced dynamic is related to blood flow or to neural dynamics, or due to a combination of both. While the mere increase of blood volume due to activation cannot explain the observed increase in the cortical diffusion coefficient, vasodilation could lead to a large increase in blood flow velocity, even for small increases in blood volume due to the strong decrease in hydrodynamic resistance with the vessel diameter. Interestingly, the measured increase in the cortical diffusion coefficient of 38.6% for contralateral stimulation is comparable to positron emission tonography (PET) measurements^{51, 52} of the functional blood flow increase in motor cortical areas slightly smaller than the ones probed in our DWS experiments. In addition, the functional increase of $D_3$ measured in our experiments is in good agreement with the corrected relative blood flow rate increase measured by Durduran using DWS.

On the other hand, the high turbidity of cortical tissue indicates that DWS is sensitive not only to motions of erythrocytes, but also to motions of and within neura. At a hematocrit of 5%, bulk blood has a transport mean free path length $l^{*}=3.8\mathrm{mm}$ at $\lambda =800\mathrm{nm}$ .⁵³ Extrapolation to a physiological hematocrit of 50% and a cerebral blood volume of $4\mathrm{g}∕100\mathrm{m}\ell$ yields $l^{*}=9.6\mathrm{mm}$ , which is more than an order of magnitude larger than the transport mean free path length of typical cortical tissue with $l_3^{*}\approx 0.4\mathrm{mm}$ ,⁴⁴ which is dominated by the scattering from neura. Photons diffusing in cortical tissue should thus experience many more scattering events from neura than from erythrocytes. If both neural organelles and erythrocytes move mainly by Brownian motion (which is a reasonable approximation for the short length scales probed by DWS), the effective diffusion coefficient $D_3$ measured by DWS is the sum of the individual diffusion coefficients, each weighted with the reduced scattering coefficient of the respective component.⁴⁸ Assuming vesicles with diameter $200\mathrm{nm}$ , their contribution to $D_3$ is expected to dominate over the one from the erythrocytes at a vesicle volume fraction of about 18%. ⁵ Although this is a rough estimate, this value is consistent with values measured by electron microscopy in chromaffin cells.⁵⁵ While the contribution of subcellular organelle motion to the cortical diffusion coefficient could be investigated by local ischemia on arteriolar occlusion in animal models (similar to the experiments in Ref. 31), such experiments are not permissible in the human brain. However, laser Doppler measurements of biological zero signals (e.g., under cuff compression) in other types of tissue show that subcellular organelle motion not related to perfusion indeed contributes to the speckle fluctuations from perfused tissue.⁵⁶

Similar to the situation at rest, cortical diffusion coefficients measured during activation could be increased not only due to enhanced blood flow rate, but also due to increased mobility of synaptic vesicles associated with primary neural signals.^{57, 58, 59} While for the current motor protocol and data collection scheme the perfusion-related dynamics cannot be separated from primary neural signals to the short-time decay of the DWS autocorrelation function, one possible way to look for primary neural signals using DWS is by measuring the temporal intensity autocorrelation function with a shorter integration time. If the DWS signal is due to cerebral hemodynamics, the maximum difference of $D_3$ (compared with the baseline measurement) should appear with a delay of several seconds with respect to the onset of stimulation. In contrast, if synaptic vesicle motion also makes a significant contribution to the DWS signal, one should detect changes in the cortical diffusion coefficient earlier, from several tens to hundreds of milliseconds after the beginning of neural activation. To address this question experimentally, visual stimulation protocols with rapid visual serial presentation could be used.⁶⁰ The high temporal resolution of EEG and magnetoencephalography (MEG) revealed that the brain is capable of switching rapidly (within less than $100\mathrm{ms}$ ) between stimuli, recognizing the emotional character of each input picture. Obviously the cerebral blood flow cannot follow such rapid changes, so if another fast dynamical process, such as synaptic vesicle motion, also contributes significantly to the DWS signal, distinct signal changes for high-arousing versus low-arousing affective stimuli over parietal (preferably right-hemispheric) areas should be seen with this technique.

5.

Conclusions

We measure diffusion coefficients in the motor cortex of humans during motor stimulation, by analyzing near-infrared diffusing-wave spectroscopy data with the analytical solution of the correlation-diffusion equation for a three-layer geometry modeling the human head. The three-layer analysis yields good agreement with the experimental data over the entire range of lag times spanning more than three decades, and yields a quantitative estimate for the cortical diffusion coefficient. Analyzing the sensitivity of the model to uncertainties in tissue optical parameters, we find that the cortical diffusion coefficient can be extracted from our data to within a few percent if we allow tissue optical parameters to vary within realistic bounds. Performing motor stimulation experiments on 11 right-handed subjects, the autocorrelation function from the C3 region is found to decay faster during the stimulation period. This enhancement of the dynamics is shown to arise from changes in the cortical dynamics and not from changes in the optical properties of the tissue; our results thus support the interpretation of earlier DWS results for somato-motor stimulation in terms of functionally enhanced blood flow rate.³² The marked hemispheric asymmetry between the measured cortical diffusion coefficients in the C3 region for left- and right-hand stimulation is consistent with the dominant activation of the contralateral motor and sensory regions in right-handed subjects by unimanual motor tasks. While large hemodynamic changes dominate the cortical diffusion coefficient for the present motor stimulation protocol, detection of speckle fluctuations with shorter integration times might allow us to use DWS for the detection of primary neural signals associated with functionally enhanced motion of subcellular organelles, such as synaptic vesicles.