1.

Introduction

Skeletal muscles not only act as force transducers allowing for macroscopic displacements but also play an important role for the stabilization of vertebra and joints. Neuromuscular disorders such as those due to peripheral nerve injury can thus severely affect the overall systemic function and health condition. Fine-wire electromyography (EMG), the standard diagnostic technique in the clinic, allows to monitor the electrical activity of single motor units. However, the pain associated with electrode introduction and the sensitivity to spasm-induced artifacts prevent the method from being used for chronical monitoring of patients, e.g., during pharmaceutical treatment or rehabilitation. Tissue Doppler imaging, on the other hand, is a noninvasive method sensitive to tissue velocity, but it suffers from susceptibility to artifacts from global muscular displacements and relatively low sensitivity.¹ Early ex vivo experiments on frog muscle have suggested that quasi-elastic light scattering (QELS) is a sensitive probe for muscle contraction.² Attempts at using the laser Doppler method, the frequency-domain analog of QELS, for in vivo measurements have been shown to differentiate denervated muscles in the human hand.³ Unfortunately, the small Doppler widths of typically several $10\mathrm{Hz}\text{to}100\mathrm{Hz}$ that are dictated by the short source–receiver distance of a few mm do not allow to resolve the contraction process.⁴

Here, we show that the contraction of skeletal muscles can be measured noninvasively and with millisecond temporal resolution using near-infrared diffusing-wave spectroscopy (DWS). DWS (also called diffuse correlation spectroscopy, or DCS) is the extension of QELS to the regime of optically dense, opaque media.^{5, 6} Being highly sensitive to subwavelength displacements of scatterers within tissue, e.g., from erythrocytes moving within the vasculature, DWS has recently been introduced as a method for the noninvasive monitoring of blood flow in tumors^{7, 8, 9} and for the noninvasive measurement of functional brain activity.^{10, 11, 12, 13} DWS signals correlate strongly with tissue blood flow rates, as was demonstrated for skeletal muscle microvasculature.¹⁴ Its high sensitivity to tissue perfusion changes and the possibility for bedside measurements makes DWS an attractive method for monitoring blood flow, in particular in situations where conventional methods such as PET or Xe CT cannot be used, such as in neurointensive care, neonatology, or rehabilitation. Indeed, a recent study demonstrated that skeletal muscle perfusion and oxygenation can be monitored simultaneously with a compact hybrid instrument combining DWS and near-infrared spectroscopy.¹⁵ While it has recently been recognized that characterization of muscular blood flow by DWS during exercise does require correction for artifacts arising from muscle fiber motion,¹⁶ the potential of DWS for the explicit investigation of muscle contraction has so far not been explored.

Our results show that DWS signals from contracting skeletal muscles may be dominated by tissue shearing due to sarcomere shortening and advection of erythrocytes within the microvasculature. Using a multispeckle detection setup,¹⁷ we measure DWS autocorrelation functions with an integration time of $6.5\mathrm{ms}$ , which allows us to resolve the initial stages of the contraction and to differentiate between fast- and slow-twitch muscles. The analysis of DWS autocorrelation functions using tissue optical properties measured with time-resolved reflectance spectroscopy yields quantitative values for the evolution of the muscle strain $e\left(t\right)$ .

2.

Experiment

2.1.

Stimulation Protocol

All experiments were performed on a single, healthy, male volunteer (age $25\text{years}$ ). Muscle contraction was induced by applying bipolar voltage pulses ( $0.1\mathrm{ms}$ width per phase) generated by a function generator, at a frequency $1∕T=2.0\mathrm{Hz}$ (see Fig. 1 ), via two electrode pads $(48\times 68{\mathrm{mm}}^2)$ placed symmetrically about the muscle belly at a distance of about $10\mathrm{cm}$ . Measurements on the M. biceps brachii were carried out with the volunteer’s arm placed horizontally on a table at the height of the heart. Stimulation blocks lasted typically $160\mathrm{s}$ . Forces were measured with a piezoelectric transducer (Kistler 9321B) placed on the wrist. The force sensor was preloaded in order to guarantee good mechanical coupling during both contraction and relaxation. Recording the stimulus trigger signal allowed to synchronize the stimulation with the optical measurements. Optical measurements were started about $60\mathrm{s}$ after the start of a stimulation block, allowing the muscles to reach a physiological steady state. Heart beat rate was monitored with a pulse oxymeter at the index finger tip.

Fig. 1

Experimental setup and stimulation protocol. The excitation and detection paths for DWS and time-resolved reflectance spectrosocpy (TRS) are shown in blue and red, respectively. The optical sensors for DWS and TRS, consisting of source and receiver fibers, are oriented parallel to the muscle fiber and are firmly pressed to the skin near the belly of the muscle in order to minimize motion artifacts. DWS and TRS measurements, respectively, are triggered by the stimulus, which is controlled by the PC. TCSPC: time-correlated single photon counting card; PMT: photomultiplier tube.

3.

Results

3.1.

Force Measurements and Physiology

The forces measured on the wrist during contraction of the M. biceps brachii show, for a stimulation amplitude between $18\mathrm{V}$ and $22\mathrm{V}$ , a rapid increase to a peak and a slower decay during relaxation (see Fig. 2 ). The minimal stimulation voltage at which an increase in the force after stimulation can be observed varies between about $15\mathrm{V}$ and $23\mathrm{V}$ , depending on the mechanical contact between force transducer and the wrist. The latencies for the peak force vary between about $56\mathrm{ms}$ for low and $72\mathrm{ms}$ for the highest stimulation voltages, in good agreement with the literature value $(66\pm 9)\mathrm{ms}$ (Ref. 23. The slight undershoot of the force during relaxation at high stimulation amplitude is due to slight rotations of the wrist. Within experimental error, the heart beat rate measured with the pulse oxymeter did not show any difference between baseline and stimulation periods.

Fig. 2

Force $F\left(t\right)$ measured at the wrist as a function of time $t$ after stimulation of the human M. biceps brachii, measured for different stimulation voltages. The data are averages over 200 trials.

3.2.

Tissue Optical Properties

Figure 3 shows that representative single-trial data for $R\left(\tau \right)$ are well described by the theory for the homogeneous semi-infinite geometry. Upon electrical stimulation, the reduced scattering coefficient, ${\mu }_{\text{s}}^{\prime }$ , increases from its baseline value at about $(5.3\pm 0.1){\mathrm{cm}}^{-1}$ to a maximum whose height increases with stimulation voltage (see Fig. 4 ). For the lowest voltage at which a change in ${\mu }_{\text{s}}^{\prime }$ can be observed $\left(27.6\mathrm{V}\right)$ , ${\mu }_{\text{s}}^{\prime }$ reaches its maximum at $6.3{\mathrm{cm}}^{-1}$ , about 17% above the baseline value, within about $80\mathrm{ms}$ . The peak value of ${\mu }_{\text{s}}^{\prime }$ increases at higher stimulation amplitudes; the latency of the peak is, to within the temporal resolution of $20\mathrm{ms}$ , independent of stimulation amplitude.

Fig. 3

Single-trial time-resolved reflectance $R\left(\tau \right)$ measured on the M. biceps brachii at $\lambda =830\mathrm{nm}$ as a function of time $\tau$ (black). The dashed line is the fit of Eq. 1, including the instrument response function (gray), to the data. The width of the instrument response function of about $500\mathrm{ps}$ is due to pulse broadening in the optical fibers and timing jitter in the detector. Integration time: $20\mathrm{ms}$ . The fitting range is marked by the vertical lines. Best-fit values for ${\mu }_{\text{a}}$ and ${\mu }_{\text{s}}^{\prime }$ are $0.232{\mathrm{cm}}^{-1}$ and $5.38{\mathrm{cm}}^{-1}$ , respectively. The rms relative error between data and model is 4.1%.

Fig. 4

(a) Reduced scattering coefficient ${\mu }_{\text{s}}^{\prime }\left(t\right)$ and (b) absorption coefficient ${\mu }_{\text{a}}\left(t\right)$ at $\lambda =802\mathrm{nm}$ for different stimulation amplitudes, determined by linear interpolation of the TRS data measured at $690\mathrm{nm}$ and $830\mathrm{nm}$ . The data are averages over 200 stimuli. Orientation of the optical sensor is parallel to the muscle fiber.

In contrast to the behavior of ${\mu }_{\text{s}}^{\prime }$ , the absorption coefficient, ${\mu }_{\text{a}}$ , is observed to decrease rapidly from its baseline value of about $0.2{\mathrm{cm}}^{-1}$ to a minimum whose depth increases with the stimulation voltage (see Fig. 4). For the lowest stimulation voltage at which a change in ${\mu }_{\text{a}}$ is observed $\left(22.8\mathrm{V}\right)$ , the minimum value is by about 5% lower than the baseline value. The latency of the peak in ${\mu }_{\text{a}}$ is about $80\mathrm{ms}$ and decreases slightly with increasing stimulation amplitude.

The variations of the baseline values ${\mu }_{\text{s}0}^{\prime }$ and ${\mu }_{\text{a}0}$ of about 5.6% and 5.4%, respectively, over the range of stimulation amplitudes are likely due to variations in the optical coupling between the time-resolved reflectance spectroscopy (TRS) sensor and the tissue. We estimate that the average number of scattering events for ${\mu }_{\text{a}}\approx 0.2{\mathrm{cm}}^{-1}$ and $5.0{\mathrm{cm}}^{-1}⩽{\mu }_{\text{s}}^{\prime }⩽6.0{\mathrm{cm}}^{-1}$ varies between 45 and 60, which is consistent with the good agreement between measured $R\left(\tau \right)$ and the diffusion theory Eq. 1.

3.3.

Tissue Dynamics

Figure 5 shows that the average decay time ${\tau }_{\text{d}}\left(t\right)$ of the field autocorrelation function follows, in contrast to the evolution of ${\mu }_{\text{a}}\left(t\right)$ and ${\mu }_{\text{s}}^{\prime }\left(t\right)$ , a biphasic pattern. For the lowest stimulation voltage $\left(18.2\mathrm{V}\right)$ showing a functional signal, ${\tau }_{\text{d}}\left(t\right)$ drops within about $20\mathrm{ms}$ from the baseline value of about $80\mu \mathrm{s}$ to about $10\mu \mathrm{s}$ . At later times, ${\tau }_{\text{d}}\left(t\right)$ increases to a maximum whose height decreases with stimulation voltage. After this maximum, ${\tau }_{\text{d}}\left(t\right)$ recovers the baseline value via a second, more shallow minimum. The average decay times ${\tau }_{\text{d}}$ in the first and in the second minimum reach plateau values of about $5\mu \mathrm{s}$ and $12\mu \mathrm{s}$ for stimulation voltages exceeding $25\mathrm{V}$ . The position of the intermediate maximum in ${\tau }_{\text{d}}\left(t\right)$ at about $80\mathrm{ms}$ shows very weak dependence on the stimulation amplitude.

Fig. 5

Average decay time ${\tau }_{\text{d}}\left(t\right)=\int g^{\left(1\right)}(\tau ,t)\mathrm{d}\tau$ measured at the M. biceps brachii as a function of time $t$ after stimulation, for different stimulation amplitudes. The data are averages over 200 trials.

The reduced intensity autocorrelation functions $g^{\left(2\right)}(\tau ,t)-1$ show a pronounced nonexponential decay that, in accordance with the variations of ${\tau }_{\text{d}}\left(t\right)$ , strongly varies during the contraction–relaxation cycle (see Fig. 6 ).

Fig. 6

Reduced intensity autocorrelation functions $g^{\left(2\right)}(\tau ,t)-1$ as a function of lag time $\tau$ for different times $t$ after stimulus application. Error bars denote the standard deviation over 200 trials.

Mean square phase fluctuations $⟨\Delta {\phi }^2(\tau ,t)⟩$ extracted from measured $g^{\left(1\right)}(\tau ,t)$ cover the range of ${10}^{-3}$ to ${10}^{-1}$ for lag times $200\mathrm{ns}⩽\tau ⩽100\mu \mathrm{s}$ . Clearly, $⟨\Delta {\phi }^2(\tau ,t)⟩$ is not described by a single power law in $\tau$ , indicating that different mechanisms are responsible for the scatterer displacements (see Fig. 7 ).

Fig. 7

Mean square phase fluctuations $⟨\Delta {\phi }^2(\tau ,t)⟩$ as a function of lag time $\tau$ for different times $t$ after stimulus application. Dashed lines are linear least-squares fits of Eq. 2 to the data.

By fitting the mixed diffusion-shear model Eq. 2 to the measured mean square phase fluctuations, we obtain, using measured optical properties ${\mu }_{\text{a}}\left(t\right)$ and ${\mu }_{\text{s}}^{\prime }\left(t\right)$ in Eqs. 3, 4, time-dependent diffusion coefficients $D\left(t\right)$ , and the modulus of the strain rate $\left|\dot{e}\left(t\right)\right|$ (see Fig. 8 ). Both quantities follow a biphasic pattern similar to the one of ${\tau }_{\text{d}}\left(t\right)$ , indicating that the temporal evolution of the latter is dominated by temporal variations of the scatterer dynamics, while variations of the optical properties have a minor impact on the DWS signal. The two maxima of $\left|\dot{e}\left(t\right)\right|$ reflect the points of maximal elongation rate during contraction and relaxation, respectively, while the minimum separating the peaks reflects maximal contraction. The minimum of $\left|\dot{e}\left(t\right)\right|$ is reached at $t_0=83\mathrm{ms}$ for the lowest stimulation voltage showing a functional DWS signal $\left(18.2\mathrm{V}\right)$ , in very good agreement with the peak time of ${\tau }_d\left(t\right)$ . The peak amplitude of $\left|\dot{e}\left(t\right)\right|$ is observed to increase with stimulation voltage, likely due to an increased number of activated motor units.

Fig. 8

(a) Particle diffusion coefficient $D\left(t\right)$ and (b) modulus of the strain rate, $\left|\dot{e}\left(t\right)\right|$ , measured on the M. biceps brachii as a function of time $t$ after stimulation. A constant offset in $\left|\dot{e}\left(t\right)\right|$ was subtracted from the raw data. Inset: $\left|\dot{e}\left(t\right)\right|$ trace for stimulation at $13.6\mathrm{V}$ .

The tissue strain $e\left(t\right)$ is obtained using

Eq. 5

$$e\left(t\right)={\int }_{-0.2T}^t\mathrm{d}{t}^{\prime }\left|\dot{e}\left(t^{\prime }\right)\right|\left[\vartheta \right(t^{\prime }-0.2T\left)\vartheta \right(-t^{\prime })-\vartheta (t^{\prime }\left)\vartheta \right(t_0-t^{\prime })+\vartheta (t^{\prime }-t_0\left)\vartheta \right(0.8T-t^{\prime }\left)\right],$$

by where $\vartheta \left(t\right)$ is Heaviside's unit step function. The minus sign in the second term accounts for the fact that the DWS signal does not allow us to determine the sign of the strain rate. This procedure results in effective strains $e\left(t\right)$ with a single, asymmetric peak at the position $t_0$ of the minimum of $\left|\dot{e}\left(t\right)\right|$ (see Fig. 9 ). The height increases with stimulation amplitude, indicating the recruitment of an increasing number of muscle fibers. For the lowest supra-threshold stimulation amplitude of $18.2\mathrm{V}$ , we find that $e\left(t\right)$ returns, to within experimental error, back to the baseline level $e=0$ after relaxation. The peak latencies of $e\left(t\right)$ are practically independent of stimulation amplitude and agree very well with the ones of ${\mu }_{\text{a}}\left(t\right)$ and ${\mu }_{\text{s}}^{\prime }\left(t\right)$ .

Fig. 9

Effective strain $e\left(t\right)$ measured by DWS on the M. biceps brachii as a function of time $t$ after stimulation for different stimulation amplitudes. For stimulation voltages of $18.2\mathrm{V}$ , $22.8\mathrm{V}$ , and $27.6\mathrm{V}$ , the peak strains have a latency of $83\mathrm{ms}$ , $80\mathrm{ms}$ , and $77\mathrm{ms}$ , respectively.

4.

Discussion

DWS and TRS were measured to estimate tissue strains associated with contraction and relaxation of the biceps muscle. For the same source–receiver distance, stimulation-induced changes of ${\mu }_{\text{a}}$ and ${\mu }_{\text{s}}^{\prime }$ are observed for voltages beyond $22.8\mathrm{V}$ and $27.6\mathrm{V}$ , respectively. At a stimulation voltage of $27.6\mathrm{V}$ where both ${\mu }_{\text{s}}^{\prime }$ and ${\mu }_{\text{a}}$ first show a change with stimulation, the one of ${\mu }_{\text{s}}^{\prime }$ (17%) is considerably higher than the one of ${\mu }_{\text{a}}$ $(-5\mathrm{\%})$ . The increasing functional reduction of ${\mu }_{\text{a}}$ at higher stimulation amplitudes suggests reduction of vascular volume during the repetitive stimulation. Our results contrast the TRS data from $60\mathrm{s}$ plantar flexion of the M. gastrocnemius,²⁴ where ${\mu }_{\text{s}}^{\prime }\left(t\right)$ and ${\mu }_{\text{a}}\left(t\right)$ change in the opposite direction. The differences are likely due to the different orientation of the optodes in Ref. 24, and due to the muscle that was able to move under the loosely attached TRS probe.

In our experiments, we observed that the peaks of $\Delta {\mu }_{\text{a}}\left(t\right)$ and $\Delta {\mu }_{\text{s}}^{\prime }\left(t\right)$ reverse their sign when the optical sensor is oriented perpendicular to the muscle fiber. While the sign reversal of $\Delta {\mu }_{\text{s}}^{\prime }$ might reflect the anisotropy of light diffusion arising the structural anisotropy of the muscle fiber,²⁵ the sign reversal of $\Delta {\mu }_{\text{a}}$ could have various origins: (1) cross talk with ${\mu }_{\text{s}}^{\prime }$ related to our analysis of the TRS data with an optically isotropic model or (2) orientation-dependent weighting of contracting and noncontracting muscle regions characterized by reduced and increased vascular volume, respectively. However, a detailed elucidation of these issues is beyond the scope of this work. ³

Interestingly, the analysis of the mean square phase fluctuation including the measured optical properties is able to resolve the first peak of $\left|\dot{e}\left(t\right)\right|$ for a stimulation voltage of $13.6\mathrm{V}$ (see inset in Fig. 8). At this voltage, ${\mu }_{\text{a}}\left(t\right)$ and ${\mu }_{\text{s}}^{\prime }\left(t\right)$ are constant, and no functional DWS signal can be seen in ${\tau }_{\text{d}}\left(t\right)$ . This shows that the decomposition of $⟨\Delta {\phi }^2(\tau ,t)⟩$ into the diffusive and shear modes provides a sensitivity superior to that of ${\tau }_{\text{d}}\left(t\right)$ . The reason for this is that ${\tau }_{\text{d}}$ strongly weights the contributions of long lag times to $g^{\left(1\right)}(\tau ,t)$ that arise from dynamics within superficial layers. Indeed, a functional signal can be seen already at $13.6\mathrm{V}$ by computing $\int d\tau {\tau }^{\gamma }{g}^{\left(1\right)}(\tau ,t)$ , which, for $\gamma <0$ , provides a stronger weighting for the long photon paths. This increase of sensitivity provides direct evidence that DWS and TRS are indeed probing the muscle.

The observation that the mean square phase fluctuations $⟨\Delta {\phi }^2(\tau ,t)⟩$ can be explained only by including a term scaling quadratically with lag time $\tau$ during contraction and relaxation indicates that a substantial part of the DWS signal has its origin in deterministic scatterer displacements. In perfused striated muscle, spatial fluctuations of the dielectric constant leading to light scattering arise mainly from the periodic arrangement of sarcomeres and from the erythrocytes within the microvasculature. At $\lambda =800\mathrm{nm}$ , blood at a physiological hematocrit has ${\mu }_{\text{s}}^{\prime }\approx 20{\mathrm{cm}}^{-1}$ (Ref. 26). At a typical vascular volume fraction of 3% (Ref. 27) and the value ${\mu }_{\text{s}}^{\prime }\approx 5{\mathrm{cm}}^{-1}$ for muscle from our experiments, we estimate that the erythrocytes account for about 12% of the total number of scattering events. We expect that during the rapid contraction shortly after stimulation, the erythrocytes will follow the strain field within the muscle. By virtue of their higher scattering cross section, the erythrocytes within the microvasculature could thus enhance the contrast for the strain rate signal from the muscle.

The origin of the diffusive mode in $⟨\Delta {\phi }^2(\tau ,t)⟩$ that closely parallels the strain mode is not entirely clear. It has been noted earlier¹¹ that DWS autocorrelation functions from perfused tissue are described by diffusion rather than by random shear flow for which $⟨\Delta {\phi }^2\left(\tau \right)⟩\sim {\tau }^2$ . Whether the diffusion coefficient measured for muscle tissue is related to residual microvascular blood flow or whether it reflects other diffusive processes, such as crossbridge dynamics, remains to be clarified.

Recently, Yu reported DWS experiments on the M. gastrocnemius during a plantar flexion exercise where reduced DWS decay times were interpreted in terms of enhanced blood flow.²⁸ Our experiments using the same protocol, but with high temporal resolution show that the reduced decay time coincides with the transition from toe-up to toe-down, while the pulsatile background of ${\tau }_{\text{d}}$ due to the muscular blood flow is unaffected by the exercise (data not shown).

The question arises whether the increase of the peak effective strain $e\left(t_0\right)$ with stimulation voltage shown in Fig. 9 is due to increasing recruitment of motor units or due to increasing peak contraction, i.e., decreasing sarcomere length at maximal contraction. Looking at the force recordings in Fig. 2, we observe that the peak force increases with stimulation voltage, indicating enhanced motor unit recruitment resulting in increased number of contracting fibers within the volume probed by the DWS experiment. Note that in order to extract true volume-averaged muscle strains, the effective strains measured with DWS need to be normalized with the fraction of fibers that are actually contracting. For stimulation voltages above $27\mathrm{V}$ , we observe a saturation of the peak effective strain at $e\left(t_0\right)\approx 0.18$ , indicating the recruitment of all excitable fibers in the muscle. We thus estimate a value for the true volume-averaged strain of about 0.18, which is in good agreement with fiber length measurements on the cat gastrocnemius.²⁹

In order to explore whether DWS is able to discriminate between fast- and slow-twitch fibers, we performed experiments on the M. gastrocnemius of the same subject, using the same electrical stimulation protocol as in the experiments on the biceps muscle. The subject was in prone position on a table, and the DWS probe was positioned on the muscle belly. In contrast to the biceps muscle, which consists mainly of fast-twitch fibers, the gastrocnemius muscle predominantly contains slow-twitch fibers, which is reflected by a larger latency of the peak force of $(118\pm 7)\mathrm{ms}$ (Ref. 30). The temporal evolution of the DWS decay time ${\tau }_{\text{d}}\left(t\right)$ is qualitatively similar to the one in the M. biceps brachii; in contrast, the latency of $96.7\mathrm{ms}$ of the peak in ${\tau }_{\text{d}}$ is significantly larger (see Fig. 10 ). Clearly, this difference is not explained by possible differences in excited muscle volume between biceps and gastrocnemius muscle, but rather suggests that it reflects the physiological difference between fast- and slow-twitch fibers. More detailed experiments are required to verify this hypothesis.

Fig. 10

Evolution of the DWS decay time ${\tau }_{\text{d}}\left(t\right)$ measured on the M. biceps brachii and on the M. gastrocnemius. The stimulation amplitude was $44.8\mathrm{V}$ for both measurements. Note the large difference in the latency of the intermediate maximum corresponding to maximal contraction (denoted by the dashed lines) in the two muscles.

5.

Conclusions

In conclusion, we show on a single subject that DWS sensitively and noninvasively measures skeletal muscle contraction with millisecond temporal resolution. Single-photon counting and multispeckle detection allow to place source and receiver at a distance that is sufficient to overcome adipose tissue layers several millimeters thick. Further improvement of signal-to-noise ratio and spatial resolution could thus open the way for quantitative, noninvasive characterization of muscle function in sports medicine and clinical applications.